geoschem / geos-chem

GEOS-Chem "Science Codebase" repository. Contains GEOS-Chem science routines, run directory generation scripts, and interface code. This repository is used as a submodule within the GCClassic and GCHP wrappers, as well as in other modeling contexts (external ESMs).
http://geos-chem.org
Other
167 stars 162 forks source link

[QUESTION] question about sulfur chemistry via KPP #1139

Closed lumosmoon closed 2 years ago

lumosmoon commented 2 years ago

Hi Geos-Chem team,

I noticed that you have implemented sulfur in-cloud and on aerosols reactions via KPP. I have some questions after reviewing the code.

  1. Both produce SO4, so why are some reactants written as SO2 (e.g., SO2+ H2O2 = SO4 )while others are written as HSO3m/SO3mm (e.g., HOCl + HSO3m = SO4 + HCl)? Is it true that the latter reaction produces SO4 without consuming SO2? It instead consumes HSO3m. This brings me to my second question.

  2. It seems that HSO3m/SO3mm are calculated in SUBROUTINE AQCHEM_SO2 in KPP/fullchem/fullchem_SulfurChemFuncs.F90. At each timestep, they are calculated from the dissolving equilibrium of SO2, which does not appear to be affected by its consumption in KPP?

Any help you could provide would be much appreciated

yantosca commented 2 years ago

Thanks for writing. I am tagging @msl3v, who can give a better explanation than I can.

msl3v commented 2 years ago

Hello lumosmoon, See below...

  1. Both produce SO4, so why are some reactants written as SO2 (e.g., |SO2+ H2O2 = SO4| )while others are written as HSO3m/SO3mm (e.g., |HOCl + HSO3m = SO4 + HCl|)? Is it true that the latter reaction produces SO4 without consuming SO2? It instead consumes HSO3m. This brings me to my second question.

