geoschem / geos-chem

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

Missing tags in ProdLoss output #2276

Open alli-moon opened 6 months ago

alli-moon commented 6 months ago

Name and Institution (Required)

Name: Alli Moon Institution: University of Washington

Confirm you have reviewed the following documentation

Hi there, I am experiencing an issue obtaining the output for heterogeneous and photolysis rxn tags in ProdLoss. The gas-phase reaction tags are working correctly. I have added around 260 tags (Halogens). I turned on the verbose output. I can see when I initialize the run I don't see the whole ProdLoss family for heterogeneous and photolysis reactions like I do the gas phase (see log file)

I have tried increasing the number of allowed families in KPP, so I don't think that is the issue. I am attaching the files I modified in KPP (custom.kpp and custom.eqn), the species database file , the geos_config file, and the HISTORY file. I can attach the Restart and Prodloss output files if that is needed, I just left them out since they are so large.

Thanks! Alli

GC.log custom.eqn.txt custom.kpp.txt GC.log geoschem_config.yml.txt HISTORY.rc.txt

yantosca commented 6 months ago

Thanks for writing @alli-moon. I checked the KPP code, it should allow you to have 300 families, so that isn't the issue. I will see if I can replicate the issue with your files. Stay tuned.

yantosca commented 6 months ago

@alli-moon: I built your mechanism but i get these errors:

  455 |   RCONST(24) = (GC_OHCO_a(1.50d-13))
      |                1
Error: Function ‘gc_ohco_a’ at (1) has no IMPLICIT type; did you mean ‘gc_epo_a’?
/n/home09/ryantosca/T/for_alli/test_std/CodeDir/src/GEOS-Chem/KPP/custom/gckpp_Rates.F90:469:16:

  469 |   RCONST(38) = (GC_OHHNO3_acacac(2.41d-14,460.0d0,2.69d-17,2199.0d0,6.51d-34,1335.0d0))
      |                1
Error: Function ‘gc_ohhno3_acacac’ at (1) has no IMPLICIT type; did you mean ‘gc_ho2ho2_acac’?
/n/home09/ryantosca/T/for_alli/test_std/CodeDir/src/GEOS-Chem/KPP/custom/gckpp_Rates.F90:493:16:

  493 |   RCONST(62) = (GC_TBRANCH_2_acabc(7.60d-12,-585.0d0,5.87d0,0.64d0,-816.0d0))
      |                1
Error: Function ‘gc_tbranch_2_acabc’ at (1) has no IMPLICIT type; did you mean ‘gc_tbranch_1_acac’?
/n/home09/ryantosca/T/for_alli/test_std/CodeDir/src/GEOS-Chem/KPP/custom/gckpp_Rates.F90:494:16:

  494 |   RCONST(63) = (GC_TBRANCH_2_acabc(7.60d-12,-585.0d0,1.7d-1,-0.64d0,816.0d0))
      |                1
Error: Function ‘gc_tbranch_2_acabc’ at (1) has no IMPLICIT type; did you mean ‘gc_tbranch_1_acac’?
/n/home09/ryantosca/T/for_alli/test_std/CodeDir/src/GEOS-Chem/KPP/custom/gckpp_Rates.F90:822:17:

  822 |   RCONST(391) = (GC_OHHNO3_acacac(2.41d-14,460.0d0,2.69d-17,2199.0d0,6.51d-34,1335.0d0))
      |  

Could you attach the rate law functions file (fullchem_RateLawFuncs.F90) that you used. I think these are extra functions that you have that I don't. I'm using 14.3.1.

alli-moon commented 6 months ago

Hi @yantosca, Thanks for looking into it! That is an interesting error. The only thing I can think of is that I added patch 2177 back when I was trying to get Obspack to work. Here is the file in my custom KPP folder! Thanks for your help! Alli

fullchem_RateLawFuncs.F90.txt

yantosca commented 6 months ago

Thanks @alli-moon. I may not get to this until Monday (5/6) at this point. But I will try it!

yantosca commented 6 months ago

Hi @alli-moon. Thanks for the rate law function. I downloaded and compiled it and it built the executable. When I tried to run the code (14.3.1) I get this error:

