-
Notifications
You must be signed in to change notification settings - Fork 170
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Added a new inorganic nucleation subroutine to tomas_mod.F90. #2528
base: main
Are you sure you want to change the base?
Conversation
The nucleation scheme is based on this paper: https://www.science.org/doi/10.1126/science.aaf2649. We added a 1000x multiplier to the total inorganic nucleation rate in order to match observations better.
Also tagging @BettyCroft who has worked with TOMAS more recently |
@@ -140,7 +140,8 @@ MODULE TOMAS_MOD | |||
REAL(fp), SAVE, ALLOCATABLE :: ECSCALE30(:) | |||
REAL(fp), SAVE, ALLOCATABLE :: ECSCALE100(:) | |||
|
|||
INTEGER :: bin_nuc = 1, tern_nuc = 1 ! Switches for nucleation type. | |||
INTEGER :: bin_nuc = 0, tern_nuc = 0 ! Switches for nucleation type. | |||
INTEGER :: dunn_nuc = 1 |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
If these are intended to be constant values, consider declaring these with the PARAMETER attribute (i.e. INTEGER, PARAMETER :: bin_nuc =0
, etc. That will provide some extra protection in case someone else tries to modify the value elsewhere.
@@ -1979,8 +1994,13 @@ SUBROUTINE getNucRate(Nk, Mk, Gci,fn,mnuc,nflg, ionrate,surf_area, & | |||
h2so4 = Gci(srtso4)/boxvol*1000.e+0_fp/98.e+0_fp*6.022e+23_fp | |||
nh3ppt = Gci(srtnh4)/17.e+0_fp/(boxmass/29.e+0_fp)*1e+12_fp* & | |||
PRES/101325.*273./TEMPTMS ! corrected for pressure (because this should be concentration) | |||
|
|||
nh3moleccm3 = Gci(srtnh4)/boxvol*1000.e+0_fp/17e+0_fp*6.022e+23_fp ! Changed by SamO |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
For consistency, please use the MW of 17.04 for NH3 as defined in the species_database.yml
file here:
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Also please use the AIRMW
parameter from Headers/physconstants.F90
for the MW of air. This is defined as 28.9644 g/mol.
https://github.com/samuelod-atmos/geos-chem/blob/main/Headers/physconstants.F90
|
||
Mair = 2.69E19*273.15/TEMPTMS*PRES/101325. |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Use _fp
to prevent loss of precision:
Mair = 2.69E19_fp * 273.15_fp / TEMPTMS * PRES / 101325.0_fp
@@ -2787,12 +2822,16 @@ SUBROUTINE NUCLEATION(Nki,Mki,Gci,Nkf,Mkf,Gcf,fn,fn1,totsulf, & | |||
|
|||
errorswitch = .false. | |||
|
|||
Mair = 2.69E19*273.15/TEMPTMS*PRES/101325. |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Use _fp
to prevent loss of precision:
Mair = 2.69E19_fp * 273.15_fp / TEMPTMS * PRES / 101325.0_fp
h2so4 = Gci(srtso4)/boxvol*1000.e+0_fp/98.e+0_fp*6.022e+23_fp | ||
nh3ppt = Gci(srtnh4)/17.e+0_fp/(boxmass/29.e+0_fp)*1e+12_fp* & | ||
PRES/101325.*273./TEMPTMS ! corrected for pressure (because this should be concentration) | ||
nh3moleccm3 = Gci(srtnh4)/boxvol*1000.e+0_fp/17e+0_fp*6.022e+23_fp ! Changed by SamO |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Use 17.04 for the MW of NH3 (see above comment)
parameter(pbn=3.95, ubn=9.70, vbn=12.6, wbn=-0.00707, ptn=2.89) | ||
parameter(utn=182., vtn=1.20, wtn=-4.19, pAn=8.00, an=1.6d-6) | ||
parameter(pbi=3.37, ubi=-11.5, vbi=25.5, wbi=0.181, pti=3.14) | ||
parameter(uti=-23.8, vti=37.0, wti=0.227, pAi=3.07, ai=0.00485) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Each of the parameters should be followed by d0
to prevent a loss of precision
temp=tempi !SamO | ||
!temp=278.0 | ||
nh3=nh3i*1E-6 ! I think they want it in units of 1E6 molec cm-3 | ||
cna=cnai*1E-6 |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Use 1d-6
to prevent loss of precision
!fion = 75.0 !SamO | ||
|
||
! CALCULATE ION CONCENTRATION | ||
alpha_ion = 6d-8*sqrt(300./temp) + 6d-26*Mair*(300./temp)**4 ! Need Mair in air molec per cm3 |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Use 300.0d0
to prevent a loss of precision.
kbn = exp(ubn - exp(vbn*(temp/1000. - wbn))) | ||
ktn = exp(utn - exp(vtn*(temp/1000. - wtn))) | ||
kbi = exp(ubi - exp(vbi*(temp/1000. - wbi))) | ||
kti = exp(uti - exp(vti*(temp/1000. - wti))) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Use 1000.0d0 to prevent a loss of precision.
ffi = nh3*cna**pti/(ai + (cna**pti)/(nh3**pAi)) | ||
else | ||
ffn = 0. | ||
ffi = 0. |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Use 0.0d0
to prevent a loss of precision
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Thanks for this PR @samuelod-atmos. I've noted several minor issues.
Also would add a note to the CHANGELOG.md file under the Added
section Something
- Added new inorganic nucleation subroutine to `tomas_mod.F90` (cf. 10.1126/science.aaf2649)
The nucleation scheme is based on this paper:
https://www.science.org/doi/10.1126/science.aaf2649.
Name and Institution (Required)
Name: Samuel O'Donnell
Institution: Colorado State University
Describe the update
We have added a new inorganic nucleation mechanism to the TOMAS aerosol package. The nucleation mechanisms are based on Dunne et al. (2016; https://www.science.org/doi/10.1126/science.aaf2649), and was originally written by Jeff Pierce at CSU. We have not included the organic nucleation mechanism outlined in Riccobono et al. (2014) which is part of the Dunne et al. (2016) paper.
A new switch has been added (dunn_nuc; 0/1) in order to turn the mechanism off/on, consistent with the existing inorganic nucleation mechanisms. The updated scheme uses sulfuric acid (molecules/cm^3) and ammonia (molecules/cm^3) in a simple power-law function to give the formation rate of small particles. In order to better match observations of new particle formation events at the Southern Great Plains (SGP) observatory, we had to scale up the nucleation mechanism by 1000x for the 15 bin simulation. This updated is relatively untested with TOMAS40, so the 1000x multiplier should be tested with a TOMAS40 simulation.
Expected changes
Overall, nucleation and particle growth is a buffered system, so the nucleation rates should result in little changes to bulk aerosol mass. Ideally, the updated mechanisms should lead to better agreement with observations in terms of new particle formation event frequency, duration, and strength throughout the globe, but we have only tested this at the SGP site so far.
Reference(s)
The work with updates to TOMAS are not published yet.
Related Github Issue
Please link to the corresponding Github issue(s) here. If fixing a bug, there should be an issue describing it with steps to reproduce.