The original code for hypohalous acids HOBr and HOCl branched the reactions into two, accounting for HSO3m and SO3mm pathways, though it did so a little differently than done here. (Note you'll find the same operations in sulfate_mod.F90 that you find in fullchem_SulfurChemFuncs.F90) The goal with the transition to KPP was to maintain as much of the original structure as possible and only move the chemical reaction computation into KPP. That is reflected here. I wanted to put HSO3m and SO3mm in explicitly to lay the groundwork for proper dissolved-phase partitioning when/if it becomes possible or necessary. I hope it's not too confusing.

  1. It seems that HSO3m/SO3mm are calculated in |SUBROUTINE AQCHEM_SO2| in |KPP/fullchem/fullchem_SulfurChemFuncs.F90|. At each timestep, they are calculated from the dissolving equilibrium of SO2, which does not appear to be affected by its consumption in KPP?

This was also meant to be consistent with the previous implementation. You point out a real potential issue here regarding accounting and mass balance. The way I justified keeping the original behavior was that consumption of any dissolved S(IV) would in any case be replenished by more uptake from gas phase, and that gas/particle would always be in equilibrium. If we can assume that little of the S(IV) dissolved in a timestep remains as as S(IV), then we can also assume mass balance error is minimized. Unless I am overlooking anything in the sulfur code, there is no way for dissolved S(IV) to return to the gas phase (right?).

Ultimately, the sum of both of your questions seems to highlight some inconsistency in the sulfur code with regards to gas/particle partitioning. If GCSC or GCST don't like this implementation, it is possible to make the reactions more like the H2O2 pathway. The HSO3m and SO3mm partitioning calculation would still have to be performed in either case.

I hope I answered your questions. Let me know if you have any more, or if you see any options for improvement. I'm eager to see this part of the code get more scrutiny.

Sincerely, Michael Long

lumosmoon commented 2 years ago

Hello @msl3v . Thanks a lot for your detailed response. I still have some questions.

SO2+H2O2=SO4 can deduct SO2 directly and the remaining SO2 will be passed to the next time step; but when it comes to HOCl+HSO3m=SO4+HCl, although HSO3m is consumed, the code does not see how the effect of HSO3m consumption is passed to the next time step.

As at the next time step, HSO3m is calculated by the gas-liquid rebalance of SO2 using this formula (see below). Unless the consumption of HSO3m/SO3mm also affects the concentration of hydrogen ions, which is linked to the thermodynamic process (ISORROPIA)? image

We use pSO2(the partial pressure of SO2) and H+ to calculate HSO3m/SO3mm, but we do not subtract these parts from SO2(g), right? So if we write other reactions e.g. H2O2+ HSO3m/SO3mm = SO4similar to HOCl, will the SO2(g) concentration never drop due to chemistry? Perhaps I am missing something in the code.

Maybe some tests could be done: ①HOCl + HSO3m = SO4 + HClHOCl + SO2 = SO4 + HCl. Then we can see how much impact these treatments would have on SO2 and sulfate abundance.

Thanks for your kind help.

msl3v commented 2 years ago

@lumosmoon Since your original questions, the decision was made to rewrite the reactions to use SO2 directly. The HSO3m/SO3mm calculation will only be used as a scaling coefficient in the rate constant calculation. This should resolve your concerns. Please let me know if you have others.

On 3/1/22 07:10, lumosmoon wrote:

Hello @msl3v https://github.com/msl3v . Thanks a lot for your detailed response. I still have some questions.

|SO2+H2O2=SO4| can deduct SO2 directly and the remaining SO2 will be passed to the next time step; but when it comes to |HOCl+HSO3m=SO4+HCl|, although HSO3m is consumed, the code does not see how the effect of HSO3m consumption is passed to the next time step.

As at the next time step, HSO3m is calculated by the gas-liquid rebalance of SO2 using this formula (see below). Unless the consumption of HSO3m/SO3mm also affects the concentration of hydrogen ions, which is linked to the thermodynamic process (ISORROPIA)? image https://user-images.githubusercontent.com/61012072/156165553-238d5d76-635e-44c2-89da-94d58c9c30da.png

We use pSO2(the partial pressure of SO2) and H+ to calculate HSO3m/SO3mm, but we do not subtract these parts from SO2(g), right? So if we write other reactions |e.g. H2O2+ HSO3m/SO3mm = SO4 |similar to HOCl, will the SO2(g) concentration never drop due to chemistry? Perhaps I am missing something in the code.

Maybe some tests could be done: ①|HOCl + HSO3m = SO4 + HCl| ②|HOCl + SO2 = SO4 + HCl|. Then we can see how much impact these treatments would have on SO2 and sulfate abundance.

Thanks for your kind help.

— Reply to this email directly, view it on GitHub https://github.com/geoschem/geos-chem/issues/1139#issuecomment-1055377983, or unsubscribe https://github.com/notifications/unsubscribe-auth/ABVE23ZXPNQC5AFVKEJ7GBLU5YCKDANCNFSM5PJVNNWA. Triage notifications on the go with GitHub Mobile for iOS https://apps.apple.com/app/apple-store/id1477376905?ct=notification-email&mt=8&pt=524675 or Android https://play.google.com/store/apps/details?id=com.github.android&referrer=utm_campaign%3Dnotification-email%26utm_medium%3Demail%26utm_source%3Dgithub.

You are receiving this because you were mentioned.Message ID: @.***>

djxjacob commented 2 years ago

Mike, will you make that change in 13.4 or does it require more testing? As you know, we need to redo the benchmarking of 13.4 to disable SSA debromination so there’s an opportunity to slip it in there. Daniel

msl3v commented 2 years ago

The code is done, and I'll do some short-term testing on it this week. Otherwise, it's not a complex change, and it can definitely make it into the re-benchmark.

msulprizio commented 2 years ago

@msl3v Is there a status update on this fix? Can you submit a PR when it's ready? As Daniel mentioned above, we'd like to include this in 13.4.0, and I hope to run the benchmarks for that version (again) next week.

msulprizio commented 2 years ago

The corresponding PR for this update is https://github.com/geoschem/geos-chem/pull/1170.

msulprizio commented 2 years ago

The update described above has now been merged into dev and will be included in GEOS-Chem 13.4.0.