Short-lived halogen chemistry by fvitt · Pull Request #1470 · ESCOMP/CAM

fvitt · 2026-01-07T15:54:26Z

new file: cime_config/testdefs/testmods_dirs/cam/ts_slh_mam4_outfrq9s/shell_commands new file: cime_config/testdefs/testmods_dirs/cam/ts_slh_mam4_outfrq9s/user_nl_cam new file: cime_config/testdefs/testmods_dirs/cam/ts_slh_mam4_outfrq9s/user_nl_clm new file: cime_config/testdefs/testmods_dirs/cam/wa_slh_mam4_outfrq9s/shell_commands new file: cime_config/testdefs/testmods_dirs/cam/wa_slh_mam4_outfrq9s/user_nl_cam new file: cime_config/testdefs/testmods_dirs/cam/wa_slh_mam4_outfrq9s/user_nl_clm

fvitt · 2026-01-09T19:47:42Z

From @RafaPedroFernandez (CC'd @CAlbertoCuevas):

Please note that I have just made a push for 2 F90 routines in my following GIT:
https://github.com/RafaPedroFernandez/CAM/tree/slh_chemistry_n04_cam6_4_097_rpf

Which is also updated local in derecho:
/glade/u/home/rpfernan/cesm_versions_rpf/derecho_cesm/git_fetch/cesm3_cam6dev_francis/slh_chemistry_n04_cam6_4_097_rpf

Routines are
chemistry.F90:
For SLH, we must call iodine_emissions_srf before the loop to outfld( sflxnam(m), cam_in%cflx(:ncol,m), ncol,lchnk ) in order for SFHOI and SFI2 to be properly written in the output file. If we do not so, the SFHOI and SFI2 are zero in the output file but are somehow included in cam_in%cflx array that is passed to other routines

mo_gas_phase_chmdr.F90
For Iodine washout ... we have always used SAD_ICETROP, which is computed based on CLDICE. This CLDICE was before computed as only ice
cldice(:ncol,k) = q(:ncol,k,cldice_ndx)
and now considers also SNOW
cldice(:ncol,k) = q(:ncol,k,cldice_ndx) + q(:ncol,k,snow_ndx)

Including snow in the variable makes the iodine deposition too efficient, and then we created a new array (ONLYICE) just to compute iodine wet-depositions without affecting other model developments.

modified: src/chemistry/mozart/chemistry.F90 modified: src/chemistry/mozart/mo_gas_phase_chemdr.F90

modified: src/chemistry/mozart/mo_usrrxt.F90

cacraigucar

Partial review - @jimmielin will be reviewing the chemistry section of this PR. My review was on the setup side.

I noticed that there are testing use_cases setup, but I do not see any regression tests being added. Is this an oversight?

cacraigucar · 2026-01-15T23:59:14Z

bld/namelist_files/namelist_defaults_cam.xml

+<xs_short_file>atm/waccm/phot/xs_short_c221005mm_slh_c250930.nc</xs_short_file>
+<xs_long_file >atm/waccm/phot/temp_prs_GT200nm_c221005mm_slh_c250930.nc</xs_long_file>


These two file still need to be "rimport"ed

@fvitt
The files are located at:
/glade/campaign/cesm/cesmdata/cseg/inputdata/atm/waccm/phot
I understand that both files should be moved to
/glade/campaign/collections/gdex/data/d651077/cesmdata/inputdata/atm/waccm/phot
and then include a symbolic link (ln -s) in the original path.

Francis, I assume you are the responsible for doing so

@cacraigucar We have now other new inputs. Is using the rimport script sufficient for adding new inputs to the input data repo?

Yes, the script performs the copy and link steps. As it has always been, once a file is moved there, it exists "forever", so files should only be rimported once they are "cast-in-stone".

cacraigucar · 2026-01-16T00:00:52Z

bld/namelist_files/namelist_definition.xml

+<entry id="SSAhno3_ScalingFactor" type="real" category="cam_chem"
+       group="slh_nl" valid_values="" >
+Scaling-factor for SLH Halocarbon surface emissions
+Default: set by build-namelist
+</entry>
+
+<entry id="SSAn2o5_ScalingFactor" type="real" category="cam_chem"
+       group="slh_nl" valid_values="" >
+Scaling-factor for SLH Halocarbon surface emissions
+Default: set by build-namelist
+</entry>
+
+<entry id="LIQfraprx_ScalingFactor_I" type="real" category="cam_chem"
+       group="slh_nl" valid_values="" >
+Scaling-factor for SLH Halocarbon surface emissions
+Default: set by build-namelist
+</entry>
+
+<entry id="ICEfraprx_ScalingFactor_I" type="real" category="cam_chem"
+       group="slh_nl" valid_values="" >
+Scaling-factor for SLH Halocarbon surface emissions
+Default: set by build-namelist
+</entry>
+
+<entry id="ICEfraprx_ScalingFactor_Br" type="real" category="cam_chem"
+       group="slh_nl" valid_values="" >
+Scaling-factor for SLH Halocarbon surface emissions
+Default: set by build-namelist
+</entry>


Comments probably need updating (since they are all the same)

@fvitt
Please update the "long-name" description as follows (highlighted in bold):
(Note that there are 6 Scaling Factors ... I added a complete new block for SSAdehal_ScalingFactor that should be just a few lines above in the same file)

SSAdehal_ScalingFactor
Scaling-factor for SLH dehalogenation reaction (XONO2,XNO2,HOX->X2,XY ... X,Y=Br,Cl,I) on Sea-Salt Aerosols
Default: set by build-namelist

Scaling-factor for SLH acid displacement reaction (HNO3->HCl) on Sea-Salt Aerosols
Default: set by build-namelist

Scaling-factor for SLH dinitrogen pentoxide reaction (N2O5->HNO3+CLNO2) on Sea-Salt Aerosols
Default: set by build-namelist

Scaling-factor for SLH washout (Iodine) on liquid droplets based on Free-Regime Approx.
Default: set by build-namelist

Scaling-factor for SLH washout (Iodine) on ice crystals based on Free-Regime Approx.
Default: set by build-namelist

Scaling-factor for SLH washout (Bromine) on ice crystals based on Free-Regime Approx.
Default: set by build-namelist

@RafaPedroFernandez
See my updates to these descriptions here: fvitt@b060e32

cacraigucar · 2026-01-16T00:05:56Z

cime_config/config_compsets.xml

 	<value  compset="HIST_CAM60%CFIRE">1995-01-01</value>
 	<value  compset="RCP[2468]_CAM\d+">2005-01-01</value>
- 	<value  compset="_CAM.*%NUDG"                   >2010-01-01</value>
+ 	<value  compset="_CAM.*%NUDG"                   >2000-01-01</value>


Why is this PR changing the date for nudging runs?

When I initially ported the code I just copy/paste the most similar PR ... and I shifted the date to 2000 because in the description paper we only validated SLH simulations in the nudgged setup for the period 2000-2005.

I suggest keeping the default value use for year 2010 and in any case the user can modify this using user_nl_cam.

Reverted this back

jimmielin

Thanks @fvitt, I have a few suggested changes but nothing major.

jimmielin · 2026-02-03T16:10:46Z

src/chemistry/mozart/clybryiy_fam.F90

+       wrk(:,:) = cloy( mmr, pcols, ncol )
+
+       factor(:ncol,:) = mmr(:ncol,:,id_cly) / wrk(:ncol,:)


Maybe a safeguard against divide by zero here if wrk(:ncol,:) are close to zero (when cloy initial conditions are very low?)

I do not know what to say here about dividing by zero. Francis, do you have any suggestion? Note I just copy/paste and updated those lines as they were initially implemented in the old clybry.F90 file when Iodine was not considered. And we never found an issue with this.

In my oppinion, the most important thing to control here are lines 101-116: In case the chem_mech.in DOES NOT include IODINE, then the "original" ids_clybry species should be considered ... so any user working with a mechanism that does not include SLH can compute the total Cl and Br abundance properly.

!rpf_CESM2_SLH
ids_clybry = (/ id_cly,id_bry, &
id_cl,id_clo,id_hocl,id_cl2,id_cl2o2,id_oclo,id_hcl,id_clono2, &
id_br,id_bro,id_hbr,id_brono2,id_brcl,id_hobr /)

ids_clybryiy = (/ id_cly,id_bry,id_iy, & id_cl,id_clo,id_hocl,id_cl2,id_cl2o2,id_oclo,id_hcl,id_clono2,id_clno2, & id_br,id_bro,id_hbr,id_brono2,id_brcl,id_hobr,id_br2,id_brno2, & id_i,id_i2,id_io,id_oio,id_hi,id_hoi,id_ino,id_ino2,id_iono2,id_ibr,id_icl,id_i2o2,id_i2o3,id_i2o4 /) if ( id_iy>0 ) then has_clybryiy = all( ids_clybryiy(:) > 0 ) else has_clybryiy = all( ids_clybry(:) > 0 ) endif

!rpf_CESM2_SLH

Thanks @RafaPedroFernandez, looks like @fvitt added a check

if ( id_cly>0 ) then wrk(:,:) = cloy( mmr, pcols, ncol ) where( wrk(:ncol,:)>0._r8 ) factor(:ncol,:) = mmr(:ncol,:,id_cly) / wrk(:ncol,:) elsewhere factor(:ncol,:) = 0._r8 end where

I am happy with this after the added checks

@fvitt
Francis ... I believe that in the extrange case thet wrk = 0.0 (we never found that issue) ... then factor should be equal to 1, and not to zero. Don't you think?
elsewhere
factor(:ncol,:) = 1._r8
end where
I always understood this portion of the code as a "double checking" that the sum of the advected individual Cly species was equivalent to the CLOY variable.

In either case, please note that in the same routine, there are also equivalent "factor" computations for id_bry (line 257 of the original file) and id_iy (line 278). All three elsewhere conditions should be implemented for Chlorine, Bromine and Iodine.

src/chemistry/mozart/mo_slh_routines.F90

jimmielin · 2026-02-03T16:12:24Z

src/chemistry/mozart/mo_slh_routines.F90

+!rpf_CESM2_SLH
+!rpf and cac: This threshold is required to aviod huge iodine emissions for low wind peed
+    where ( ws_mask(:ncol) < 3._r8 )
+      ws_mask(:ncol) = 3._r8
+    endwhere
+!rpf_CESM2_SLH


Maybe clean up this comment a bit and fix typos:

Suggested change

!rpf_CESM2_SLH

!rpf and cac: This threshold is required to aviod huge iodine emissions for low wind peed

where ( ws_mask(:ncol) < 3._r8 )

ws_mask(:ncol) = 3._r8

endwhere

!rpf_CESM2_SLH

!rpf and cac: This threshold is required to avoid huge iodine emissions for low wind speeds

where ( ws_mask(:ncol) < 3._r8 )

ws_mask(:ncol) = 3._r8

endwhere

Thanks for finding the typo. I agree with the change. Please replace "huge" by "large"

jimmielin · 2026-02-03T16:14:08Z

src/chemistry/mozart/mo_slh_routines.F90

+       !!! Be CAREFULL that if SST = 0 then the division is not possible. Should it be forced to zero?
+       if (ocnfrac(i) == 0) then ! Actually it should be landfrac = 1 ... but we do not want to bring landfrac as an additional agrument
+          iodide_aq(i) = 0 
+       else
+!         iodide_aq(i) = 1.46e6_r8 * exp( -9134._r8 / sst(i) ) * ocnfrac(i)
+          iodide_aq(i) = 1.46e6_r8 * exp( -9134._r8 / sst(i) )
+       end if


Since SST is in Kelvin it shouldn't get close to zero but I would check if it is intended that iodine_aq is computed when sst is below freezing? This may be a science question

You are right ... SST is in Kelvin and therefore no division by zero will occur here.
The "scientific" condigion about considereing or not iodide_aq fields for iodine emissions is considered above in the logical condition "if (ocnfrac(i) == 0) then"

@RafaPedroFernandez Is this block of code okay as it is? Should the "Be CAREFULL ..." be removed?

@fvitt
Yes francis. You can remove the "Be carefull ..." line completely.

In addition, please also notice that we keep the following line 272 commented as this was the original computation in CESM1.
! iodide_aq(i) = 1.46e6_r8 * exp( -9134._r8 / sst(i) ) * ocnfrac(i)

But while porting the code to CESM2, we recognized that the multiplication by ocnfrac(i) was incorrect and replace it by line 273
iodide_aq(i) = 1.46e6_r8 * exp( -9134._r8 / sst(i) )

You know better than me if both (commented and uncommented) lines should remain in the trunk for "trazability" of the implementation. In case it is better to keep only one, please remove line 272 (commented) and keep only line 273

src/chemistry/mozart/mo_slh_routines.F90

src/chemistry/mozart/mo_strato_rates.F90

jimmielin · 2026-02-03T16:20:53Z

src/chemistry/mozart/mo_strato_rates.F90

+                     if ( rid_slf_str_i_2 > 0 ) then
+                        if( hclvmr > small_div .and. hoivmr > small_div ) then
+! CAC 2020 RPF START OF implmentation of the HOI + HCl(liq) =  ICl + H2O  Sulfate Aerosol Reaction
+! NOTE: gprob for HOI species is mapped to gprob_hobr_hcl 
+                           if ( rid_slf_str_i_2 > 0 ) then


The check in line 918 is a duplicate of the check in line 914 for rid_slf_str_i_2 > 0.

You are correct. This is a duplicate condition and the second one could be removed.
Note the same happens for rid_slf_str_i_3 in line 932 and also for rid_slf_str_i_4 in line 948.

Idem for rid_nat_str_i_2, rid_nat_str_i_3, rid_nat_str_i_4 in lines 1131, 1147, 1163.

Idem for rid_ice_str_i_2, rid_ice_str_i_3, rid_ice_str_i_4 in lines 1337, 1353, 1369.

jimmielin · 2026-02-03T16:21:02Z

src/chemistry/mozart/mo_strato_rates.F90

+                     if ( rid_slf_str_i_3 > 0 ) then
+                        if( hbrvmr > small_div .and. hoivmr > small_div ) then
+! CAC 2020 RPF START OF implmentation of the HOI + HBr(liq) =  IBr + H2O  Sulfate Aerosol Reaction
+! NOTE: gprob for HOI species is mapped to gprob_hobr_hcl 
+                           if ( rid_slf_str_i_3 > 0 ) then


Same duplicate checks here

jimmielin · 2026-02-03T16:21:11Z

src/chemistry/mozart/mo_strato_rates.F90

+                     if ( rid_slf_str_i_4 > 0 ) then
+                        if( hivmr > small_div .and. hoivmr > small_div ) then
+! CAC 2020 RPF START OF implmentation of the HOI + HI (liq) =  I2  + H2O  Sulfate Aerosol Reaction
+! NOTE: gprob for HOI species is mapped to gprob_hobr_hcl 
+                           if ( rid_slf_str_i_4 > 0 ) then


Same duplicate checks here

jimmielin · 2026-02-03T16:23:51Z

src/chemistry/mozart/mo_usrrxt.F90

+! Is there a typo in the pre-exponential factor? (1.72 instead of 1.73?) Which one is correct?
+!            rxt(:,k,usr_NO2XNO3_M_ndx) = rxt(:,k,tag_NO2_NO3_ndx) * 1.734138e26_r8 * exp_fac(:)
+             rxt(:,k,usr_NO2XNO3_M_ndx) = rxt(:,k,tag_NO2_NO3_ndx) * 1.724138e26_r8 * exp_fac(:)
+!rpf_CESM2_SLH


This comment is leftover and should be checked for correctness

@RafaPedroFernandez Can you comment on this? Do you know what the correct value is?

@fvitt @lkemmons @tilmes
Francis. I have not developed not used those reactions for our SLH scheme. I just noticed that the 1.734138e26_r8 value used in original line 1282 was different to the equivalent expresions used in lines 1264 ad 1272 for equivalent reactions (that used the value 1.724138e26_r8). As numbers are very similar and expresions are equivalent I thought this could be an unnoticed typo. I believe Louisa or Simone can tell you which is the correct value (or who should we ask for).

@fvitt Please correct this to 1.724.
Thank you, Rafa, for finding this.

Thanks for this. Cleaned up the comments here and accepted the correction.

modified: bld/namelist_files/namelist_defaults_cam.xml modified: bld/namelist_files/use_cases/2000_trop_strat_cam6_slh.xml modified: bld/namelist_files/use_cases/hist_trop_strat_cam6_slh.xml modified: bld/namelist_files/use_cases/hist_trop_strat_nudged_cam6_slh.xml modified: bld/namelist_files/use_cases/waccm_tsmlt_2000_cam6_slh.xml modified: bld/namelist_files/use_cases/waccm_tsmlt_hist_cam6_slh.xml modified: bld/namelist_files/use_cases/waccm_tsmlt_nudged_cam6_slh.xml

modified: src/chemistry/mozart/clybryiy_fam.F90 modified: src/chemistry/mozart/mo_slh_routines.F90 modified: src/chemistry/mozart/mo_srf_emissions.F90 modified: src/chemistry/mozart/mo_strato_rates.F90

modified: cime_config/config_compsets.xml modified: src/chemistry/mozart/mo_usrrxt.F90

fvitt added 2 commits January 6, 2026 16:06

mods from slh_chemistry branch

628491f

fvitt self-assigned this Jan 7, 2026

cacraigucar requested review from cacraigucar and jimmielin January 7, 2026 17:39

fvitt marked this pull request as draft January 8, 2026 00:09

more changes from Rafa

55f61c9

modified: src/chemistry/mozart/chemistry.F90 modified: src/chemistry/mozart/mo_gas_phase_chemdr.F90

fvitt marked this pull request as ready for review January 13, 2026 23:11

add r8 real kind specifiers

48d7fbd

modified: src/chemistry/mozart/mo_usrrxt.F90

cacraigucar added this to CAM Development Jan 15, 2026

cacraigucar requested changes Jan 16, 2026

View reviewed changes

jimmielin requested changes Feb 3, 2026

View reviewed changes

fvitt added 4 commits February 10, 2026 19:13

code clean up

7209e28

modified: src/chemistry/mozart/clybryiy_fam.F90 modified: src/chemistry/mozart/mo_slh_routines.F90 modified: src/chemistry/mozart/mo_srf_emissions.F90 modified: src/chemistry/mozart/mo_strato_rates.F90

update discriptions of SLH scaling namelist params

b060e32

revert start data; clean up comments in usrxt (accept correction)

0144e21

modified: cime_config/config_compsets.xml modified: src/chemistry/mozart/mo_usrrxt.F90

		<xs_short_file>atm/waccm/phot/xs_short_c221005mm_slh_c250930.nc</xs_short_file>
		<xs_long_file >atm/waccm/phot/temp_prs_GT200nm_c221005mm_slh_c250930.nc</xs_long_file>

		wrk(:,:) = cloy( mmr, pcols, ncol )

		factor(:ncol,:) = mmr(:ncol,:,id_cly) / wrk(:ncol,:)

Conversation

fvitt commented Jan 7, 2026

Uh oh!

fvitt commented Jan 9, 2026

Uh oh!

cacraigucar left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

RafaPedroFernandez Feb 12, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

RafaPedroFernandez Feb 12, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

jimmielin left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Uh oh!

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

RafaPedroFernandez Feb 12, 2026 •

edited

Loading

RafaPedroFernandez Feb 12, 2026 •

edited

Loading