===============================================================================
GEOS-Chem ERROR: MW_g for species RXN_s_SO2_SALAAL is undefined!
 -> at Init_Species_Database (in module Headers/species_database_mod.F90
===============================================================================

===============================================================================
GEOS-Chem ERROR: Error encountered in routine "Init_Species_Database"!
 -> at Init_State_Chm (in module Headers/state_chm_mod.F90)
===============================================================================

===============================================================================
GEOS-Chem ERROR: Error encountered within call to "Init_State_Chm"!
 -> at GC_Init_StateObj (in GeosCore/gc_environment_mod.F90)
===============================================================================

===============================================================================
GEOS-CHEM ERROR: Error encountered in "GC_Init_StateObj!"!
STOP at  -> at GEOS-Chem (in GeosCore/main.F90)
===============================================================================

Do you have the molecular weight specified in your own species_database.yml file? If so, could you please post that here? I'll try again.

alli-moon commented 6 months ago

Hi @yantosca , Here is the file I'm currently using where the S reactions have a MW_g specified. I just started a run with it this morning. Does the MW_g matter for ProdLoss as long as there is a value there? Thanks! Alli species_database.yml.txt

yantosca commented 6 months ago

Hi @alli-moon. I think it just needs a default value (1.0 would work). Thanks!

yantosca commented 6 months ago

Hi again @alli-moon. Thanks for the species_database.yml file. I found and fixed a few issues in it, see this file: species_database.for_alli.txt.

I think one of the issues may have been that several of the prod/loss variables in the species_database.yml had a space between the end of the variable name and the semicolon. Maybe that caused the code to think that the space was part of the species name.

In any case, with the fixed file attached here, I was able to get all of the PROD_RXN_* tags saved to the GEOSChem.ProdLoss collection file.

$ ncdump -cts GEOSChem.ProdLoss.20190701_0000.nc4 | grep "float PROD_RXN"
128:    float AREA(lat, lon) ;

167:    float Prod_RXN_p_ICl(time, lev, lat, lon) ;
177:    float Prod_RXN_p_IBr(time, lev, lat, lon) ;
. . .
2707:   float Prod_RXN_s_HNO3_SALCAL(time, lev, lat, lon) ;
2717:   float Prod_RXN_s_HCl_SALCAL(time, lev, lat, lon) ;
2727:   float Prod_RXN_s_SO2_SALCAL(time, lev, lat, lon) ;
2737:   float Prod_RXN_s_HNO3_SALAAL(time, lev, lat, lon) ;
2747:   float Prod_RXN_s_HCl_SALAAL(time, lev, lat, lon) ;
2757:   float Prod_RXN_s_SO2_SALAAL(time, lev, lat, lon) ;

and then I counted the number of PROD_RXN variables in the file with this command:

$ ncd GEOSChem.ProdLoss.20190701_0000z.nc4 | g "float PROD_RXN"  | wc -l
260

so it looks like it works. Let me know if you can replicate this on your end.

alli-moon commented 6 months ago

Hi @yantosca, Thank you so much for the new file! I don't think I could have caught the space issue without help. I will add this file and re-run now and report back. Thanks!! Alli

alli-moon commented 6 months ago

Hi @yantosca, I added the new species_database file and also removed some spaces I found behind the tags in the custom.kpp file. I recompiled it in KPP and my run directory. I'm still not seeing the heterogeneous and photolysis production reactions. Are there any other changes I need to add? Thanks! Alli

yantosca commented 6 months ago

Hi @alli-moon. I used the 14.3.1 "out-of-the-box" code and created a fullchem/standard run directory. Then I added your custom.kpp, custom.eqn,KPP/custom/fullchem_RateLawFuncs.F90and the fixedspecies_database.yml. When I compiled and built with-DMECH=custom` everything just worked. Which base version do you use? Also I built with GCC 10.2.0 compilers (although I don't think the compiler version would be the issue).

alli-moon commented 6 months ago

Hi @yantosca, This is with v14.2.3.

We are currently thinking of using version 14.4 since we are excited about the HETP addition. Perhaps this issue will resolve itself once I download the new version (and maybe my obspack bug will also be fixed).

In the meantime, I may try what you just described in 14.3 as a sanity check.

Thank you for your help! Best, Alli

yantosca commented 6 months ago

Thanks @alli-moon. We have 14.4.0 in benchmarking now. Hopefully to be released by IGC11. Stay tuned for release news.

alli-moon commented 4 months ago

Hi @yantosca, I implemented the same tags in version 14.4 and I see the photolysis and heterogeneous output in ProdLoss now! Thank you for your help! Also if anyone else needs to add halogen tags I am happy to share my KPP and species database files with them :) Best, Alli

yantosca commented 4 months ago

HI @alli-moon, nice to have met you at IGC11. Thanks for the update! Also note, we will be releasing a 14.4.1 version soon due to a bug that caused a slowdown (see the workaround in PR #2347).

At IGC11, some folks still reported issues in the ProdLoss so I will have to take a deeper dive (when I get a chance).