diff --git a/.gitmodules b/.gitmodules index ffb20365f1..c6ca4f5c77 100644 --- a/.gitmodules +++ b/.gitmodules @@ -51,13 +51,27 @@ path = src/chemistry/geoschem/geoschem_src url = https://github.com/geoschem/geos-chem.git fxrequired = AlwaysRequired - fxtag = 14.1.2 + fxtag = 14.4.3 fxDONOTUSEurl = https://github.com/geoschem/geos-chem.git +[submodule "cloud_j"] + path = src/chemistry/cloud_j + url = https://github.com/geoschem/cloud-j.git + fxrequired = AlwaysRequired + fxtag = 7.7.3 + fxDONOTUSEurl = https://github.com/geoschem/cloud-j.git + +[submodule "hetp"] + path = src/chemistry/hetp + url = https://github.com/geoschem/heterogeneous-vectorized-or-parallel.git + fxrequired = AlwaysRequired + fxtag = geoschem_hetp_1.0 + fxDONOTUSEurl = https://github.com/geoschem/heterogeneous-vectorized-or-parallel.git + [submodule "hemco"] path = src/hemco url = https://github.com/ESCOMP/HEMCO_CESM.git - fxtag = hemco-cesm1_2_1_hemco3_6_3_cesm_rme01 + fxtag = hemco-cesm2_0_hemco3_9_0 fxrequired = AlwaysRequired fxDONOTUSEurl = https://github.com/ESCOMP/HEMCO_CESM.git diff --git a/bld/build-namelist b/bld/build-namelist index d0a1c8c8d0..66c3574a62 100755 --- a/bld/build-namelist +++ b/bld/build-namelist @@ -2038,10 +2038,13 @@ if (($chem =~ /trop_mozart/ or $chem =~ /trop_strat/ or $chem =~ /waccm_tsmlt/) if ($chem =~ /geoschem/) { - my $val; + # Input directories + add_default($nl, 'geoschem_chem_inputs'); + add_default($nl, 'geoschem_aeropt_inputs'); + add_default($nl, 'geoschem_photol_inputs'); # Species with fixed lower boundary - $val = "'CCL4','CH4','N2O','CO2','CFC11','CFC12','CH3BR','CH3CCL3','CH3CL'" + my $val = "'CCL4','CH4','N2O','CO2','CFC11','CFC12','CH3BR','CH3CCL3','CH3CL'" .",'HCFC22','CFC114','CFC115','HCFC141B','HCFC142B','CH2BR2','CHBR3','H2402'"; if ($chem_has_ocs) { diff --git a/bld/configure b/bld/configure index 5168b0402c..e64e78b6a6 100755 --- a/bld/configure +++ b/bld/configure @@ -1450,9 +1450,9 @@ if ($chem_pkg =~ '_mam3') { # Customize GEOS-Chem advected species and chemistry CPP definitions if ($chem_pkg =~ 'geoschem') { - $chem_cppdefs .= ' -DEXTERNAL_GRID -DEXTERNAL_FORCING '; + $chem_cppdefs .= ' -DEXTERNAL_GRID -DEXTERNAL_FORCING -DMODEL_CESM -DMODEL_GEOSCHEM'; if ($chem_pkg =~ '_mam4') { - $chem_nadv = 267; # includes GC advected species (233), CO2 (1), and MAM aerosols (33) + $chem_nadv = 269; # includes all GEOS-Chem advected species, CO2, and MAM aerosols } } # CARMA sectional microphysics @@ -2219,8 +2219,10 @@ sub write_filepath print $fh "$chem_src_dir\n"; } - # GEOS-Chem must be prior to Mozart + # GEOS-Chem and dependencies; GEOS-Chem must be prior to Mozart if ($chem_pkg =~ 'geoschem') { + print $fh "$camsrcdir/src/chemistry/cloud_j/src/Core\n"; + print $fh "$camsrcdir/src/chemistry/hetp/src/Core\n"; print $fh "$chem_src_dir/geoschem_src/GeosCore\n"; print $fh "$chem_src_dir/geoschem_src/GeosUtil\n"; print $fh "$chem_src_dir/geoschem_src/Headers\n"; diff --git a/bld/namelist_files/geoschem_master_gas_drydep_list.xml b/bld/namelist_files/geoschem_master_gas_drydep_list.xml index eebafa33a7..dcf1ed67f8 100644 --- a/bld/namelist_files/geoschem_master_gas_drydep_list.xml +++ b/bld/namelist_files/geoschem_master_gas_drydep_list.xml @@ -30,6 +30,7 @@ ETHP ETNO3 ETP + FURA GLYC GLYX H2O2 diff --git a/bld/namelist_files/geoschem_master_gas_wetdep_list.xml b/bld/namelist_files/geoschem_master_gas_wetdep_list.xml index 419f518c32..61c4ddd3a2 100644 --- a/bld/namelist_files/geoschem_master_gas_wetdep_list.xml +++ b/bld/namelist_files/geoschem_master_gas_wetdep_list.xml @@ -20,6 +20,7 @@ ETHN ETHP ETP + FURA GLYC GLYX H2O2 diff --git a/bld/namelist_files/namelist_defaults_cam.xml b/bld/namelist_files/namelist_defaults_cam.xml index f55d17a59c..3dfd7dbf5b 100644 --- a/bld/namelist_files/namelist_defaults_cam.xml +++ b/bld/namelist_files/namelist_defaults_cam.xml @@ -2026,7 +2026,7 @@ atm/cam/chem/trop_mozart/dvel/dep_data_c20221208.nc -atm/cam/chem/geoschem/dvel/dep_data_file_geoschem_c230417.nc +atm/cam/chem/geoschem/dvel/dep_data_file_geoschem_v14.3_2023Jan08.nc atm/waccm/phot/effxstex.txt @@ -3700,4 +3700,9 @@ 144 91 + +atm/cam/geoschem/emis/ExtData/CHEM_INPUTS +atm/cam/geoschem/emis/ExtData/CHEM_INPUTS/FAST_JX/v2024-05/ +atm/cam/geoschem/emis/ExtData/CHEM_INPUTS/CLOUD_J/v2023-05/ + diff --git a/bld/namelist_files/namelist_definition.xml b/bld/namelist_files/namelist_definition.xml index fad758cad1..bd003c779a 100644 --- a/bld/namelist_files/namelist_definition.xml +++ b/bld/namelist_files/namelist_definition.xml @@ -5241,12 +5241,23 @@ Default: set by build-namelist for climo cases, otherwise -1 to use model clock. - -Full pathname to GEOS-Chem chemistry inputs directory +Full pathname of the GEOS-Chem chemistry inputs directory. Default: set by build-namelist. + +Full pathname of the aerosols optical property inputs directory used in GEOS-Chem outside of Cloud-J. +Default: set by build-namelist. + + + +Full pathname of the Cloud-J photolysis inputs directory used in GEOS-Chem. +Default: set by build-namelist. + diff --git a/bld/namelist_files/use_cases/2000_geoschem.xml b/bld/namelist_files/use_cases/2000_geoschem.xml index 7463a49361..3d2c5507b5 100644 --- a/bld/namelist_files/use_cases/2000_geoschem.xml +++ b/bld/namelist_files/use_cases/2000_geoschem.xml @@ -6,8 +6,6 @@ -atm/cam/geoschem/ExtData/CHEM_INPUTS/ - atm/cam/geoschem/initial_conditions/f.e20.FC2010.f09_f09.144.GC_vbsext.001.cam.i.0007-01-01-00000.nc atm/cam/geoschem/initial_conditions/f.e20.FC2010.f19_f19.144.GC_vbsext.001.cam.i.0007-01-01-00000.nc diff --git a/bld/namelist_files/use_cases/2010_geoschem.xml b/bld/namelist_files/use_cases/2010_geoschem.xml index 2d3ee5db95..b1b0f9f2eb 100644 --- a/bld/namelist_files/use_cases/2010_geoschem.xml +++ b/bld/namelist_files/use_cases/2010_geoschem.xml @@ -4,8 +4,6 @@ -atm/cam/geoschem/ExtData/CHEM_INPUTS/ - atm/cam/geoschem/initial_conditions/f.e20.FC2010.f09_f09.144.GC_vbsext.001.cam.i.0007-01-01-00000.nc atm/cam/geoschem/initial_conditions/f.e20.FC2010.f19_f19.144.GC_vbsext.001.cam.i.0007-01-01-00000.nc diff --git a/bld/namelist_files/use_cases/hist_geoschem.xml b/bld/namelist_files/use_cases/hist_geoschem.xml index 78b681e572..1cfff4a8a9 100644 --- a/bld/namelist_files/use_cases/hist_geoschem.xml +++ b/bld/namelist_files/use_cases/hist_geoschem.xml @@ -6,8 +6,6 @@ -atm/cam/geoschem/ExtData/CHEM_INPUTS/ - atm/cam/geoschem/initial_conditions/f.e20.FC2010.f09_f09.144.GC_vbsext.001.cam.i.0007-01-01-00000.nc atm/cam/geoschem/initial_conditions/f.e20.FC2010.f19_f19.144.GC_vbsext.001.cam.i.0007-01-01-00000.nc @@ -29,7 +27,7 @@ 0.25D0 SERIAL -atm/waccm/lb/LBC_1750-2014_CMIP6_0p5degLat_c170126.nc +atm/waccm/lb/LBC_17500116-25001216_CMIP6_SSP585_0p5degLat_c20200824.nc SERIAL diff --git a/bld/namelist_files/use_cases/hist_geoschem_nudged.xml b/bld/namelist_files/use_cases/hist_geoschem_nudged.xml index 0550880d80..3c87bcb4fc 100644 --- a/bld/namelist_files/use_cases/hist_geoschem_nudged.xml +++ b/bld/namelist_files/use_cases/hist_geoschem_nudged.xml @@ -6,8 +6,6 @@ -atm/cam/geoschem/ExtData/CHEM_INPUTS/ - atm/cam/geoschem/initial_conditions/f.e20.FC2010.f09_f09.144.GC_vbsext.001.cam.i.0007-01-01-00000.nc atm/cam/geoschem/initial_conditions//f.e20.FC2010.f19_f19.144.GC_vbsext.001.cam.i.0007-01-01-00000.nc @@ -29,7 +27,7 @@ 0.25D0 SERIAL -atm/waccm/lb/LBC_1750-2014_CMIP6_0p5degLat_c170126.nc +atm/waccm/lb/LBC_17500116-25001216_CMIP6_SSP585_0p5degLat_c20200824.nc SERIAL diff --git a/cime_config/config_component.xml b/cime_config/config_component.xml index 4a64730bc1..b4d45d98e7 100644 --- a/cime_config/config_component.xml +++ b/cime_config/config_component.xml @@ -10,8 +10,8 @@ --> CAM cam7 physics: CAM cam6 physics: - CAM cam5 physics: - CAM cam4 physics: + CAM cam5 physics: + CAM cam4 physics: CAM simplified and non-versioned physics : diff --git a/cime_config/testdefs/testlist_cam.xml b/cime_config/testdefs/testlist_cam.xml index 2d736a2d51..b2c62fd5e7 100644 --- a/cime_config/testdefs/testlist_cam.xml +++ b/cime_config/testdefs/testlist_cam.xml @@ -17,14 +17,6 @@ - - - - - - - - @@ -42,14 +34,6 @@ - - - - - - - - @@ -1890,6 +1874,14 @@ + + + + + + + + @@ -1900,6 +1892,11 @@ + + + + + diff --git a/doc/ChangeLog b/doc/ChangeLog index 4cc0bfd27c..c576afc24e 100644 --- a/doc/ChangeLog +++ b/doc/ChangeLog @@ -1,3 +1,137 @@ +Tag name: cam6_4_038 +Originator(s): lizziel +Date: 1 Oct 2024 +One-line Summary: Update to GEOS-Chem 14.4.3 and HEMCO 3.9.0, and add Cloud-J 7.7.3 and HETP 1.0 +Github PR URL: https://github.com/ESCOMP/CAM/pull/1093 + +Purpose of changes (include the issue number and title text for each relevant GitHub issue): +Update GEOS-Chem and dependencies to newer versions. (issue #959) + +Describe any changes made to build system: + +Describe any changes made to the namelist: + +List any changes to the defaults for the boundary datasets: + +Describe any substantial timing or memory changes: + +Code reviewed by: fvitt cacraigucar + +List all files eliminated: N/A + +List all files added and what they do: N/A + +List all existing files that have been modified, and describe the changes: + +M .gitmodules + - Update GEOS-Chem and HEMCO_CESM versions + - Add new dependencies cloud_j and hetp + +M bld/build-namelist + - Add defaults for GEOS-Chem namelist input directories + +M bld/configure + - Add MODEL_CESM and MODEL_GEOSCHEM to geoschem chem_cppdefs for use in Cloud-J + - Update hard-coded chem_nadv used for GEOS-Chem + +M bld/namelist_files/geoschem_master_gas_drydep_list.xml +M bld/namelist_files/geoschem_master_gas_wetdep_list.xml + - Add new GEOS-Chem deposited species FURA + - Add SO4 which is a non-MAM aerosol species in GEOS-Chem + +M bld/namelist_files/namelist_defaults_cam.xml + - Update default dep_data_file used for GEOS-Chem + - Add two input directories for GEOS-Chem + +M bld/namelist_files/namelist_definition.xml + - Change GEOS-Chem chem inputs directory name + - Add new GEOS-Chem input directory used for Cloud-J + +M bld/namelist_files/use_cases/2000_geoschem.xml +M bld/namelist_files/use_cases/2010_geoschem.xml +M bld/namelist_files/use_cases/hist_geoschem.xml +M bld/namelist_files/use_cases/hist_geoschem_nudged.xml + - Delete GEOS-Chem chem inputs directory since now in defaults + +M cime_config/buildnml + - No-diff cleanup + +M cime_config/config_component.xml + - Remove GEOS-chem from CAM40 and CAM50 options + +M cime_config/config_compsets.xml + - Fix bug in FCnudged_GC where nudging incorrectly specified in compset definition + +M testdefs/testlist_cam.xml + - Move GEOS-Chem tests to be grouped with similar compsets + +A src/chemistry/cloud_j + - tagged version 7.7.3 + +M src/chemistry/geoschem/chem_mods.F90 + - Update hard-coded nTracersMax, gas_pcnst, and nslvd used for + GEOS-Chem compsets + +M src/chemistry/geoschem/chemistry.F90 + - Updates for compatibility with GEOS-Chem 14.3.1 + - Replace hard-coded photolysis inputs directory with namelist value + - Add aerosol distribution fix to deposit SO4 using Neu as gas instead of aerosol + - Update precision modifiers + +M src/chemistry/geoschem/geoschem_diagnostics_mod.F90 + - Update precision modifiers + +M src/chemistry/geoschem/geoschem_emissions_mod.F90 + - Update precision modifiers + +M src/chemistry/geoschem/geoschem_src + - tagged version 14.4.3 + +M src/chemistry/geoschem/mo_sim_dat.F90 + - Add new species to solsym: BUTDI (butendedial) and FURA (furan) + +A src/chemistry/hetp + - tagged version 1.0 + +M src/chemistry/mozart/mo_neu_wetdep.F90 + - Restrict henrys law coefficient assignment error print to masterproc + - Add GEOS-Chem bulk sulfates SO4 and SO4S + +M src/hemco + - tag hemco-cesm2_0_hemco3_9_0 + +M test/system/TR8.sh + - Skip over GEOS-Chem and dependency modules + +If there were any failures reported from running test_driver.sh on any test +platform, and checkin with these failures has been OK'd by the gatekeeper, +then copy the lines from the td.*.status files for the failed tests to the +appropriate machine below. All failed tests must be justified. + +derecho/intel/aux_cam: + + PEND SMS_D_Ln9_P1280x1.ne0CONUSne30x8_ne0CONUSne30x8_mt12.FCHIST.derecho_intel.cam-outfrq9s + - pre-existing failures -- need fix in CLM external + + FAIL ERP_Ln9.f09_f09_mg17.FCSD_HCO.derecho_intel.cam-outfrq9s COMPARE_base_rest + DIFF SMS_Lh12.f09_f09_mg17.FCSD_HCO.derecho_intel.cam-outfrq3h + - pre-existing failure due to HEMCO not having reproducible results issues #1018 and #856 + + DIFF SMS_Ld1.f09_f09_mg17.FCHIST_GC.derecho_intel.cam-outfrq1d + - expected baseline failure due to updates to GEOS-Chem + +derecho/nvhpc/aux_cam: All PASS + +izumi/nag/aux_cam: + FAIL DAE.f45_f45_mg37.FHS94.izumi_nag.cam-dae + - pre-existing failure - issue #670 + +izumi/gnu/aux_cam: ALL PASS + +Summarize any changes to answers: bit-for-bit unchanged + +=============================================================== + Tag name: cam6_4_037 Originator(s): jimmielin Date: Sep 30, 2024 @@ -144,7 +278,6 @@ izumi/gnu/aux_cam: ALL PASS Summarize any changes to answers: bit-for-bit -=============================================================== =============================================================== Tag name: cam6_4_035 @@ -4964,7 +5097,6 @@ izumi/nag/aux_cam: all B4B, except: izumi/gnu/aux_cam: all BFB Summarize any changes to answers: bit-for-bit unchanged except GEOS-Chem and HEMCO tests described in issue #1018 - =============================================================== Tag name: cam6_3_158 diff --git a/src/chemistry/cloud_j b/src/chemistry/cloud_j new file mode 160000 index 0000000000..d20050f1ef --- /dev/null +++ b/src/chemistry/cloud_j @@ -0,0 +1 @@ +Subproject commit d20050f1ef9e3895d58f3041efd2da59ce1ed421 diff --git a/src/chemistry/geoschem/chem_mods.F90 b/src/chemistry/geoschem/chem_mods.F90 index 2d8a500253..7ce98d2ab5 100644 --- a/src/chemistry/geoschem/chem_mods.F90 +++ b/src/chemistry/geoschem/chem_mods.F90 @@ -7,7 +7,7 @@ module chem_mods implicit none save - INTEGER, PARAMETER :: nTracersMax = 267 ! Must be equal to chem_nadv + INTEGER, PARAMETER :: nTracersMax = 269 ! = chem_nadv in cam/bld/configure INTEGER :: nTracers REAL(r8) :: ref_MMR(pcnst) @@ -62,10 +62,10 @@ module chem_mods rxntot = 212, & ! number of total reactions gascnt = 172, & ! number of gas phase reactions nabscol = 2, & ! number of absorbing column densities - gas_pcnst = 269, & ! number of "gas phase" species (same as solsym length) - ! Includes GC advected species (233), MAM aerosols (33), - ! and CO2 (1), as well as any non-advected species added - ! to solsym and mo_sim_dat.F90. + gas_pcnst = 357, & ! number of "gas phase" species (=solsym length) + ! Includes GC advected species, MAM aerosols, + ! and CO2, and any non-advected species added + ! to solsym within mo_sim_dat.F90. nfs = 6, & ! number of "fixed" species relcnt = 0, & ! number of relationship species grpcnt = 0, & ! number of group members @@ -81,7 +81,7 @@ module chem_mods clsze = 1, & ! loop length for implicit chemistry rxt_tag_cnt = 0, & ! number of tagged reactions (unused in GEOS-Chem) enthalpy_cnt = 0, & - nslvd = 86 ! number of short-lived (non-advected) species + nslvd = 88 ! number of short-lived (non-advected) species minus CO2 integer :: clscnt(5) = 0 integer :: cls_rxt_cnt(4,5) = 0 integer :: clsmap(gas_pcnst,5) = 0 diff --git a/src/chemistry/geoschem/chemistry.F90 b/src/chemistry/geoschem/chemistry.F90 index 945346f263..fbb99e4b8f 100644 --- a/src/chemistry/geoschem/chemistry.F90 +++ b/src/chemistry/geoschem/chemistry.F90 @@ -60,17 +60,48 @@ module chemistry public :: chem_emissions public :: chem_timestep_init + ! + ! Private routines: + ! + private :: sect02_mam4 + private :: erfc_num_recipes + ! Location of valid geoschem_config.yml and species_database.yml ! Use local files in run folder CHARACTER(LEN=500) :: gcConfig = 'geoschem_config.yml' CHARACTER(LEN=500) :: speciesDB = 'species_database.yml' - ! Location of chemistry input - CHARACTER(LEN=shr_kind_cl) :: geoschem_cheminputs + CHARACTER(LEN=shr_kind_cl) :: geoschem_chem_inputs + CHARACTER(LEN=shr_kind_cl) :: geoschem_aeropt_inputs + CHARACTER(LEN=shr_kind_cl) :: geoschem_photol_inputs ! Debugging LOGICAL :: debug = .TRUE. + ! Compile-time logical controls. These options are usually not expected + ! to change from run to run and significantly affect the behavior + ! of the model and thus are not read through the namelist. + + ! Use SOA initial conditions from MAM4 (soaX_aY) or GEOS-Chem (SOA*) + ! in the initial conditions (ncdata) / restart file? + LOGICAL :: useSOAICfromMAM4 = .TRUE. + + ! Map back SOAs from MAM4 at the beginning of every chemistry timestep? + ! There are several implications: + ! - MAM4 will perform deposition of SOAs, changing the bulk mass; + ! if disabled, only one-way mapping of GC aerosols to MAM4 is done. + ! deposition of SOAs will still be performed but based on GEOS-Chem species. + ! Either approach is scientifically valid. + LOGICAL :: useMAM4mapBackSOA = .FALSE. + + ! Prescribe aerosol size distributions based on Feng et al. (2021) GMD? + ! This is intended to stabilize the model if only for gas-phase chemistry + ! purposes and will provide a more reasonable radiative/cloud properties. + ! However, it will complicate climate/geoengineering simulations as MAM4 + ! will lose control of the sulfate size distribution. + LOGICAL :: usePrescribedAerDistribution = .FALSE. + + ! Derived type objects TYPE(OptInput) :: Input_Opt ! Input Options object TYPE(ChmState),ALLOCATABLE :: State_Chm(:) ! Chemistry State object @@ -252,7 +283,7 @@ subroutine chem_register ! Options needed by Init_State_Chm IO%ITS_A_FULLCHEM_SIM = .True. IO%LLinoz = .True. - IO%LPRT = .False. + IO%Verbose = .False. IO%N_Advect = nTracers DO I = 1, nTracers IO%AdvectSpc_Name(I) = TRIM(tracerNames(I)) @@ -971,7 +1002,6 @@ subroutine chem_init(phys_state, pbuf2d) use geoschem_history_mod, only : HistoryExports_SetServices ! GEOS-Chem modules - use Chemistry_Mod, only : Init_Chemistry use DiagList_Mod, only : Init_DiagList, Print_DiagList use Drydep_Mod, only : depName, Ndvzind use Error_Mod, only : Init_Error @@ -980,9 +1010,10 @@ subroutine chem_init(phys_state, pbuf2d) use GC_Grid_Mod, only : SetGridFromCtrEdges use Input_Mod, only : Read_Input_File, Validate_Directories use Input_Opt_Mod, only : Set_Input_Opt - use isorropiaII_Mod, only : Init_IsorropiaII + use Aerosol_Thermodynamics_Mod, only : Init_ATE use Linear_Chem_Mod, only : Init_Linear_Chem use Linoz_Mod, only : Linoz_Read + use Photolysis_Mod, only : Init_Photolysis use PhysConstants, only : PI, PI_180, Re use Pressure_Mod, only : Accept_External_ApBp use State_Chm_Mod, only : Ind_ @@ -990,6 +1021,7 @@ subroutine chem_init(phys_state, pbuf2d) use TaggedDiagList_Mod, only : Init_TaggedDiagList, Print_TaggedDiagList use Time_Mod, only : Accept_External_Date_Time use Ucx_Mod, only : Init_Ucx + use Unitconv_Mod, only : MOLES_SPECIES_PER_MOLES_DRY_AIR use Vdiff_Mod, only : Max_PblHt_For_Vdiff TYPE(physics_state), INTENT(IN ) :: phys_state(BEGCHUNK:ENDCHUNK) @@ -1096,12 +1128,12 @@ subroutine chem_init(phys_state, pbuf2d) ! Find maximum tropopause level, set at 40 hPa (based on GEOS-Chem 72 and 47 ! layer grids) nTrop = nZ - DO WHILE ( hyam(nZ+1-nTrop) * ps0 < 4000.0 ) + DO WHILE ( hyam(nZ+1-nTrop) * ps0 < 4000.0_r8 ) nTrop = nTrop-1 ENDDO ! Find stratopause level, defined at 1 hPa nStrat = nZ - DO WHILE ( hyam(nZ+1-nStrat) * ps0 < 100.0 ) + DO WHILE ( hyam(nZ+1-nStrat) * ps0 < 100.0_r8 ) nStrat = nStrat-1 ENDDO @@ -1132,10 +1164,12 @@ subroutine chem_init(phys_state, pbuf2d) State_Grid = maxGrid, & RC = RC ) - ! First setup directories - Input_Opt%Chem_Inputs_Dir = TRIM(geoschem_cheminputs) + ! First setup directories. FAST-JX directory is still used for + ! optical properties of aerosols outside of Cloud-J. + Input_Opt%Chem_Inputs_Dir = TRIM(geoschem_chem_inputs) Input_Opt%SpcDatabaseFile = TRIM(speciesDB) - Input_Opt%FAST_JX_DIR = TRIM(geoschem_cheminputs)//'FAST_JX/v2020-02/' + Input_Opt%FAST_JX_DIR = TRIM(geoschem_aeropt_inputs) + Input_Opt%CLOUDJ_DIR = TRIM(geoschem_photol_inputs) !---------------------------------------------------------- ! CESM-specific input flags @@ -1158,7 +1192,7 @@ subroutine chem_init(phys_state, pbuf2d) CALL Validate_Directories( Input_Opt, RC ) IF ( RC /= GC_SUCCESS ) THEN - ErrMsg = 'Error encountered in "Validation_Directories"!' + ErrMsg = 'Error encountered in "Validate_Directories"!' CALL Error_Stop( ErrMsg, ThisLoc ) ENDIF @@ -1168,7 +1202,7 @@ subroutine chem_init(phys_state, pbuf2d) RC = RC ) IF ( RC /= GC_SUCCESS ) THEN - ErrMsg = 'Error encountered within call to "GC_Init_Grid"!' + ErrMsg = 'Error encountered within call to "GC_Init_Grid" (1 - maxGrid)!' CALL Error_Stop( ErrMsg, ThisLoc ) ENDIF @@ -1189,6 +1223,10 @@ subroutine chem_init(phys_state, pbuf2d) CALL Error_Stop( ErrMsg, ThisLoc ) ENDIF + ! Set grid metadata. This has to be after State_Grid is initialized. + State_Grid(I)%CPU_Subdomain_ID = I + State_Grid(I)%CPU_Subdomain_FirstID = BEGCHUNK + State_Grid(I)%NX = nX State_Grid(I)%NY = NCOL(I) State_Grid(I)%NZ = nZ @@ -1199,7 +1237,7 @@ subroutine chem_init(phys_state, pbuf2d) RC = RC ) IF ( RC /= GC_SUCCESS ) THEN - ErrMsg = 'Error encountered within call to "GC_Init_Grid"!' + ErrMsg = 'Error encountered within call to "GC_Init_Grid" (2 - chunk)!' CALL Error_Stop( ErrMsg, ThisLoc ) ENDIF @@ -1356,7 +1394,7 @@ subroutine chem_init(phys_state, pbuf2d) ENDIF ! Set a flag to denote if we should print ND70 debug output - prtDebug = ( Input_Opt%LPRT .and. MasterProc ) + prtDebug = ( Input_Opt%Verbose .and. MasterProc ) historyConfigFile = 'HISTORY.rc' ! This requires geoschem_config.yml and HISTORY.rc to be in the run directory @@ -1400,6 +1438,7 @@ subroutine chem_init(phys_state, pbuf2d) ENDIF DO I = BEGCHUNK, ENDCHUNK + ! Restrict prints to one thread only Input_Opt%amIRoot = (MasterProc .AND. (I == BEGCHUNK)) CALL GC_Init_StateObj( Diag_List = Diag_List, & ! Diagnostic list obj @@ -1418,7 +1457,9 @@ subroutine chem_init(phys_state, pbuf2d) ENDIF ! Start with v/v dry (CAM standard) - State_Chm(I)%Spc_Units = 'v/v dry' + DO N = 1, State_Chm(I)%nSpecies + State_Chm(I)%Species(N)%Units = MOLES_SPECIES_PER_MOLES_DRY_AIR + ENDDO ENDDO Input_Opt%amIRoot = MasterProc @@ -1572,17 +1613,27 @@ subroutine chem_init(phys_state, pbuf2d) CALL Error_Stop( ErrMsg, ThisLoc ) ENDIF - IF ( Input_Opt%Its_A_FullChem_Sim .OR. & - Input_Opt%Its_An_Aerosol_Sim ) THEN - ! This also initializes Fast-JX - CALL Init_Chemistry( Input_Opt = Input_Opt, & - State_Chm = State_Chm(BEGCHUNK), & - State_Diag = State_Diag(BEGCHUNK), & - State_Grid = State_Grid(BEGCHUNK), & - RC = RC ) + ! Initialize photolysis, including reading files for optical properties. + IF ( Input_Opt%ITS_A_FULLCHEM_SIM .or. & + Input_Opt%ITS_AN_AEROSOL_SIM ) THEN + DO I = BEGCHUNK, ENDCHUNK + CALL Init_Photolysis( Input_Opt = Input_Opt, & + State_Grid = State_Grid(I), & + State_Chm = State_Chm(I), & + State_Diag = State_Diag(I), & + RC = RC ) + + ! Only the root chunk (on all CPUs) should be reading the data + ! in State_Chm%Phot%OREF and State_Chm%Phot%TREF, and the rest should be copied. + ! This fixes a hang condition in the ne30 (SE dycore) compsets. (hplin, 7/3/24) + IF( I .ne. BEGCHUNK ) THEN + State_Chm(I)%Phot%TREF = State_Chm(BEGCHUNK)%Phot%TREF + State_Chm(I)%Phot%OREF = State_Chm(BEGCHUNK)%Phot%OREF + ENDIF + ENDDO IF ( RC /= GC_SUCCESS ) THEN - ErrMsg = 'Error encountered in "Init_Chemistry"!' + ErrMsg = 'Error encountered in "Init_Photolysis"!' CALL Error_Stop( ErrMsg, ThisLoc ) ENDIF ENDIF @@ -1607,7 +1658,7 @@ subroutine chem_init(phys_state, pbuf2d) CALL mpi_bcast( State_Chm(I)%NOXCOEFF, size(State_Chm(I)%NOXCOEFF), mpi_real8, masterprocid, mpicom, ierr ) IF ( ierr /= mpi_success ) CALL endrun('Error in mpi_bcast of NOXCOEFF in first chunk') ELSE - State_CHM(I)%NOXCOEFF = State_Chm(BEGCHUNK)%NOXCOEFF + State_Chm(I)%NOXCOEFF = State_Chm(BEGCHUNK)%NOXCOEFF ENDIF ENDDO ENDIF @@ -1626,7 +1677,7 @@ subroutine chem_init(phys_state, pbuf2d) ENDIF IF ( Input_Opt%LSSalt ) THEN - CALL INIT_ISORROPIAII( State_Grid = maxGrid ) + CALL INIT_ATE( State_Grid = maxGrid ) ENDIF ! Get some indices @@ -1770,7 +1821,9 @@ subroutine geoschem_readnl(nlfile) integer :: unitn, ierr character(len=*), parameter :: subname = 'geoschem_readnl' - namelist /geoschem_nl/ geoschem_cheminputs + namelist /geoschem_nl/ geoschem_chem_inputs + namelist /geoschem_nl/ geoschem_aeropt_inputs + namelist /geoschem_nl/ geoschem_photol_inputs ! Read namelist IF ( MasterProc ) THEN @@ -1788,9 +1841,19 @@ subroutine geoschem_readnl(nlfile) ENDIF ! Broadcast namelist variables - CALL mpi_bcast(geoschem_cheminputs, LEN(geoschem_cheminputs), mpi_character, masterprocid, mpicom, ierr) + CALL mpi_bcast(geoschem_chem_inputs, LEN(geoschem_chem_inputs), mpi_character, masterprocid, mpicom, ierr) + IF ( ierr /= mpi_success ) then + CALL endrun(subname//': MPI_BCAST ERROR: geoschem_chem_inputs') + ENDIF + + CALL mpi_bcast(geoschem_aeropt_inputs, LEN(geoschem_aeropt_inputs), mpi_character, masterprocid, mpicom, ierr) IF ( ierr /= mpi_success ) then - CALL endrun(subname//': MPI_BCAST ERROR: geoschem_cheminputs') + CALL endrun(subname//': MPI_BCAST ERROR: geoschem_aeropt_inputs') + ENDIF + + CALL mpi_bcast(geoschem_photol_inputs, LEN(geoschem_photol_inputs), mpi_character, masterprocid, mpicom, ierr) + IF ( ierr /= mpi_success ) then + CALL endrun(subname//': MPI_BCAST ERROR: geoschem_photol_inputs') ENDIF end subroutine geoschem_readnl @@ -1840,14 +1903,12 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) use GeosChem_History_Mod, only : HistoryExports_SetDataPointers, CopyGCStates2Exports ! GEOS-Chem modules - use Aerosol_Mod, only : Set_AerMass_Diagnostic use Calc_Met_Mod, only : Set_Dry_Surface_Pressure, AirQnt use Chemistry_Mod, only : Do_Chemistry - use CMN_FJX_MOD, only : ZPJ use CMN_Size_Mod, only : NSURFTYPE, PTop use Diagnostics_Mod, only : Zero_Diagnostics_StartOfTimestep, Set_Diagnostics_EndofTimestep + use Diagnostics_Mod, only : Set_AerMass_Diagnostic use Drydep_Mod, only : Do_Drydep, DEPNAME, NDVZIND, Update_DryDepFreq - use FAST_JX_MOD, only : RXN_NO2, RXN_O3_1 use GC_Grid_Mod, only : SetGridFromCtr use HCO_Interface_GC_Mod,only : Compute_Sflx_For_Vdiff use Linear_Chem_Mod, only : TrID_GC, GC_Bry_TrID, NSCHEM @@ -1862,7 +1923,8 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) use Time_Mod, only : Accept_External_Date_Time use Toms_Mod, only : Compute_Overhead_O3 use UCX_Mod, only : Set_H2O_Trac - use Unitconv_Mod, only : Convert_Spc_Units + use Unitconv_Mod, only : Convert_Spc_Units, UNIT_STR + use Unitconv_Mod, only : KG_SPECIES_PER_KG_DRY_AIR, KG_SPECIES_PER_M2 use Wetscav_Mod, only : Setup_Wetscav REAL(r8), INTENT(IN) :: dT ! Time step @@ -1939,6 +2001,14 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) INTEGER :: iMap, nMapping, iBin, binSOA_1, binSOA_2 INTEGER :: K1, K2, K3, K4 LOGICAL :: isSOA_aerosol + CHARACTER(LEN=64) :: aerName + + ! For prescribed aerosol distributions. + REAL(r8) :: prescr_aer_xnum(3) + REAL(r8) :: prescr_aer_xmas(3) + REAL(r8) :: vmr_so4_sum(state%ncol, pver) + REAL(r8) :: prescr_aer_lbnd, prescr_aer_abnd, prescr_aer_cbnd, prescr_aer_ubnd + REAL(r8), POINTER :: dgncur_a(:,:,:) #endif @@ -1990,7 +2060,7 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) TYPE(Species), POINTER :: SpcInfo TYPE(SfcMrObj), POINTER :: iSfcMrObj - CHARACTER(LEN=63) :: OrigUnit + INTEGER :: previous_units REAL(r8) :: SlsData(PCOLS, PVER, nSls) @@ -2012,6 +2082,7 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) LOGICAL :: lastChunk INTEGER :: RC + call t_startf('GEOSChem_All_Tend') ! Initialize pointers SpcInfo => NULL() @@ -2028,6 +2099,7 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) pbuf_ik => NULL() pbuf_i => NULL() + call t_startf('GEOSChem_CAM_Interfacing') ! LCHNK: which chunk we have on this process LCHNK = state%LCHNK ! NCOL: number of atmospheric columns on this chunk @@ -2098,9 +2170,9 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) ! 2. Copy tracers into State_Chm ! Data was received in kg/kg dry - State_Chm(LCHNK)%Spc_Units = 'kg/kg dry' - ! Initialize ALL State_Chm species data to zero, not just tracers + ! Initialize ALL State_Chm species data to zero, not just tracer DO N = 1, State_Chm(LCHNK)%nSpecies + State_Chm(LCHNK)%Species(N)%Units = KG_SPECIES_PER_KG_DRY_AIR State_Chm(LCHNK)%Species(N)%Conc = 0.0e+0_fp ENDDO @@ -2196,6 +2268,11 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) State_Chm(LCHNK)%Species(M)%Conc(1,:nY,:nZ) = REAL(SlsData(:nY,nZ:1:-1,N),fp) ENDDO + call t_stopf('GEOSChem_CAM_Interfacing') + + ! We want to put t_startf timers outside of C-preprocessor flags + ! in order to always have these timers present even if zero. (hplin, 4/30/24) + call t_startf('GEOSChem_MAM_Interfacing') #if defined( MODAL_AERO ) ! NOTE: GEOS-Chem bulk aerosol concentrations (BCPI, BCPO, SO4, ...) are ZEROED OUT ! here in order to be reconstructed from the modal concentrations. @@ -2208,6 +2285,12 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) DO SM = 1, nspec_amode(M) P = map2MAM4(SM,M) ! Constituent index for GEOS-Chem IF ( P > 0 ) K = map2GC(P) ! Index in State_Chm + + ! do not zero out sulfate aerosol here since aerosol distribution for sulfate + ! will be prescribed (hplin, 5/9/23) + call rad_cnst_get_info(0,M,SM,spec_name=aerName) + IF ( to_upper(aerName(:3)) == "SO4" ) CYCLE + IF ( K > 0 ) State_Chm(LCHNK)%Species(K)%Conc(1,:nY,:nZ) = 0.0e+00_fp ENDDO ENDDO @@ -2223,6 +2306,11 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) ! species (with cnst index P, which corresponds to index K in ! State_Chm) + ! do not zero out sulfate aerosol here since aerosol distribution for sulfate + ! will be prescribed (hplin, 5/9/23) + call rad_cnst_get_info(0,M,SM,spec_name=aerName) + IF ( to_upper(aerName(:3)) == "SO4" ) CYCLE + ! Multiple MAM4 bins are mapped to same GEOS-Chem species State_Chm(LCHNK)%Species(K)%Conc(1,:nY,:nZ) = State_Chm(LCHNK)%Species(K)%Conc(1,:nY,:nZ) & + REAL(state%q(:nY,nZ:1:-1,N),fp) * & @@ -2292,76 +2380,81 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) ! TSOG2 + ASOG2 <- SOAG3 ! TSOG3 + ASOG3 <- SOAG4 - IF ( iStep > 1 ) THEN - ! Do not perform this mapping on initialization as we first want to - ! overwrite soa*_a* with the GEOS-Chem SOAs. - nMapping = 8 - DO iMap = 1, nMapping - speciesName_1 = '' - speciesName_2 = '' - speciesName_3 = '' - speciesName_4 = '' - IF ( iMap == 1 ) THEN - binSOA_1 = 1 - binSOA_2 = 2 - speciesName_1 = 'TSOA0' - speciesName_2 = 'ASOAN' - speciesName_3 = 'SOAIE' - speciesName_4 = 'SOAGX' - ELSEIF ( iMap == 2 ) THEN - binSOA_1 = 3 - binSOA_2 = 3 - speciesName_1 = 'TSOA1' - speciesName_2 = 'ASOA1' - ELSEIF ( iMap == 3 ) THEN - binSOA_1 = 4 - binSOA_2 = 4 - speciesName_1 = 'TSOA2' - speciesName_2 = 'ASOA2' - ELSEIF ( iMap == 4 ) THEN - binSOA_1 = 5 - binSOA_2 = 5 - speciesName_1 = 'TSOA3' - speciesName_2 = 'ASOA3' - ELSEIF ( iMap == 5 ) THEN - binSOA_1 = 1 - binSOA_2 = 2 - speciesName_1 = 'TSOG0' - speciesName_2 = 'TSOG0' - ELSEIF ( iMap == 6 ) THEN - binSOA_1 = 3 - binSOA_2 = 3 - speciesName_1 = 'TSOG1' - speciesName_2 = 'ASOG1' - ELSEIF ( iMap == 7 ) THEN - binSOA_1 = 4 - binSOA_2 = 4 - speciesName_1 = 'TSOG2' - speciesName_2 = 'ASOG2' - ELSEIF ( iMap == 8 ) THEN - binSOA_1 = 5 - binSOA_2 = 5 - speciesName_1 = 'TSOG3' - speciesName_2 = 'ASOG3' - ELSE - CALL ENDRUN('Unknown SOA mapping!') - ENDIF - isSOA_aerosol = .False. - IF ( iMap <= 4 ) isSOA_aerosol = .True. - - ! Compute total mass from GEOS-Chem species. This sets the ratio between - ! speciesId_1 and speciesId_2 - totMass(:nY,:nZ) = 0.0e+00_r8 - - CALL cnst_get_ind( speciesName_1, speciesId_1, abort=.True. ) - CALL cnst_get_ind( speciesName_2, speciesId_2, abort=.False. ) - CALL cnst_get_ind( speciesName_3, speciesId_3, abort=.False. ) - CALL cnst_get_ind( speciesName_4, speciesId_4, abort=.False. ) - IF ( speciesId_1 > 0 ) totMass(:nY,:nZ) = totMass(:nY,:nZ) + state%q(:nY,:nZ,speciesId_1) - IF ( speciesId_2 > 0 ) totMass(:nY,:nZ) = totMass(:nY,:nZ) + state%q(:nY,:nZ,speciesId_2) - IF ( speciesId_3 > 0 ) totMass(:nY,:nZ) = totMass(:nY,:nZ) + state%q(:nY,:nZ,speciesId_3) - IF ( speciesId_4 > 0 ) totMass(:nY,:nZ) = totMass(:nY,:nZ) + state%q(:nY,:nZ,speciesId_4) - + nMapping = 8 + DO iMap = 1, nMapping + speciesName_1 = '' + speciesName_2 = '' + speciesName_3 = '' + speciesName_4 = '' + IF ( iMap == 1 ) THEN + binSOA_1 = 1 + binSOA_2 = 2 + speciesName_1 = 'TSOA0' + speciesName_2 = 'ASOAN' + speciesName_3 = 'SOAIE' + speciesName_4 = 'SOAGX' + ELSEIF ( iMap == 2 ) THEN + binSOA_1 = 3 + binSOA_2 = 3 + speciesName_1 = 'TSOA1' + speciesName_2 = 'ASOA1' + ELSEIF ( iMap == 3 ) THEN + binSOA_1 = 4 + binSOA_2 = 4 + speciesName_1 = 'TSOA2' + speciesName_2 = 'ASOA2' + ELSEIF ( iMap == 4 ) THEN + binSOA_1 = 5 + binSOA_2 = 5 + speciesName_1 = 'TSOA3' + speciesName_2 = 'ASOA3' + ELSEIF ( iMap == 5 ) THEN + binSOA_1 = 1 + binSOA_2 = 2 + speciesName_1 = 'TSOG0' + speciesName_2 = 'TSOG0' + ELSEIF ( iMap == 6 ) THEN + binSOA_1 = 3 + binSOA_2 = 3 + speciesName_1 = 'TSOG1' + speciesName_2 = 'ASOG1' + ELSEIF ( iMap == 7 ) THEN + binSOA_1 = 4 + binSOA_2 = 4 + speciesName_1 = 'TSOG2' + speciesName_2 = 'ASOG2' + ELSEIF ( iMap == 8 ) THEN + binSOA_1 = 5 + binSOA_2 = 5 + speciesName_1 = 'TSOG3' + speciesName_2 = 'ASOG3' + ELSE + CALL ENDRUN('Unknown SOA mapping!') + ENDIF + isSOA_aerosol = .False. + IF ( iMap <= 4 ) isSOA_aerosol = .True. + + ! Compute total mass from GEOS-Chem species. This sets the ratio between + ! speciesId_1 and speciesId_2 + totMass(:nY,:nZ) = 0.0e+00_r8 + + CALL cnst_get_ind( speciesName_1, speciesId_1, abort=.True. ) + CALL cnst_get_ind( speciesName_2, speciesId_2, abort=.False. ) + CALL cnst_get_ind( speciesName_3, speciesId_3, abort=.False. ) + CALL cnst_get_ind( speciesName_4, speciesId_4, abort=.False. ) + IF ( speciesId_1 > 0 ) totMass(:nY,:nZ) = totMass(:nY,:nZ) + state%q(:nY,:nZ,speciesId_1) + IF ( speciesId_2 > 0 ) totMass(:nY,:nZ) = totMass(:nY,:nZ) + state%q(:nY,:nZ,speciesId_2) + IF ( speciesId_3 > 0 ) totMass(:nY,:nZ) = totMass(:nY,:nZ) + state%q(:nY,:nZ,speciesId_3) + IF ( speciesId_4 > 0 ) totMass(:nY,:nZ) = totMass(:nY,:nZ) + state%q(:nY,:nZ,speciesId_4) + + K1 = Ind_(speciesName_1) + K2 = Ind_(speciesName_2) + K3 = Ind_(speciesName_3) + K4 = Ind_(speciesName_4) + + ! Check whether to overwrite GEOS-Chem SOAs using concentrations from MAM4. + IF ( (useMAM4mapBackSOA .or. (iStep == 1 .and. useSOAICfromMAM4)) .and. & ! If MAM4 should map back SOAs, then + (useSOAICfromMAM4 .or. iStep > 1) ) THEN ! If use IC, run at all times; otherwise, only run after 1st step to overwrite GC SOAs from soaX_aY ! Compute total bulk mass from MAM bulkMass(:nY,:nZ) = 0.0e+00_r8 IF ( isSOA_aerosol ) THEN @@ -2382,10 +2475,6 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) ENDDO ENDIF - K1 = Ind_(speciesName_1) - K2 = Ind_(speciesName_2) - K3 = Ind_(speciesName_3) - K4 = Ind_(speciesName_4) DO J = 1, nY DO L = 1, nZ ! Total SOA aerosol masses from GC are available. Partition according to the ratio given in speciesId_N to totMass summed above. @@ -2406,12 +2495,19 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) ENDIF ENDDO ENDDO - IF ( K1 > 0 ) MMR_Beg(:nY,:nZ,K1) = State_Chm(LCHNK)%Species(K1)%Conc(1,:nY,:nZ) - IF ( K2 > 0 ) MMR_Beg(:nY,:nZ,K2) = State_Chm(LCHNK)%Species(K2)%Conc(1,:nY,:nZ) - IF ( K3 > 0 ) MMR_Beg(:nY,:nZ,K3) = State_Chm(LCHNK)%Species(K3)%Conc(1,:nY,:nZ) - IF ( K4 > 0 ) MMR_Beg(:nY,:nZ,K4) = State_Chm(LCHNK)%Species(K4)%Conc(1,:nY,:nZ) - ENDDO - ENDIF + ENDIF + + ! Regardless of whether MAM4 will overwrite GEOS-Chem SOA, this part must run, as MMR_Beg is used + ! for computing the flux. If this step is skipped in the first time step, then MMR_Beg is taken + ! as zero and this will result in the entire mass to be provided to the GEOS-Chem species as flux, + ! doubling the species MMRs. + ! + ! Thus, the short-circuiting of the MAM4 to GEOS-Chem mapping must only be done above. (hplin, 5/11/23) + IF ( K1 > 0 ) MMR_Beg(:nY,:nZ,K1) = State_Chm(LCHNK)%Species(K1)%Conc(1,:nY,:nZ) + IF ( K2 > 0 ) MMR_Beg(:nY,:nZ,K2) = State_Chm(LCHNK)%Species(K2)%Conc(1,:nY,:nZ) + IF ( K3 > 0 ) MMR_Beg(:nY,:nZ,K3) = State_Chm(LCHNK)%Species(K3)%Conc(1,:nY,:nZ) + IF ( K4 > 0 ) MMR_Beg(:nY,:nZ,K4) = State_Chm(LCHNK)%Species(K4)%Conc(1,:nY,:nZ) + ENDDO ! Add gas-phase H2SO4 to GEOS-Chem SO4 (which lumps SO4 aerosol and gaseous) K = iSO4 @@ -2446,7 +2542,9 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) mmr_tend(:nY,:nZ,N) = state%q(:nY,:nZ,-M) ENDIF ENDDO + call t_stopf('GEOSChem_MAM_Interfacing') + call t_startf('GEOSChem_CAM_Interfacing') ! If H2O tendencies are propagated to specific humidity, then make sure ! that Q actually applies tendencies IF ( Input_Opt%applyQtend ) lq(cQ) = .True. @@ -2664,7 +2762,7 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) ENDIF ENDDO - ! Field : FRCLND, FRLAND, FROCEAN, FRSEAICE, FRLAKE, FRLANDIC + ! Field : FRCLND, FRLAND, FROCEAN, FRSEAICE, FRLAKE, FRLANDICE ! Description: Olson land fraction ! Fraction of land ! Fraction of ocean @@ -2680,8 +2778,8 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) State_Met(LCHNK)%FROCEAN (1,:nY) = cam_in%ocnFrac(:nY) + cam_in%iceFrac(:nY) State_Met(LCHNK)%FRSEAICE (1,:nY) = cam_in%iceFrac(:nY) State_Met(LCHNK)%FRLAKE (1,:nY) = 0.0e+0_fp - State_Met(LCHNK)%FRLANDIC (1,:nY) = 0.0e+0_fp - State_Met(LCHNK)%FRSNO (1,:nY) = 0.0e+0_fp + State_Met(LCHNK)%FRLANDICE (1,:nY) = 0.0e+0_fp + State_Met(LCHNK)%FRSNOW (1,:nY) = 0.0e+0_fp ! Field : GWETROOT, GWETTOP ! Description: Root and top soil moisture @@ -3188,11 +3286,11 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) ! isIce, which are based on albedo. Rather, we use CLM landFranc, ocnFrac ! and iceFrac. We also compute isSnow DO J = 1, nY - iMaxLoc = MAXLOC( (/ State_Met(LCHNK)%FRLAND(1,J) + & - State_Met(LCHNK)%FRLANDIC(1,J) + & - State_Met(LCHNK)%FRLAKE(1,J), & - State_Met(LCHNK)%FRSEAICE(1,J), & - State_Met(LCHNK)%FROCEAN(1,J) - & + iMaxLoc = MAXLOC( (/ State_Met(LCHNK)%FRLAND(1,J) + & + State_Met(LCHNK)%FRLANDICE(1,J) + & + State_Met(LCHNK)%FRLAKE(1,J), & + State_Met(LCHNK)%FRSEAICE(1,J), & + State_Met(LCHNK)%FROCEAN(1,J) - & State_Met(LCHNK)%FRSEAICE(1,J) /) ) IF ( iMaxLoc(1) == 3 ) iMaxLoc(1) = 0 ! reset ocean to 0 @@ -3217,7 +3315,7 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) State_Met(LCHNK)%isSnow(1,J) = & ( State_Met(LCHNK)%FRSEAICE(1,J) > 0.0e+0_fp & - .or. State_Met(LCHNK)%SNODP(1,J) > 0.01 ) + .or. State_Met(LCHNK)%SNODP(1,J) > 0.01_fp ) ENDDO @@ -3388,6 +3486,7 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) ENDIF ENDIF + call t_stopf('GEOSChem_CAM_Interfacing') ! This is not necessary as we prescribe CH4 surface mixing ratios ! through CAM. @@ -3483,6 +3582,8 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) ! Thibaud M. Fritz - 27 Feb 2020 !================================================================== + call t_startf('GEOSChem_DryDep') + IF ( Input_Opt%LDryD ) THEN ! Compute the Olson landmap fields of State_Met ! (e.g. State_Met%IREG, State_Met%ILAND, etc.) @@ -3580,7 +3681,11 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) ENDIF !=========================================================== - ! ***** M I X E D L A Y E R M I X I N G ***** + ! ***** S U R F A C E F L U X E S ***** + ! Note: Turbulence (PBL mixing) is NOT done by GEOS-Chem routines + ! and is handled by CAM. But we reuse GEOS-Chem code here to compute + ! the surface *deposition-only* fluxes (-dflx) to merge with the CAM + ! fluxes passed to turbulence. (hplin, 4/30/24) !=========================================================== ! Updates from Bob Yantosca, 06/2020 @@ -3605,23 +3710,30 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) ENDIF ENDIF + ! This dry deposition timer intentionally ends after Compute_Sflx_For_Vdiff + ! because the SurfaceFlux is only the GEOS-Chem deposition flux. (hplin, 4/30/24) + call t_stopf('GEOSChem_DryDep') + !----------------------------------------------------------------------- ! Get emissions from HEMCO + Lightning + Fire ! Add surface emissions to cam_in !----------------------------------------------------------------------- + call t_startf('GEOSChem_Emissions') CALL GC_Emissions_Calc( state = state, & hco_pbuf2d = hco_pbuf2d, & State_Met = State_Met(LCHNK), & cam_in = cam_in, & eflx = eflx, & iStep = iStep ) + call t_stopf('GEOSChem_Emissions') !----------------------------------------------------------------------- - ! Add dry deposition flux + ! Add dry deposition flux from GEOS-Chem State_Chm%SurfaceFlux ! (stored as SurfaceFlux = -dflx) !----------------------------------------------------------------------- + call t_startf('GEOSChem_DryDep') IF ( Input_Opt%LDryD ) THEN DO ND = 1, State_Chm(BEGCHUNK)%nDryDep ! Get the species ID from the drydep ID @@ -3635,19 +3747,23 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) + State_Chm(LCHNK)%SurfaceFlux(1,1:nY,N) ENDDO ENDIF + call t_stopf('GEOSChem_DryDep') !----------------------------------------------------------------------- ! Add non-surface emissions !----------------------------------------------------------------------- + call t_startf('GEOSChem_Emissions') + ! Use units of kg/m2 as State_Chm%Species to add emissions fluxes - CALL Convert_Spc_Units( Input_Opt = Input_Opt, & - State_Chm = State_Chm(LCHNK), & - State_Grid = State_Grid(LCHNK), & - State_Met = State_Met(LCHNK), & - OutUnit = 'kg/m2', & - RC = RC, & - OrigUnit = OrigUnit ) + CALL Convert_Spc_Units( Input_Opt = Input_Opt, & + State_Chm = State_Chm(LCHNK), & + State_Grid = State_Grid(LCHNK), & + State_Met = State_Met(LCHNK), & + new_units = KG_SPECIES_PER_M2, & + previous_units = previous_units, & + RC = RC ) + IF ( RC /= GC_SUCCESS ) THEN ErrMsg = 'Error encountered in "Convert_Spc_Units"!' @@ -3674,7 +3790,7 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) State_Chm = State_Chm(LCHNK), & State_Grid = State_Grid(LCHNK), & State_Met = State_Met(LCHNK), & - OutUnit = OrigUnit, & + new_units = previous_units, & RC = RC ) ! Convert State_Chm%Species back to original units @@ -3683,13 +3799,15 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) CALL Error_Stop( ErrMsg, ThisLoc ) ENDIF + call t_stopf('GEOSChem_Emissions') + !============================================================== ! ***** C H E M I S T R Y ***** !============================================================== call t_startf( 'chemdr' ) - ! Get the overhead column O3 for use with FAST-J + ! Get the overhead column O3 for computing J-values IF ( Input_Opt%Its_A_FullChem_Sim .OR. & Input_Opt%Its_An_Aerosol_Sim ) THEN @@ -3757,8 +3875,10 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) iSfcMrObj => iSfcMrObj%Next ENDDO + call t_startf('GEOSChem_Fullchem') + ! Reset photolysis rates - ZPJ = 0.0e+0_r8 + State_Chm(LCHNK)%Phot%ZPJ = 0.0e+0_r8 ! Perform chemistry CALL Do_Chemistry( Input_Opt = Input_Opt, & @@ -3773,6 +3893,8 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) CALL Error_Stop( ErrMsg, ThisLoc ) ENDIF + call t_stopf('GEOSChem_Fullchem') + ! GEOS-Chem considers CO2 as a dead species and resets its concentration ! internally. Right after the call to `Do_Chemistry`, State_Chm%Species(iCO2) ! corresponds to the chemically-produced CO2. The real CO2 concentration @@ -3780,10 +3902,10 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) State_Chm(LCHNK)%Species(iCO2)%Conc(1,:nY,:nZ) = State_Chm(LCHNK)%Species(iCO2)%Conc(1,:nY,:nZ) & + MMR_Beg(:nY,:nZ,iCO2) - ! Make sure State_Chm(LCHNK) is back in kg/kg dry! - IF ( TRIM(State_Chm(LCHNK)%Spc_Units) /= 'kg/kg dry' ) THEN - Write(iulog,*) 'Current unit = ', TRIM(State_Chm(LCHNK)%Spc_Units) - Write(iulog,*) 'Expected unit = kg/ kg dry' + ! Make sure State_Chm(LCHNK) is back in kg/kg dry! Just check first species. + IF ( State_Chm(LCHNK)%Species(1)%Units /= KG_SPECIES_PER_KG_DRY_AIR ) THEN + Write(iulog,*) 'Current unit = ', TRIM(UNIT_STR(State_Chm(LCHNK)%Species(1)%Units)) + Write(iulog,*) 'Expected unit = ', TRIM(UNIT_STR(KG_SPECIES_PER_KG_DRY_AIR)) CALL ENDRUN('Incorrect unit in GEOS-Chem State_Chm%Species') ENDIF @@ -3800,7 +3922,7 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) CALL pbuf_get_field(pbuf_chnk, tmpIdx, pbuf_i) ! RXN_NO2: NO2 + hv --> NO + O - pbuf_i(:nY) = ZPJ(1,RXN_NO2,1,:nY) + pbuf_i(:nY) = State_Chm(LCHNK)%Phot%ZPJ(1,State_Chm(LCHNK)%Phot%RXN_NO2,1,:nY) pbuf_chnk => NULL() pbuf_i => NULL() @@ -3815,7 +3937,7 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) CALL pbuf_get_field(pbuf_chnk, tmpIdx, pbuf_i) ! RXN_O3_1: O3 + hv --> O2 + O - pbuf_i(:nY) = ZPJ(1,RXN_O3_1,1,:nY) + pbuf_i(:nY) = State_Chm(LCHNK)%Phot%ZPJ(1,State_Chm(LCHNK)%Phot%RXN_O3_1,1,:nY) pbuf_chnk => NULL() pbuf_i => NULL() ENDIF @@ -3837,25 +3959,134 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) !============================================================== #if defined( MODAL_AERO ) - ! Repartition SO4 into H2SO4 and so4_a* - IF ( l_H2SO4 > 0 .AND. l_SO4 > 0 ) THEN - P = l_H2SO4 - ! SO4_gasRatio is mol(SO4) (gaseous) / mol(SO4) (gaseous+aerosol) - vmr1(:nY,:nZ,P) = SO4_gasRatio(:nY,:nZ) * vmr1(:nY,:nZ,l_SO4) - ! binRatio is mol(SO4) (current bin) / mol(SO4) (all bins) - DO M = 1, ntot_amode + call t_startf('GEOSChem_MAM_Interfacing') + ! Construct dgncur_a array for the dry geometric mean diameter [m] + ! of given number distribution. (hplin, 3/6/23) + ! Requires a pbuf field DGNUM + call pbuf_get_field(pbuf, pbuf_get_index('DGNUM'), dgncur_a) + + ! Amount of chemically-produced H2SO4 (mol/mol) + ! This is archived from fullchem_mod.F90 using SO2 + OH rate from KPP (hplin, 1/25/23) + del_h2so4_gasprod(:nY,:nZ) = State_Chm(LCHNK)%H2SO4_PRDR(1,:nY,nZ:1:-1) + + ! Prescribe aerosol size distribution using the method in Feng et al., 2021? + ! If yes, then vmr1 is overwritten with prescribed values for sulfate distribution + ! rather than using values from within MAM4. + ! + ! There are two approaches to where to put in this prescribed value. + ! Either before aero_model_gasaerexch, simulating the change in distribution as chem fluxes (1) + ! or after, directly overwriting any distribution information coming out of aero_model_gasaerexch (2). + ! + ! Testing needs to be done to see if 1 or 2 is 'better' (more reasonable vs. CC and obs.) + if(usePrescribedAerDistribution) then + ! Assume all chemically-produced SO4 is in the gas-phase + ! and use MAM4 to partition into a1, a2, a3. Later, the bins can be prescribed. (hplin, 3/16/23) + ! This will allow total sulfur to be correctly conserved and H2SO4 partitioning is still held. + vmr1(:nY,:nZ,l_H2SO4) = vmr0(:nY,:nZ,l_H2SO4) + State_Chm(LCHNK)%H2SO4_PRDR(1,:nY,nZ:1:-1) + DO M = 1, ntot_amode N = lptr_so4_a_amode(M) IF ( N <= 0 ) CYCLE P = mapCnst(N) - vmr1(:nY,:nZ,P) = vmr1(:nY,:nZ,l_SO4) & - * ( 1.0_r8 - SO4_gasRatio(:nY,:nZ) ) & - * binRatio(iSulf(M),M,:nY,:nZ) - ENDDO + vmr1(:nY,:nZ,P) = vmr0(:nY,:nZ,P) + ENDDO + + ! sect02_new only takes dlo and dhi and auto-computes bin spacing. + ! it may be necessary to prescribe the bin spacing as well to fit a1-3 definitions + ! D < 0.05 is aitken mode, 0.05 to 2 is accumulation mode, D > 2um is coarse. + + prescr_aer_lbnd = 0.0_r8 + prescr_aer_ubnd = 10.0_r8 + + ! First, sum all the available VMRs + vmr_so4_sum(:nY,:nZ) = vmr0(:nY,:nZ,mapCnst(lptr_so4_a_amode(1))) + & + vmr0(:nY,:nZ,mapCnst(lptr_so4_a_amode(2))) + & + vmr0(:nY,:nZ,mapCnst(lptr_so4_a_amode(3)) ) + + ! Loop through chunks as geometric mean dia is different... + DO J = 1, nY + DO L = 1, nZ + ! Get geometric mean diameters of all aerosol bins before further MAM4 calculation + prescr_aer_abnd = (dgncur_a(J, L, 1) + dgncur_a(J, L, 2)) * 1e6_r8 / 2.0_r8 + prescr_aer_cbnd = (dgncur_a(J, L, 1) + dgncur_a(J, L, 3)) * 1e6_r8 / 2.0_r8 + + ! thus, sect02_mam4 is developed for this use. + call sect02_mam4(dgnum_um = 0.14_r8, & ! SO4 geometric mean dia. of log-normal distr [um] + sigmag = 1.6_r8, & ! sigma + duma = 1.0_r8, & ! (unknown scaling factor) + nbin = 3, & ! put into a1, a2, a3 three 'bins' + dlo_sect = (/prescr_aer_lbnd, prescr_aer_abnd, prescr_aer_cbnd/), & ! diameter bounds [um] + dhi_sect = (/prescr_aer_abnd, prescr_aer_cbnd, prescr_aer_ubnd/), & ! diameter bounds [um] + !dlo_sect = (/0.0390625_r8, 0.1_r8, 2.0_r8/), & ! diameter bounds [um] + !dhi_sect = (/0.1_r8, 2.0_r8, 10.0_r8/), & ! diameter bounds [um] + !dlo_um = 0.0390625_r8, & ! lower bound um + !dhi_um = 10.0_r8, & ! upper bound um + xnum_sect= prescr_aer_xnum, & ! prescribed aerosol number ratios + xmas_sect= prescr_aer_xmas) ! prescribed aerosol mass ratios + + ! apply the ratios into the distribution + ! (currently there is no good way to allocate the numbers - thus, this process + ! must be done as "chemical fluxes" and ask gasaerexch to do two hard things here: + ! - partition H2SO4 into the total aer-phase sulfate + ! - readjust num_aX according to changes in SO4 fluxes?) TBD hplin 5/8/23 + + ! use prescr_aer_xnum to scale since we are dealing with mol/mol and not mass quantities. + + ! so4_a3 (coarse mode) + ! ensure that total num molec is conserved. otherwise, this will be a silent sulfate sink. + ! note that prescr_xnum(2) is large, so for maximum precision, use 1.0 to minus it first. + prescr_aer_xnum(3) = (1.0_r8 - prescr_aer_xnum(2)) - prescr_aer_xnum(1) + prescr_aer_xmas(3) = (1.0_r8 - prescr_aer_xmas(2)) - prescr_aer_xmas(1) + + ! there may also be the case that presc_aer_x(3) is lower than 0 (!) which would be unphysical. + ! for safety sake, ensure that coarse sulfate ratio is no lower than 0.01. + ! we can compensate from x(2) (accumulation mode) which is generally the greatest. + if(prescr_aer_xnum(3) .lt. 0.01_r8) then + prescr_aer_xnum(3) = 0.01_r8 + prescr_aer_xnum(2) = 1.0_r8 - prescr_aer_xnum(3) - prescr_aer_xnum(1) + endif + + if(prescr_aer_xmas(3) .lt. 0.01_r8) then + prescr_aer_xmas(3) = 0.01_r8 + prescr_aer_xmas(2) = 1.0_r8 - prescr_aer_xmas(3) - prescr_aer_xmas(1) + endif + + ! so4_a2 (aitken mode) -- note this is smallest even though its a2! + vmr1(:nY,:nZ,mapCnst(lptr_so4_a_amode(2))) = vmr_so4_sum(:nY,:nZ) * prescr_aer_xnum(1) + + ! so4_a1 (accumulation mode) + vmr1(:nY,:nZ,mapCnst(lptr_so4_a_amode(1))) = vmr_so4_sum(:nY,:nZ) * prescr_aer_xnum(2) + + vmr1(:nY,:nZ,mapCnst(lptr_so4_a_amode(3))) = vmr_so4_sum(:nY,:nZ) * prescr_aer_xnum(3) + + ! write out? + ! if(masterproc .and. J .eq. 1 .and. L .eq. nZ) then + ! write(iulog,*) "prescribe aer (so4): L = ", L, " dgncur_a (um) = ", dgncur_a(J, L, :) * 1e6_r8 + ! write(iulog,*) "prescribe aer (so4): L = ", L, " prescr_aer_xnum = ", prescr_aer_xnum + ! write(iulog,*) "prescribe aer (so4): L = ", L, " prescr_aer_xmas = ", prescr_aer_xmas + ! endif + ENDDO + ENDDO + ELSE + ! Original approach: no prescribing aerosol size distribution. + ! Repartition SO4 into H2SO4 and so4_a* + IF ( l_H2SO4 > 0 .AND. l_SO4 > 0 ) THEN + P = l_H2SO4 + ! SO4_gasRatio is mol(SO4) (gaseous) / mol(SO4) (gaseous+aerosol) + vmr1(:nY,:nZ,P) = SO4_gasRatio(:nY,:nZ) * vmr1(:nY,:nZ,l_SO4) + ! binRatio is mol(SO4) (current bin) / mol(SO4) (all bins) + DO M = 1, ntot_amode + N = lptr_so4_a_amode(M) + IF ( N <= 0 ) CYCLE + P = mapCnst(N) + vmr1(:nY,:nZ,P) = vmr1(:nY,:nZ,l_SO4) & + * ( 1.0_r8 - SO4_gasRatio(:nY,:nZ) ) & + * binRatio(iSulf(M),M,:nY,:nZ) + ENDDO + ENDIF ENDIF - ! Amount of chemically-produced H2SO4 (mol/mol) - ! This is archived from fullchem_mod.F90 using SO2 + OH rate from KPP (hplin, 1/25/23) - del_h2so4_gasprod(:nY,:nZ) = State_Chm(LCHNK)%H2SO4_PRDR(1,:nY,nZ:1:-1) + call t_stopf('GEOSChem_MAM_Interfacing') + call t_startf('GEOSChem_MAM_GasAerExch') call aero_model_gasaerexch( loffset = iFirstCnst - 1, & ncol = NCOL, & @@ -3880,6 +4111,9 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) vmr = vmr1, & pbuf = pbuf ) + call t_stopf('GEOSChem_MAM_GasAerExch') + call t_startf('GEOSChem_MAM_Interfacing') + ! Repartition MAM SOAs following mapping: ! TSOA0 + ASOAN + SOAIE + SOAGX -> soa1_a* + soa2_a* ! TSOA1 + ASOA1 -> soa3_a* @@ -3987,10 +4221,12 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) ENDDO #endif + call t_stopf('GEOSChem_MAM_Interfacing') !============================================================== ! ***** W E T D E P O S I T I O N (rainout + washout) ***** !============================================================== + call t_startf('GEOSChem_Neu_Wetdep') IF ( Input_Opt%LWetD ) THEN IF ( gas_wetdep_method == 'NEU' ) THEN @@ -4014,6 +4250,7 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) ENDIF ENDIF + call t_stopf('GEOSChem_Neu_Wetdep') !============================================================== ! ***** B O U N D A R Y C O N D I T I O N S ***** @@ -4036,10 +4273,10 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) adv_mass(N) / MWDry ENDDO - ! Make sure State_Chm(LCHNK) is back in kg/kg dry! - IF ( TRIM(State_Chm(LCHNK)%Spc_Units) /= 'kg/kg dry' ) THEN - Write(iulog,*) 'Current unit = ', TRIM(State_Chm(LCHNK)%Spc_Units) - Write(iulog,*) 'Expected unit = kg/ kg dry' + ! Make sure State_Chm(LCHNK) is back in kg/kg dry! Only check first species. + IF ( State_Chm(LCHNK)%Species(1)%Units /= KG_SPECIES_PER_KG_DRY_AIR ) THEN + Write(iulog,*) 'Current unit = ', TRIM(UNIT_STR(State_Chm(LCHNK)%Species(1)%Units)) + Write(iulog,*) 'Expected unit = ', TRIM(UNIT_STR(KG_SPECIES_PER_KG_DRY_AIR)) CALL ENDRUN('Incorrect unit in GEOS-Chem State_Chm%Species') ENDIF @@ -4068,6 +4305,7 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) ENDDO #if defined( MODAL_AERO ) + call t_startf('GEOSChem_MAM_Interfacing') ! Here apply tendencies to MAM aerosols ! Initial mass in bin SM is stored as state%q(N) ! Final mass in bin SM is stored as binRatio(SM,M) * State_Chm(P) @@ -4118,6 +4356,7 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) * adv_mass(P) / MWDry ENDIF ENDDO + call t_stopf('GEOSChem_MAM_Interfacing') #endif DO N = 1, gas_pcnst @@ -4136,6 +4375,7 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) ptend%q(:,:,cQ) = ptend%q(:,:,cH2O) ENDIF + call t_startf('GEOSChem_Diagnostics') CALL GC_Diagnostics_Calc( Input_Opt = Input_Opt, & State_Chm = State_Chm(LCHNK), & State_Diag = State_Diag(LCHNK), & @@ -4179,6 +4419,7 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) HistoryConfig = HistoryConfig, & LCHNK = LCHNK, & RC = RC ) + call t_stopf('GEOSChem_Diagnostics') IF ( ghg_chem ) THEN ptend%lq(1) = .True. @@ -4216,6 +4457,8 @@ subroutine chem_timestep_tend( state, ptend, cam_in, cam_out, dT, pbuf, fh2o ) FIRST = .false. ENDIF + call t_stopf('GEOSChem_All_Tend') + end subroutine chem_timestep_tend !================================================================================================ @@ -4264,9 +4507,7 @@ subroutine chem_final use geoschem_history_mod, only : Destroy_HistoryConfig ! GEOS-Chem modules - use Aerosol_Mod, only : Cleanup_Aerosol use Carbon_Mod, only : Cleanup_Carbon - use CMN_FJX_Mod, only : Cleanup_CMN_FJX use Drydep_Mod, only : Cleanup_Drydep use Dust_Mod, only : Cleanup_Dust use Error_Mod, only : Cleanup_Error @@ -4289,7 +4530,6 @@ subroutine chem_final ! Finalize GEOS-Chem - CALL Cleanup_Aerosol CALL Cleanup_Carbon CALL Cleanup_Drydep CALL Cleanup_Dust @@ -4308,12 +4548,6 @@ subroutine chem_final CALL short_lived_species_final() - CALL Cleanup_CMN_FJX( RC ) - IF ( RC /= GC_SUCCESS ) THEN - ErrMsg = 'Error encountered in "Cleanup_CMN_FJX"!' - CALL Error_Stop( ErrMsg, ThisLoc ) - ENDIF - ! Cleanup Input_Opt CALL Cleanup_Input_Opt( Input_Opt, RC ) @@ -4435,5 +4669,111 @@ subroutine chem_emissions( state, cam_in, pbuf ) ENDDO end subroutine chem_emissions - +! +! P R E S C R I B E A E R O S O L D I S T R I B U T I O N +! +! Based on code from Feng et al., 2021 GMD (WRF-GC v2.0), by Xu Feng et al. +! in module_diag_aero_size_info.F, originally based from WRF-Chem. +! +! Reference: +! Feng, X., Lin, H., Fu, T.-M., Sulprizio, M. P., Zhuang, J., Jacob, D. J., Tian, H., Ma, Y., Zhang, L., Wang, X., Chen, Q., and Han, Z.: WRF-GC (v2.0): online two-way coupling of WRF (v3.9.1.1) and GEOS-Chem (v12.7.2) for modeling regional atmospheric chemistry–meteorology interactions, Geosci. Model Dev., 14, 3741–3768, https://doi.org/10.5194/gmd-14-3741-2021, 2021. +! + + real(8) function erfc_num_recipes( x ) + ! + ! from press et al, numerical recipes, 1990, page 164 + ! + implicit none + real(r8) :: x, erfc_dbl, dum, t, z + z = abs(x) + t = 1.0_r8/(1.0_r8 + 0.5_r8*z) + dum = ( -z*z - 1.26551223_r8 + t*(1.00002368_r8 + t*(0.37409196_r8 + & + t*(0.09678418_r8 + t*(-0.18628806_r8 + t*(0.27886807_r8 + & + t*(-1.13520398_r8 + & + t*(1.48851587_r8 + t*(-0.82215223_r8 + t*0.17087277_r8 ))))))))) + erfc_dbl = t * exp(dum) + if (x .lt. 0.0_r8) erfc_dbl = 2.0_r8 - erfc_dbl + erfc_num_recipes = erfc_dbl + return + end function erfc_num_recipes + + ! sect02_mam4 is based off sect02_new in WRF-GC, which is based off + ! sect02 in WRF-Chem chem/module_optical_averaging.F. + ! + ! user specifies a single log-normal mode and a set of section boundaries + ! prog calculates mass and number for each section. + subroutine sect02_mam4(dgnum_um, sigmag, duma, nbin, dlo_sect, dhi_sect, & + xnum_sect, xmas_sect) + ! INPUT PARAMETERS: + ! dgnum_um *diameter* geometric mean of log-normal distribution [um] + ! sigmag geometric standard deviation of log-normal dist. [unitless] + ! duma 1.0 ? + ! nbin # of target bins (wrf-gc = 4, MAM4 = 3) [count] + ! dlo_sect(nbin) low diameter limit (wrf-gc = 0.0390625) [um] + ! dhi_sect(nbin) high diameter limit (wrf-gc = 10.0) [um] + + ! OUTPUT PARAMETERS: + ! xnum_sect(nbin) aerosol number per bin, ratio of total [unitless] + ! xmas_sect(bin) aerosol mass per bin, ratio of total [unitless] + + implicit none + real(8), dimension(nbin), intent(out) :: xnum_sect, xmas_sect + integer :: n, nbin + real(8) :: dgnum, dgnum_um, dhi, & + dlo, duma, dumfrac, & + dx, sigmag, & + sx, sxroot2, thi, tlo, x0, x3, & + xhi, xlo, xmtot, xntot + real(8), intent(in) :: dlo_sect(nbin), dhi_sect(nbin) + real(8) :: my_dlo_sect(nbin), my_dhi_sect(nbin) + real(8) :: pi + parameter (pi = 3.141592653589_r8) + + xmtot = duma + xntot = duma + + ! Compute bins based on number of bins. Originally sect02_new. + ! For MAM4, we prescribe the bin ranges as well. + ! dlo = dlo_um*1.0E-4_r8 + ! dhi = dhi_um*1.0E-4_r8 + ! xlo = log( dlo ) + ! xhi = log( dhi ) + ! dx = (xhi - xlo)/nbin + ! do n = 1, nbin + ! dlo_sect(n) = exp( xlo + dx*(n-1) ) + ! dhi_sect(n) = exp( xlo + dx*n ) + ! end do + + ! dlo_sect and dhi_sect have to be scaled by 1e-4 + ! in order to fit parameters in the above calculation, if they are prescribed. + + my_dlo_sect(:) = dlo_sect(:) * 1.0e-4_r8 + my_dhi_sect(:) = dhi_sect(:) * 1.0e-4_r8 + + dgnum = dgnum_um*1.0E-4_r8 + sx = log( sigmag ) + x0 = log( dgnum ) + x3 = x0 + 3.0_r8*sx*sx + sxroot2 = sx * sqrt( 2.0_r8 ) + do n = 1, nbin + xlo = log( my_dlo_sect(n) ) + xhi = log( my_dhi_sect(n) ) + tlo = (xlo - x0)/sxroot2 + thi = (xhi - x0)/sxroot2 + if (tlo .le. 0.0_r8) then + dumfrac = 0.5_r8*( erfc_num_recipes(-thi) - erfc_num_recipes(-tlo) ) + else + dumfrac = 0.5_r8*( erfc_num_recipes(tlo) - erfc_num_recipes(thi) ) + end if + xnum_sect(n) = xntot*dumfrac + tlo = (xlo - x3)/sxroot2 + thi = (xhi - x3)/sxroot2 + if (tlo .le. 0.0_r8) then + dumfrac = 0.5_r8*( erfc_num_recipes(-thi) - erfc_num_recipes(-tlo) ) + else + dumfrac = 0.5_r8*( erfc_num_recipes(tlo) - erfc_num_recipes(thi) ) + endif + xmas_sect(n) = xmtot*dumfrac + enddo + end subroutine sect02_mam4 end module chemistry diff --git a/src/chemistry/geoschem/geoschem_diagnostics_mod.F90 b/src/chemistry/geoschem/geoschem_diagnostics_mod.F90 index 447d2c29cd..5da738038c 100644 --- a/src/chemistry/geoschem/geoschem_diagnostics_mod.F90 +++ b/src/chemistry/geoschem/geoschem_diagnostics_mod.F90 @@ -56,10 +56,10 @@ MODULE GeosChem_Diagnostics_Mod REAL(r8) :: NHx_MWs(2) REAL(r8) :: TOTH_MWs(3) - REAL(r8), PARAMETER :: MW_NIT = 62.01 - REAL(r8), PARAMETER :: MW_HNO3 = 63.01 - REAL(r8), PARAMETER :: MW_HCl = 36.45 - REAL(r8), PARAMETER :: MW_H2O = 18.02 + REAL(r8), PARAMETER :: MW_NIT = 62.01_r8 + REAL(r8), PARAMETER :: MW_HNO3 = 63.01_r8 + REAL(r8), PARAMETER :: MW_HCl = 36.45_r8 + REAL(r8), PARAMETER :: MW_H2O = 18.02_r8 ! NOx species INTEGER :: i_NO, i_NO2, i_N diff --git a/src/chemistry/geoschem/geoschem_emissions_mod.F90 b/src/chemistry/geoschem/geoschem_emissions_mod.F90 index 9d9dfc6bd1..dc674775ba 100644 --- a/src/chemistry/geoschem/geoschem_emissions_mod.F90 +++ b/src/chemistry/geoschem/geoschem_emissions_mod.F90 @@ -432,13 +432,13 @@ SUBROUTINE GC_Emissions_Calc( state, hco_pbuf2d, State_Met, cam_in, eflx, iStep ! Convert from kg/m2/s to molec/cm3/s ! Note 1: cnst_mw is in kg/kmole ! Note 2: avogad is in molecules/kmole - CALL Outfld( TRIM(SpcName), eflx(:nY,:nZ,N) / State_Met%BXHEIGHT(1,:nY,nZ:1:-1) * 1.0E-06 / cnst_mw(N) * avogad, nY, LCHNK ) + CALL Outfld( TRIM(SpcName), eflx(:nY,:nZ,N) / State_Met%BXHEIGHT(1,:nY,nZ:1:-1) * 1.0E-06_r8 / cnst_mw(N) * avogad, nY, LCHNK ) SpcName = TRIM(cnst_name(N))//'_CLXF' ! Convert from kg/m2/s to molec/cm2/s ! Note 1: cnst_mw is in kg/kmole ! Note 2: avogad is in molecules/kmole - CALL Outfld( TRIM(SpcName), SUM(eflx(:nY,:nZ-1,N), DIM=2) * 1.0E-04 / cnst_mw(N) * avogad, nY, LCHNK ) + CALL Outfld( TRIM(SpcName), SUM(eflx(:nY,:nZ-1,N), DIM=2) * 1.0E-04_r8 / cnst_mw(N) * avogad, nY, LCHNK ) SpcName = TRIM(cnst_name(N))//'_CMXF' CALL Outfld( TRIM(SpcName), SUM(eflx(:nY,:nZ-1,N), DIM=2), nY, LCHNK ) @@ -454,7 +454,7 @@ SUBROUTINE GC_Emissions_Calc( state, hco_pbuf2d, State_Met, cam_in, eflx, iStep ! Multiply by MWNO * BXHEIGHT * 1.0E+06 / AVO ! = mole/molec * kg NO/mole * m * cm^3/m^3 ! cnst_mw(N) is in g/mole - SCALFAC = cnst_mw(N) * 1.0E-03 * 1.0E+06 / AVO + SCALFAC = cnst_mw(N) * 1.0E-03_r8 * 1.0E+06_r8 / AVO DO J = 1, nY DO L = 1, nZ eflx(J,L,N) = eflx(J,L,N) & diff --git a/src/chemistry/geoschem/geoschem_src b/src/chemistry/geoschem/geoschem_src index 28345ee76e..bef56c605e 160000 --- a/src/chemistry/geoschem/geoschem_src +++ b/src/chemistry/geoschem/geoschem_src @@ -1 +1 @@ -Subproject commit 28345ee76e5631d6d14869a36dc73e9dd6e0ce1e +Subproject commit bef56c605e018eecbd91646a51ce82c7cd77f56a diff --git a/src/chemistry/geoschem/mo_sim_dat.F90 b/src/chemistry/geoschem/mo_sim_dat.F90 index 13f5740645..0d52edc2e3 100644 --- a/src/chemistry/geoschem/mo_sim_dat.F90 +++ b/src/chemistry/geoschem/mo_sim_dat.F90 @@ -36,146 +36,150 @@ subroutine set_sim_dat ! cls_rxt_cnt(:,1) = (/ 37, 61, 0, 30 /) ! cls_rxt_cnt(:,4) = (/ 23, 174, 326, 191 /) - ! GEOS-Chem tracers (advected species) are placed first along MAM - ! aerosols, as those will be constituents. MAM requires that there + ! GEOS-Chem tracers (advected species) are placed before MAM + ! aerosols which are also constituents. MAM requires that there ! is a linear mapping between solsym and constituents - ! ewl notes: added HMS (for GEOS-Chem 13.3) - ! added AONITA, AROMP4, AROMP5, BALD, BENZP, BZCO3H, - ! BZPAN, C2H2, C2H4, CSL, ETHN, ETHP, MCT, NPHEN, PHEN for 14.0 - ! Removed non-advected GEOS-Chem species for 14.0, except CO2 - ! which is a constituent, as well as OH and HO2 for diagnostic - ! output. + ! note: + ! - solsym and adv_mass here, + ! - gas_pcnst, nTracersMax, nSlvd in chem_mods.F90, + ! - chem_nadv in bld/configure (clean all needed after change) + ! and .xml files for dry and wet deposition lists, should be updated + ! programmatically using https://github.com/jimmielin/geos-chem-coupling-tools + ! from a species .yml file generated by a corresponding GEOS-Chem "Classic" simulation. ! - ! Currently include GC advected species (233), MAM aerosols (33), CO2 (1), - ! and OH and HO2 (2). - ! If changed, update to match solsym length: - ! 1. cam/bld/configure variable $chem_nadv - ! 2. cam/src/chemistry/geoschem/chem_mods.F90 vars gas_pcnst and nTracersMax - ! Also update adv_mass to store MWs for species in solsym (ewl, 8/8/22) - solsym(:269) = (/ 'ACET ', & - 'ACTA ','AERI ', & - 'ALD2 ','ALK4 ','ASOA1 ', & - 'ASOA2 ','ASOA3 ','ASOAN ', & - 'ASOG1 ','ASOG2 ','ASOG3 ', & - 'AONITA ','AROMP4 ','AROMP5 ', & - 'ATOOH ','BALD ','BCPI ', & - 'BCPO ','BENZ ','BENZP ', & - 'BR ','BR2 ','BRCL ', & - 'BRNO2 ','BRNO3 ','BRO ', & - 'BRSALA ','BRSALC ','BZCO3H ', & - 'BZPAN ','C2H2 ','C2H4 ', & - 'C2H6 ','C3H8 ','CCL4 ', & - 'CFC11 ','CFC113 ','CFC114 ', & - 'CFC115 ','CFC12 ','CH2BR2 ', & - 'CH2CL2 ','CH2I2 ','CH2IBR ', & - 'CH2ICL ','CH2O ','CH3BR ', & - 'CH3CCL3 ','CH3CL ','CH3I ', & - 'CH4 ','CHBR3 ','CHCL3 ', & - 'CL ','CL2 ','CL2O2 ', & - 'CLNO2 ','CLNO3 ','CLO ', & - 'CLOO ','CLOCK ','CO ', & - 'CSL ', & - 'DMS ','DST1 ','DST2 ', & - 'DST3 ','DST4 ','EOH ', & - 'ETHLN ','ETHN ','ETHP ', & - 'ETNO3 ','ETP ', & - 'GLYC ','GLYX ', & - 'H1211 ','H1301 ','H2402 ', & - 'H2O ','H2O2 ','HAC ', & - 'HBR ','HC5A ','HCFC123 ', & - 'HCFC141B ','HCFC142B ','HCFC22 ', & - 'HCL ','HCOOH ','HI ', & - 'HMHP ','HMML ','HMS ', & - 'HNO2 ', & - 'HNO3 ','HNO4 ','HOBR ', & - 'HOCL ','HOI ','HONIT ', & - 'HPALD1 ','HPALD2 ','HPALD3 ', & - 'HPALD4 ','HPETHNL ','I ', & - 'I2 ','I2O2 ','I2O3 ', & - 'I2O4 ','IBR ','ICHE ', & - 'ICL ','ICN ','ICPDH ', & - 'IDC ','IDCHP ','IDHDP ', & - 'IDHPE ','IDN ','IEPOXA ', & - 'IEPOXB ','IEPOXD ','IHN1 ', & - 'IHN2 ','IHN3 ','IHN4 ', & - 'INDIOL ','INO ','INPB ', & - 'INPD ','IO ','IONITA ', & - 'IONO ','IONO2 ','IPRNO3 ', & - 'ISALA ','ISALC ','ISOP ', & - 'ITCN ','ITHN ','LIMO ', & - 'LVOC ','LVOCOA ','MACR ', & - 'MACR1OOH ','MAP ','MCRDH ', & - 'MCRENOL ','MCRHN ','MCRHNB ', & - 'MCRHP ','MCT ','MEK ', & - 'MENO3 ', & - 'MGLY ','MOH ','MONITA ', & - 'MONITS ','MONITU ','MP ', & - 'MPAN ','MPN ','MSA ', & - 'MTPA ','MTPO ','MVK ', & - 'MVKDH ','MVKHC ','MVKHCB ', & - 'MVKHP ','MVKN ','MVKPC ', & - 'N2O ','N2O5 ','NH3 ', & - 'NH4 ','NIT ','NITS ', & - 'NO ','NO2 ','NO3 ', & - 'NPHEN ', & - 'NPRNO3 ','O3 ','OCLO ', & - 'OCPI ','OCPO ','OCS ', & - 'OIO ','PAN ','PFE ', & - 'PHEN ', & - 'PIP ','PP ','PPN ', & - 'PROPNN ','PRPE ','PRPN ', & - 'PYAC ','R4N2 ','R4P ', & - 'RA3P ','RB3P ','RCHO ', & - 'RIPA ','RIPB ','RIPC ', & - 'RIPD ','RP ','SALA ', & - 'SALAAL ','SALACL ','SALC ', & - 'SALCAL ','SALCCL ','SO2 ', & - 'SO4 ','SO4S ','SOAGX ', & - 'SOAIE ','TOLU ','TSOA0 ', & - 'TSOA1 ','TSOA2 ','TSOA3 ', & - 'TSOG0 ','TSOG1 ','TSOG2 ', & - 'TSOG3 ','XYLE ','bc_a1 ', & - 'bc_a4 ','dst_a1 ','dst_a2 ', & - 'dst_a3 ','ncl_a1 ','ncl_a2 ', & - 'ncl_a3 ','num_a1 ','num_a2 ', & - 'num_a3 ','num_a4 ','pom_a1 ', & - 'pom_a4 ','so4_a1 ','so4_a2 ', & - 'so4_a3 ','soa1_a1 ','soa1_a2 ', & - 'soa2_a1 ','soa2_a2 ','soa3_a1 ', & - 'soa3_a2 ','soa4_a1 ','soa4_a2 ', & - 'soa5_a1 ','soa5_a2 ','H2SO4 ', & - 'SOAG0 ','SOAG1 ','SOAG2 ', & - 'SOAG3 ','SOAG4 ','CO2 ', & - 'HO2 ','OH ' /) -!non-advected GEOS-Chem species in 14.0 (beware this includes OH and HO2 already listed above) -! 'LBRO2H ','LBRO2N ','LISOPOH ', & -! 'LISOPNO3 ','LTRO2H ','LTRO2N ', & -! 'LXRO2H ','LXRO2N ','SO4H1 ', & -! 'SO4H2 ','SO4H3 ','SO4H4 ', & -! 'POX ','LOX ','PCO ', & -! 'LCO ','PSO4 ','LCH4 ', & -! 'PH2O2 ','BRO2 ','TRO2 ', & -! 'N ','XRO2 ','HPALD2OO ', & -! 'HPALD1OO ','INA ','C4HVP1 ', & -! 'C4HVP2 ','IDNOO ','ICNOO ', & -! 'ISOPNOO2 ','ROH ','ISOPNOO1 ', & -! 'IDHNDOO1 ','IDHNDOO2 ','H ', & -! 'IHPOO2 ','IHPOO1 ','IHPOO3 ', & -! 'IHPNDOO ','ICHOO ','R4N1 ', & -! 'PRN1 ','MVKOHOO ','MCROHOO ', & -! 'MACR1OO ','PO2 ','OLNN ', & -! 'OLND ','ETO2 ','IHPNBOO ', & -! 'RCO3 ','LIMO2 ','KO2 ', & -! 'IEPOXAOO ','IEPOXBOO ','CH3CHOO ', & -! 'PIO2 ','IDHNBOO ','A3O2 ', & -! 'IHOO4 ','IHOO1 ','INO2D ', & -! 'INO2B ','MACRNO2 ','ATO2 ', & -! 'OTHRO2 ','R4O2 ','B3O2 ', & -! 'CH2OO ','MCO3 ','MO2 ', & -! 'O1D ','OH ','HO2 ', & -! 'O ','H2 ','N2 ', & -! 'O2 ','RCOOH ' /) + ! editing of these lists manually may be error prone. + ! solsym contains, in order: + ! - GEOS-Chem advected species. The order of these must match exactly as specified in + ! geoschem_config.yml. + ! - MAM aerosols (33 in MAM4; hardcoded in chem_readnl of geoschem/chemistry.F90) + ! these are the 4 modes for bc (black carbon), dst (dust), ncl (NaCl, seasalt), + ! pom (primary organic matter), SOA, SO4 (sulfate), and + ! number concentration, + ! plus H2SO4, SOAG0-4, for MAM4. + ! - CO2, and + ! - GEOS-Chem non-advected species, as specified in gckpp_Parameters.F90. (hplin, 6/24/24) + solsym(:357) = (/ & + 'ACET ', 'ACTA ', 'AERI ', & + 'ALD2 ', 'ALK4 ', 'AONITA ', & + 'AROMP4 ', 'AROMP5 ', 'ASOA1 ', & + 'ASOA2 ', 'ASOA3 ', 'ASOAN ', & + 'ASOG1 ', 'ASOG2 ', 'ASOG3 ', & + 'ATOOH ', 'BALD ', 'BCPI ', & + 'BCPO ', 'BENZ ', 'BENZP ', & + 'BR ', 'BR2 ', 'BRCL ', & + 'BRNO2 ', 'BRNO3 ', 'BRO ', & + 'BRSALA ', 'BRSALC ', 'BUTDI ', & + 'BZCO3H ', 'BZPAN ', 'C2H2 ', & + 'C2H4 ', 'C2H6 ', 'C3H8 ', & + 'CCL4 ', 'CFC11 ', 'CFC113 ', & + 'CFC114 ', 'CFC115 ', 'CFC12 ', & + 'CH2BR2 ', 'CH2CL2 ', 'CH2I2 ', & + 'CH2IBR ', 'CH2ICL ', 'CH2O ', & + 'CH3BR ', 'CH3CCL3 ', 'CH3CL ', & + 'CH3I ', 'CH4 ', 'CHBR3 ', & + 'CHCL3 ', 'CL ', 'CL2 ', & + 'CL2O2 ', 'CLNO2 ', 'CLNO3 ', & + 'CLO ', 'CLOO ', 'CLOCK ', & + 'CO ', 'CSL ', 'DMS ', & + 'DST1 ', 'DST2 ', 'DST3 ', & + 'DST4 ', 'EOH ', 'ETHLN ', & + 'ETHN ', 'ETHP ', 'ETNO3 ', & + 'ETP ', 'FURA ', 'GLYC ', & + 'GLYX ', 'H1211 ', 'H1301 ', & + 'H2402 ', 'H2O ', 'H2O2 ', & + 'HAC ', 'HBR ', 'HC5A ', & + 'HCFC123 ', 'HCFC141B ', 'HCFC142B ', & + 'HCFC22 ', 'HCL ', 'HCOOH ', & + 'HI ', 'HMHP ', 'HMML ', & + 'HMS ', 'HNO2 ', 'HNO3 ', & + 'HNO4 ', 'HOBR ', 'HOCL ', & + 'HOI ', 'HONIT ', 'HPALD1 ', & + 'HPALD2 ', 'HPALD3 ', 'HPALD4 ', & + 'HPETHNL ', 'I ', 'I2 ', & + 'I2O2 ', 'I2O3 ', 'I2O4 ', & + 'IBR ', 'ICHE ', 'ICL ', & + 'ICN ', 'ICPDH ', 'IDC ', & + 'IDCHP ', 'IDHDP ', 'IDHPE ', & + 'IDN ', 'IEPOXA ', 'IEPOXB ', & + 'IEPOXD ', 'IHN1 ', 'IHN2 ', & + 'IHN3 ', 'IHN4 ', 'INDIOL ', & + 'INO ', 'INPB ', 'INPD ', & + 'IO ', 'IONITA ', 'IONO ', & + 'IONO2 ', 'IPRNO3 ', 'ISALA ', & + 'ISALC ', 'ISOP ', 'ITCN ', & + 'ITHN ', 'LIMO ', 'LVOC ', & + 'LVOCOA ', 'MACR ', 'MACR1OOH ', & + 'MAP ', 'MCRDH ', 'MCRENOL ', & + 'MCRHN ', 'MCRHNB ', 'MCRHP ', & + 'MCT ', 'MEK ', 'MENO3 ', & + 'MGLY ', 'MOH ', 'MONITA ', & + 'MONITS ', 'MONITU ', 'MP ', & + 'MPAN ', 'MPN ', 'MSA ', & + 'MTPA ', 'MTPO ', 'MVK ', & + 'MVKDH ', 'MVKHC ', 'MVKHCB ', & + 'MVKHP ', 'MVKN ', 'MVKPC ', & + 'N2O ', 'N2O5 ', 'NH3 ', & + 'NH4 ', 'NIT ', 'NITS ', & + 'NO ', 'NO2 ', 'NO3 ', & + 'NPHEN ', 'NPRNO3 ', 'O3 ', & + 'OCLO ', 'OCPI ', 'OCPO ', & + 'OCS ', 'OIO ', 'PAN ', & + 'PFE ', 'PHEN ', 'PIP ', & + 'PP ', 'PPN ', 'PROPNN ', & + 'PRPE ', 'PRPN ', 'PYAC ', & + 'R4N2 ', 'R4P ', 'RA3P ', & + 'RB3P ', 'RCHO ', 'RIPA ', & + 'RIPB ', 'RIPC ', 'RIPD ', & + 'RP ', 'SALA ', 'SALAAL ', & + 'SALACL ', 'SALC ', 'SALCAL ', & + 'SALCCL ', 'SO2 ', 'SO4 ', & + 'SO4S ', 'SOAGX ', 'SOAIE ', & + 'TOLU ', 'TSOA0 ', 'TSOA1 ', & + 'TSOA2 ', 'TSOA3 ', 'TSOG0 ', & + 'TSOG1 ', 'TSOG2 ', 'TSOG3 ', & + 'XYLE ', 'bc_a1 ', 'bc_a4 ', & + 'dst_a1 ', 'dst_a2 ', 'dst_a3 ', & + 'ncl_a1 ', 'ncl_a2 ', 'ncl_a3 ', & + 'num_a1 ', 'num_a2 ', 'num_a3 ', & + 'num_a4 ', 'pom_a1 ', 'pom_a4 ', & + 'so4_a1 ', 'so4_a2 ', 'so4_a3 ', & + 'soa1_a1 ', 'soa1_a2 ', 'soa2_a1 ', & + 'soa2_a2 ', 'soa3_a1 ', 'soa3_a2 ', & + 'soa4_a1 ', 'soa4_a2 ', 'soa5_a1 ', & + 'soa5_a2 ', 'H2SO4 ', 'SOAG0 ', & + 'SOAG1 ', 'SOAG2 ', 'SOAG3 ', & + 'SOAG4 ', 'A3O2 ', 'AROMRO2 ', & + 'ATO2 ', 'B3O2 ', 'BENZO ', & + 'BENZO2 ', 'BRO2 ', 'BZCO3 ', & + 'C4HVP1 ', 'C4HVP2 ', 'CH2OO ', & + 'CH3CHOO ', 'CO2 ', 'ETO ', & + 'ETO2 ', 'ETOO ', 'H ', & + 'H2 ', 'HO2 ', 'HPALD1OO ', & + 'HPALD2OO ', 'ICHOO ', 'ICNOO ', & + 'IDHNBOO ', 'IDHNDOO1 ', 'IDHNDOO2 ', & + 'IDNOO ', 'IEPOXAOO ', 'IEPOXBOO ', & + 'IHOO1 ', 'IHOO4 ', 'IHPNBOO ', & + 'IHPNDOO ', 'IHPOO1 ', 'IHPOO2 ', & + 'IHPOO3 ', 'INA ', 'INO2B ', & + 'INO2D ', 'ISOPNOO1 ', 'ISOPNOO2 ', & + 'KO2 ', 'LBRO2H ', 'LBRO2N ', & + 'LCH4 ', 'LCO ', 'LIMO2 ', & + 'LISOPNO3 ', 'LISOPOH ', 'LNRO2H ', & + 'LNRO2N ', 'LOX ', 'LTRO2H ', & + 'LTRO2N ', 'LXRO2H ', 'LXRO2N ', & + 'MACR1OO ', 'MACRNO2 ', 'MCO3 ', & + 'MCROHOO ', 'MO2 ', 'MVKOHOO ', & + 'N ', 'N2 ', 'NAP ', & + 'NRO2 ', 'O ', 'O1D ', & + 'O2 ', 'OH ', 'OLND ', & + 'OLNN ', 'OTHRO2 ', 'PCO ', & + 'PH2O2 ', 'PH2SO4 ', 'PIO2 ', & + 'PO2 ', 'POX ', 'PRN1 ', & + 'PSO4 ', 'PSO4AQ ', 'R4N1 ', & + 'R4O2 ', 'RCO3 ', 'RCOOH ', & + 'ROH ', 'TRO2 ', 'XRO2 ' & + /) inv_lst(: 6) = (/ 'M ', 'N2 ', 'O2 ', & 'H2 ', 'MOH ', 'RCOOH ' /) @@ -183,60 +187,130 @@ subroutine set_sim_dat fix_mass(: 6) = (/ 0.00000000_r8, 28.0134800_r8, 31.9988000_r8, 2.020000_r8, 32.050000_r8, & 74.090000_r8 /) - adv_mass(:269) = (/ 58.090000_r8, 60.060000_r8, 126.900000_r8, 44.060000_r8, 58.120000_r8, & - 150.000000_r8, 150.000000_r8, 150.000000_r8, 150.000000_r8, 150.000000_r8, & - 150.00000_r8, 150.000000_r8, 189.12_r8, 68.08_r8, 98.10_r8, & - 90.0900000_r8, 106.12_r8, 12.010000_r8, 12.010000_r8, & - 78.120000_r8, 110.11_r8, 79.900000_r8, 159.800000_r8, 115.450000_r8, & - 125.910000_r8, 141.910000_r8, 95.900000_r8, 79.900000_r8, & - 79.900000_r8, 138.12_r8, 183.12_r8, 26.05_r8, 28.05_r8, 30.080000_r8, & - 44.110000_r8, 153.820000_r8, 137.370000_r8, 187.380000_r8, 170.920000_r8, & - 154.470000_r8, 120.910000_r8, 173.830000_r8, 84.930000_r8, 267.840000_r8, & - 220.840000_r8, 176.380000_r8, 30.030000_r8, 94.940000_r8, 133.350000_r8, & - 50.450000_r8, 141.940000_r8, 16.050000_r8, 252.730000_r8, 119.350000_r8, & - 35.450000_r8, 70.900000_r8, 102.910000_r8, 81.450000_r8, 97.450000_r8, & - 51.450000_r8, 67.450000_r8, 1.000000_r8, 28.010000_r8, 108.14_r8, 62.130000_r8, & - 29.000000_r8, 29.000000_r8, 29.000000_r8, 29.000000_r8, 46.080000_r8, & - 105.060000_r8, 107.07_r8, 78.07_r8, 91.080000_r8, 62.080000_r8, 60.060000_r8, 58.040000_r8, & - 165.360000_r8, 148.910000_r8, 259.820000_r8, 18.020000_r8, 34.020000_r8, & - 74.080000_r8, 80.910000_r8, 100.130000_r8, 152.930000_r8, 116.940000_r8, & - 100.500000_r8, 86.470000_r8, 36.450000_r8, 46.030000_r8, 127.910000_r8, & - 64.050000_r8, 102.100000_r8, 110.000000_r8, 47.010000_r8, 63.010000_r8, 79.010000_r8, & - 96.910000_r8, 52.450000_r8, 143.890000_r8, 215.000000_r8, 116.130000_r8, & - 116.130000_r8, 116.130000_r8, 116.130000_r8, 76.060000_r8, 126.900000_r8, & - 253.800000_r8, 285.800000_r8, 301.800000_r8, 317.800000_r8, 206.900000_r8, & - 116.130000_r8, 162.450000_r8, 145.130000_r8, 150.150000_r8, 98.110000_r8, & - 148.130000_r8, 168.170000_r8, 150.150000_r8, 192.150000_r8, 106.140000_r8, & - 106.140000_r8, 106.140000_r8, 147.150000_r8, 147.150000_r8, 147.150000_r8, & - 147.150000_r8, 102.000000_r8, 156.910000_r8, 163.150000_r8, 163.150000_r8, & - 142.900000_r8, 14.010000_r8, 172.910000_r8, 188.910000_r8, 105.110000_r8, & - 126.900000_r8, 126.900000_r8, 68.130000_r8, 195.150000_r8, 197.170000_r8, & - 136.260000_r8, 154.190000_r8, 154.190000_r8, 70.100000_r8, 102.100000_r8, & - 76.060000_r8, 104.120000_r8, 86.100000_r8, 149.110000_r8, 149.110000_r8, & - 120.120000_r8, 124.0_r8, 72.110000_r8, 77.050000_r8, 72.070000_r8, 32.050000_r8, & - 14.010000_r8, 215.280000_r8, 215.280000_r8, 48.050000_r8, 147.100000_r8, & - 93.050000_r8, 96.100000_r8, 136.260000_r8, 136.260000_r8, 70.090000_r8, & - 105.130000_r8, 102.100000_r8, 102.100000_r8, 120.120000_r8, 149.120000_r8, & - 118.100000_r8, 44.020000_r8, 108.020000_r8, 17.040000_r8, 18.050000_r8, & - 62.010000_r8, 31.400000_r8, 30.010000_r8, 46.010000_r8, 62.010000_r8, & - 139.11_r8, 105.110000_r8, 48.000000_r8, 67.450000_r8, 12.010000_r8, 12.010000_r8, & - 60.070000_r8, 158.900000_r8, 121.060000_r8, 55.850000_r8, 94.11_r8, 186.280000_r8, & - 92.110000_r8, 135.080000_r8, 119.080000_r8, 42.090000_r8, 137.110000_r8, & - 88.070000_r8, 119.100000_r8, 90.140000_r8, 76.110000_r8, 76.110000_r8, & - 58.090000_r8, 118.150000_r8, 118.150000_r8, 118.150000_r8, 118.150000_r8, & - 90.090000_r8, 31.400000_r8, 31.400000_r8, 35.450000_r8, 31.400000_r8, & - 31.400000_r8, 35.450000_r8, 64.040000_r8, 96.060000_r8, 31.400000_r8, & - 58.040000_r8, 118.150000_r8, 92.150000_r8, 150.000000_r8, 150.000000_r8, & - 150.000000_r8, 150.000000_r8, 150.000000_r8, 150.000000_r8, 150.000000_r8, & - 150.000000_r8, 106.180000_r8, 12.011000_r8, 12.011000_r8, 135.064039_r8, & - 135.064039_r8, 135.064039_r8, 58.442468_r8, 58.442468_r8, 58.442468_r8, & - 1.007400_r8, 1.007400_r8, 1.007400_r8, 1.007400_r8, 12.011000_r8, & - 12.011000_r8, 115.107340_r8, 115.107340_r8, 115.107340_r8, 250.445000_r8, & - 250.445000_r8, 250.445000_r8, 250.445000_r8, 250.445000_r8, 250.445000_r8, & - 250.445000_r8, 250.445000_r8, 250.445000_r8, 250.445000_r8, 98.078400_r8, & - 250.445000_r8, 250.445000_r8, 250.445000_r8, 250.445000_r8, 250.445000_r8, & - 44.010000_r8, 33.0100000_r8, 17.0100000_r8 /) + adv_mass(:357) = (/ & + 58.09_r8, 60.06_r8, 126.9_r8, & + 44.06_r8, 58.12_r8, 189.12_r8, & + 68.08_r8, 98.1_r8, 150.0_r8, & + 150.0_r8, 150.0_r8, 150.0_r8, & + 150.0_r8, 150.0_r8, 150.0_r8, & + 90.09_r8, 106.12_r8, 12.01_r8, & + 12.01_r8, 78.12_r8, 110.11_r8, & + 79.9_r8, 159.8_r8, 115.45_r8, & + 125.91_r8, 141.91_r8, 95.9_r8, & + 79.9_r8, 79.9_r8, 84.07_r8, & + 138.12_r8, 183.12_r8, 26.05_r8, & + 28.05_r8, 30.08_r8, 44.11_r8, & + 153.82_r8, 137.37_r8, 187.38_r8, & + 170.92_r8, 154.47_r8, 120.91_r8, & + 173.83_r8, 84.93_r8, 267.84_r8, & + 220.84_r8, 176.38_r8, 30.03_r8, & + 94.94_r8, 133.35_r8, 50.45_r8, & + 141.94_r8, 16.04_r8, 252.73_r8, & + 119.35_r8, 35.45_r8, 70.9_r8, & + 102.91_r8, 81.45_r8, 97.45_r8, & + 51.45_r8, 67.45_r8, 1.0_r8, & + 28.01_r8, 108.14_r8, 62.13_r8, & + 29.0_r8, 29.0_r8, 29.0_r8, & + 29.0_r8, 46.07_r8, 105.06_r8, & + 107.07_r8, 78.07_r8, 91.08_r8, & + 62.08_r8, 68.07_r8, 60.06_r8, & + 58.04_r8, 165.36_r8, 148.91_r8, & + 259.82_r8, 18.02_r8, 34.02_r8, & + 74.08_r8, 80.91_r8, 100.13_r8, & + 152.93_r8, 116.94_r8, 100.5_r8, & + 86.47_r8, 36.45_r8, 46.03_r8, & + 127.91_r8, 64.05_r8, 102.1_r8, & + 111.1_r8, 47.01_r8, 63.01_r8, & + 79.01_r8, 96.91_r8, 52.45_r8, & + 143.89_r8, 215.0_r8, 116.13_r8, & + 116.13_r8, 116.13_r8, 116.13_r8, & + 76.06_r8, 126.9_r8, 253.8_r8, & + 285.8_r8, 301.8_r8, 317.8_r8, & + 206.9_r8, 116.13_r8, 162.45_r8, & + 145.13_r8, 150.15_r8, 98.11_r8, & + 148.13_r8, 168.17_r8, 150.15_r8, & + 192.15_r8, 106.14_r8, 106.14_r8, & + 106.14_r8, 147.15_r8, 147.15_r8, & + 147.15_r8, 147.15_r8, 102.0_r8, & + 156.91_r8, 163.15_r8, 163.15_r8, & + 142.9_r8, 14.01_r8, 172.91_r8, & + 188.91_r8, 105.11_r8, 126.9_r8, & + 126.9_r8, 68.13_r8, 195.15_r8, & + 197.17_r8, 136.26_r8, 154.19_r8, & + 154.19_r8, 70.1_r8, 102.1_r8, & + 76.06_r8, 104.12_r8, 86.1_r8, & + 149.11_r8, 149.11_r8, 120.12_r8, & + 124.0_r8, 72.11_r8, 77.05_r8, & + 72.07_r8, 32.05_r8, 14.01_r8, & + 215.28_r8, 215.28_r8, 48.05_r8, & + 147.1_r8, 93.05_r8, 96.1_r8, & + 136.26_r8, 136.26_r8, 70.09_r8, & + 105.13_r8, 102.1_r8, 102.1_r8, & + 120.12_r8, 149.12_r8, 118.1_r8, & + 44.02_r8, 108.02_r8, 17.04_r8, & + 18.05_r8, 62.01_r8, 31.4_r8, & + 30.01_r8, 46.01_r8, 62.01_r8, & + 139.11_r8, 105.11_r8, 48.0_r8, & + 67.45_r8, 12.01_r8, 12.01_r8, & + 60.07_r8, 158.9_r8, 121.06_r8, & + 55.85_r8, 94.11_r8, 186.28_r8, & + 92.11_r8, 135.08_r8, 119.08_r8, & + 42.09_r8, 137.11_r8, 88.07_r8, & + 119.1_r8, 90.14_r8, 76.11_r8, & + 76.11_r8, 58.09_r8, 118.15_r8, & + 118.15_r8, 118.15_r8, 118.15_r8, & + 90.09_r8, 31.4_r8, 31.4_r8, & + 35.45_r8, 31.4_r8, 31.4_r8, & + 35.45_r8, 64.04_r8, 96.06_r8, & + 31.4_r8, 58.04_r8, 118.15_r8, & + 92.15_r8, 150.0_r8, 150.0_r8, & + 150.0_r8, 150.0_r8, 150.0_r8, & + 150.0_r8, 150.0_r8, 150.0_r8, & + 106.18_r8, 12.011_r8, 12.011_r8, & + 135.064039_r8, 135.064039_r8, 135.064039_r8, & + 58.442468_r8, 58.442468_r8, 58.442468_r8, & + 1.0074_r8, 1.0074_r8, 1.0074_r8, & + 1.0074_r8, 12.011_r8, 12.011_r8, & + 115.10734_r8, 115.10734_r8, 115.10734_r8, & + 250.445_r8, 250.445_r8, 250.445_r8, & + 250.445_r8, 250.445_r8, 250.445_r8, & + 250.445_r8, 250.445_r8, 250.445_r8, & + 250.445_r8, 98.0784_r8, 250.445_r8, & + 250.445_r8, 250.445_r8, 250.445_r8, & + 250.445_r8, 75.1_r8, 127.0_r8, & + 89.08_r8, 75.1_r8, 93.0_r8, & + 109.0_r8, 159.13_r8, 137.0_r8, & + 103.11_r8, 103.11_r8, 46.03_r8, & + 60.06_r8, 44.01_r8, 61.06_r8, & + 61.07_r8, 77.06_r8, 1.01_r8, & + 2.02_r8, 33.01_r8, 147.12_r8, & + 147.12_r8, 149.14_r8, 194.14_r8, & + 196.16_r8, 196.16_r8, 196.16_r8, & + 241.14_r8, 149.14_r8, 149.14_r8, & + 117.14_r8, 117.14_r8, 212.16_r8, & + 212.16_r8, 167.16_r8, 167.16_r8, & + 167.16_r8, 146.14_r8, 162.14_r8, & + 162.14_r8, 196.16_r8, 196.16_r8, & + 101.09_r8, 159.13_r8, 159.13_r8, & + 16.04_r8, 28.01_r8, 185.27_r8, & + 68.13_r8, 68.13_r8, 159.17_r8, & + 159.17_r8, 48.0_r8, 173.16_r8, & + 173.16_r8, 187.19_r8, 187.19_r8, & + 101.09_r8, 180.1_r8, 75.05_r8, & + 119.11_r8, 47.04_r8, 119.11_r8, & + 14.01_r8, 28.02_r8, 128.18_r8, & + 159.17_r8, 16.0_r8, 16.0_r8, & + 32.0_r8, 17.01_r8, 230.27_r8, & + 230.27_r8, 61.07_r8, 28.01_r8, & + 34.02_r8, 96.06_r8, 185.27_r8, & + 91.1_r8, 48.0_r8, 136.09_r8, & + 96.06_r8, 96.06_r8, 150.13_r8, & + 89.13_r8, 89.08_r8, 74.09_r8, & + 60.11_r8, 173.16_r8, 187.19_r8 & + /) + ! extfrc_lst contains a list of species that have 3-D emissions (distributed in the vertical by HEMCO) + ! frc_from_dataset (set to false) and extcnt in chem_mods.F90 also needs to be updated. (hplin, 6/24/24) extfrc_lst(: 34) = (/ 'NO ', 'CO ', 'SO2 ', 'SO4 ', & 'NH3 ', 'ACET ', 'ALD2 ', 'ALK4 ', & 'C2H6 ', 'C3H8 ', 'CH2O ', 'PRPE ', & diff --git a/src/chemistry/hetp b/src/chemistry/hetp new file mode 160000 index 0000000000..2a99b24625 --- /dev/null +++ b/src/chemistry/hetp @@ -0,0 +1 @@ +Subproject commit 2a99b24625ed26cf87ae88697ddd6cf8bbdec812 diff --git a/src/chemistry/mozart/mo_neu_wetdep.F90 b/src/chemistry/mozart/mo_neu_wetdep.F90 index 78b3779fe4..f6d94a6bcc 100644 --- a/src/chemistry/mozart/mo_neu_wetdep.F90 +++ b/src/chemistry/mozart/mo_neu_wetdep.F90 @@ -27,6 +27,7 @@ module mo_neu_wetdep real(r8),allocatable, dimension(:) :: mol_weight logical ,allocatable, dimension(:) :: ice_uptake integer :: index_cldice,index_cldliq,nh3_ndx,co2_ndx,so2_ndx + integer :: so4_ndx,so4s_ndx ! geos-chem logical :: debug = .false. integer :: hno3_ndx = 0 ! @@ -142,6 +143,12 @@ subroutine neu_wetdep_init if ( trim(test_name) == 'SO2' ) then so2_ndx = m end if + if ( trim(test_name) == 'SO4' ) then ! GEOS-Chem bulk sulfate + so4_ndx = m + end if + if ( trim(test_name) == 'SO4S' ) then ! GEOS-Chem bulk sulfate on surface seasalt + so4s_ndx = m + end if ! end do @@ -364,20 +371,20 @@ subroutine neu_wetdep_tend(lchnk,ncol,mmr,pmid,pdel,zint,tfld,delt, & endif ! if( dheff(5,l) /= 0._r8 ) then - if( nh3_ndx > 0 .or. co2_ndx > 0 .or. so2_ndx > 0 ) then + if( nh3_ndx > 0 .or. co2_ndx > 0 .or. so2_ndx > 0 .or. so4_ndx > 0 .or. so4s_ndx > 0 ) then e298 = dheff(3,l) dhr = dheff(4,l) dk1s(:) = e298*exp( dhr*wrk(:) ) e298 = dheff(5,l) dhr = dheff(6,l) dk2s(:) = e298*exp( dhr*wrk(:) ) - if( m == co2_ndx .or. m == so2_ndx ) then + if( m == co2_ndx .or. m == so2_ndx .or. m == so4_ndx .or. m == so4s_ndx ) then heff(:,k,m) = heff(:,k,m)*(1._r8 + dk1s(:)*ph_inv*(1._r8 + dk2s(:)*ph_inv)) else if( m == nh3_ndx ) then heff(:,k,m) = heff(:,k,m)*(1._r8 + dk1s(:)*ph/dk2s(:)) else - write(iulog,*) 'error in assigning henrys law coefficients' - write(iulog,*) 'species ',m + if ( masterproc ) write(iulog,*) 'error in assigning henrys law coefficients' + if ( masterproc ) write(iulog,*) 'species ',m end if end if end if diff --git a/src/hemco b/src/hemco index 3a6d999ab0..7bd8358229 160000 --- a/src/hemco +++ b/src/hemco @@ -1 +1 @@ -Subproject commit 3a6d999ab0dbee9f03ab6b9a13dd3b6d9670eb54 +Subproject commit 7bd8358229eefd2cbb910ff30c46dfc97e34fb6f diff --git a/test/system/TR8.sh b/test/system/TR8.sh index 498c9bb57f..cbdb400463 100755 --- a/test/system/TR8.sh +++ b/test/system/TR8.sh @@ -58,12 +58,12 @@ fi #Check Chemistry if [ -d "${CAM_ROOT}/components/cam" ]; then -ruby $ADDREALKIND_EXE -r r8 -l 1 -d $CAM_ROOT/components/cam/src/chemistry -s geoschem +ruby $ADDREALKIND_EXE -r r8 -l 1 -d $CAM_ROOT/components/cam/src/chemistry -s geoschem/geoschem_src,cloud_j,hetp rc=`expr $? + $rc` else -ruby $ADDREALKIND_EXE -r r8 -l 1 -d $CAM_ROOT/src/chemistry -s geoschem +ruby $ADDREALKIND_EXE -r r8 -l 1 -d $CAM_ROOT/src/chemistry -s geoschem/geoschem_src,cloud_j,hetp rc=`expr $? + $rc` fi