ISORROPIA II

From Geos-chem
Revision as of 17:33, 14 January 2020 by Melissa Payer (Talk | contribs) (Investigating persistent noise observed in ISORROPIA output)

Jump to: navigation, search

Overview

Brief description

The ISORROPIA II package performs aerosol thermodynamical equilibrium. It partitions nitrate (HNO3 and NIT) and ammonia (NH3 and NH4) between the gas and aerosol phases. Inputs to the partitioning routine include temperature and RH. ISORROPIA II has significant benefits over previous implementations of ISORROPIA, especially for partitioning of nitrate at low RH.

From Pye et al [2009]:

ISORROPIA II [Fountoukis and Nenes, 2007] is implemented in GEOS-Chem to compute gas-aerosol equilibrium partitioning of nitric acid and ammonia. Particles in this study are not size-resolved; however, they can be generally assumed to represent PM2.5 since formation of sulfate-nitrate-ammonium on coarse mode sea salt and dust is excluded. Submicrometer-sized particles are likely to reach gas-aerosol equilibrium on time-scales less than the 1 hour computational time step used here [Meng and Seinfeld, 1996].
Sodium and chloride from accumulation mode sea salt are considered in the gas-aerosol equilibrium along with sulfate, nitrate, and ammonium. Calcium, magnesium, and potassium concentrations are not considered in the present study because of the issues with dust emissions previously mentioned. All inorganic aerosols are assumed to exist on the upper, metastable branch of the hygroscopic hysterisis curve. Although this assumption may not hold at higher altitudes in the free troposphere [Wang et al., 2008], since the focus of this study is mainly on surface-level concentrations, where humidities reach high values on a daily basis, the metastable assumption is acceptable.

Authors and collaborators

  • Thanos Nenes (Georgia Tech) -- Principal Investigator
  • Havala O. T. Pye (Caltech)

--Bob Y. 12:48, 2 March 2010 (EST)

Implementation notes

ISORROPIA II has been implemented into GEOS-Chem v8-03-01 and later versions.

Code structure

ISORROPIA II consists of the following files:

Files in ISOROPIA/ subdirectory (prior to GEOS-Chem v11-02 (aka 12.0.0)
-----------------------------------------------------------------------
Makefile                 Makefile for ISORROPIA II code
isoropiaIIcode.F         Source code file with ISORROPIA II routines
isrpia.inc               ISOROPIA II header file with common blocks

Files in ISOROPIA/ subdirectory (in GEOS-Chem v11-02 (aka 12.0.0) and later versions)
-------------------------------------------------------------------------------------
Makefile                 Makefile for ISORROPIA II code
isoropiaII_main_mod.F90  Source code file with ISORROPIA II routines

Files in GeosCore/ subdirectory:
--------------------------------
isoropiaII_mod.f         "Interface" code between GEOS-Chem and ISORROPIA II

--Bob Yantosca (talk) 20:14, 6 June 2018 (UTC)

Additional Documentation

Havala O. T. Pye wrote:

  1. Documentation (including a user manual) for ISORROPIAII can be found on the ISORROPIA website.
  2. isoropiaIIcode.f is essentially ISOFWD.FOR and ISOCOM.FOR of the ISORROPIA II box model pasted together.
  3. For more information on ISORROPIA II, see
    • Fountoukis, C. and Nenes, A.: ISORROPIA II: a computationally efficient thermodynamic equilibrium model for K+–Ca2+–Mg2+–NH4+–Na+–SO42−–NO3−–Cl−–H2O aerosols, Atmos. Chem. Phys., 7, 4639-4659, 2007. pdf
  4. The implementation by Pye et al. 2009 JGR did not include Ca, K, or Mg since dust emissions were not used.

--Bob Y. 12:28, 17 March 2010 (EDT)

Validation

See Pye et al [2009]

References

  1. Fountoukis, C., and A. Nenes (2007), ISORROPIA II: A computationally efficient thermodynamic equilibrium model for K+-Ca2+-Mg2+-NH4+-Na+-SO42-NO3--Cl-H2O aerosols, Atmos. Chem. Phys., 7(17), 4639-4659.
  2. Meng, Z. Y., and J. H. Seinfeld, Time scales to achieve atmospheric gas-aerosol equilibrium for volatile species, Atmos. Environ., 30(16), 28892900, doi:10.1016/1352-2310(95)00493-9., 1996.
  3. Pye, H. O. T., H. Liao, S. Wu, L. J. Mickley, D. J.Jacob, D. K. Henze, and J. H. Seinfeld, Effect of changes in climate and emissions on future sulfate-nitrate-ammonium aerosol levels in the United States, J. Geophys. Res., 114, D01205, 2009. PDF
  4. Wang, J., A. A. Hoffman, R. J. Park, D. Jacob, and S. T. Martin, Global distribution of solid and aqueous sulfate aerosols: Effect of the hysteresis of particle phase transitions, J. Geophys. Res., 113, D11206, doi:10.1029/2007JD009367, 2008. PDF

--Bob Y. 16:00, 22 February 2010 (EST)

Previous issues now resolved

Bug fixes for ISORROPIA II stable mode

These two fixes were included in v11-02e (approved 24 Mar 2018).

Jingyuan Shao (UW) wrote:

In isoropiaII_mod.F (line 558) when calculate aerosol PH by isorropia, the value would be too high when we choose stable mode. After I checked the code, I think code should be changed as follows :
        IF ( AERLIQ(8) < 1e-32_fp ) THEN
           ! Aerosol is dry so HPLUSTEMP and PH_SAV are undefined
           ! We force HPLUSTEMP to 1d20 and PH_SAV to -999e+0_fp.
           ! (hotp, ccc, 12/18/09)
           HPLUSTEMP       = 1e+20_fp
           PH_SAV(I,J,L)   = -999e+0_fp
        ELSE
           HPLUSTEMP       = AERLIQ(1) / AERLIQ(8) * 1e+3_fp/18e+0_fp
to
        IF ( AERLIQ(8) < 1e-18_fp ) THEN
           ! Aerosol is dry so HPLUSTEMP and PH_SAV are undefined
           ! We force HPLUSTEMP to 1d20 and PH_SAV to -999e+0_fp.
           ! (hotp, ccc, 12/18/09)
           HPLUSTEMP       = 1e+20_fp
           PH_SAV(I,J,L)   = -999e+0_fp
        ELSE
           HPLUSTEMP       = AERLIQ(1) / AERLIQ(8) * 1e+3_fp/18e+0_fp
or
        IF ( OTHER(1) .eq. 1.d0 ) THEN
           ! Aerosol is dry so HPLUSTEMP and PH_SAV are undefined
           ! We force HPLUSTEMP to 1d20 and PH_SAV to -999e+0_fp.
           ! (hotp, ccc, 12/18/09)
           HPLUSTEMP       = 1e+20_fp
           PH_SAV(I,J,L)   = -999e+0_fp
        ELSE
           HPLUSTEMP       = AERLIQ(1) / AERLIQ(8) * 1e+3_fp/18e+0_fp
the default is metastable mode, there is no influence if we choose metastable mode because metastable means supersaturation state, water always exist on aerosol if we choose metastable mode.

Shaojie Song (Harvard) wrote:

We have identified and fixed coding errors in ISORROPIA which significantly affect aerosol water pH calculations when the stable state for aerosol phase is applied. These errors exist in source codes for the cases G2 and O2, and also affect the pH results for the cases G1 and O1. The standard ISORROPIA codes fail to account for the partitioning of ammonia between the gas and aqueous phases. More information can be seen in our recent ACPD paper: https://www.atmos-chem-phys-discuss.net/acp-2018-6/.
In order to resolve this issue, in isoropiaIIcode.F of GEOS-Chem v11-01, both line 20401 and line 25829 (which consist of the same codes):
     PSI4 = MIN(PSI5+PSI6,CHI4)

should be replaced with

     IF(W(2).GT.TINY) THEN       ! Accounts for NH3 evaporation
        BB   =-(CHI4 + PSI6 + PSI5 + 1.d0/A4)
        CC   = CHI4*(PSI5+PSI6) - 2.d0*PSI2/A4
        DD   = MAX(BB*BB-4.d0*CC,ZERO)  ! Patch proposed by Uma Shankar, 19/11/01
        PSI4 =0.5d0*(-BB - SQRT(DD))
     ELSE
        PSI4 = TINY
     ENDIF
The same coding errors probably exist in all of the previous versions of GEOS-Chem because they exist in the versions from 1.5 to 2.1 of ISORROPIA. Below is an example showing the average aerosol pH by GEOS-Chem v10-01 during Jan 1st-Jan 15th, 2015 (left: standard model; right: revised model). A big difference can be found over some land areas such as USA, Europe, and East Asia, which is mainly because their relatively higher concentrations of ammonia. The marked value represents the pH of Beijing, which is reduced from about 7 (neutral) to about 4 (moderately acidic). Note that there is no influence when using metastable mode (the default).
Map ph global.png

--Melissa Sulprizio (talk) 18:48, 30 January 2018 (UTC)

Implement ISORROPIA v2.0 as a Fortran module

This update was included in v11-02a and approved on 12 May 2017.

Seb Eastham has repackaged the existing ISORROPIA v2.0 code as a Fortran module. This finally allows us to remove files ISOROPIA/isoropiaII_code.F and ISOROPIA/isrpia.inc, which have proven to be problematic for some newer compilers (e.g. GNU Fortran v5 and v6). He writes:

  1. I’ve adapted ISORROPIA v2.0 into a Fortran module, converting COMMON blocks and BLOCK DATA statements into module-level arrays. This can safely replace the current implementation in GEOS-Chem.
  2. All REAL variables are declared using a custom KIND parameter. Implicit variables are forbidden (using IMPLICIT NONE at the module level) and all variables are declared explicitly. This allows precision to be specified easily.
  3. The number of OpenMP THREADPRIVATE variables has also been significantly reduced to minimize overhead and improve code readability.
  4. The module can be switched between the legacy version (v2.0) and the updated version (v2.2) by (un)commenting the line #define ISORROPIA_V22 at the top of the module file.
  5. When compiled as ISORROPIA v2.0, the file exactly reproduces results produced with GEOS-Chem (zero differences w/r/t prior model versions).
  6. I’ve put together a wrapper program which can take the isorropiaII_main_mod.F file exactly as written, compile it into a standalone program, run a parameter sweep and store the output into a NetCDF file. Contact me if this might be of use to you.

Update 21 Feb 2017: The ISORROPIA module now allows you to compile GEOS-Chem with GNU Fortran v4.9 and higher.

Update 26 May 2017: We have brought the ISORROPIA module update into the GEOS-Chem v11-01 public release code, in the v11-01-ISORROPIA-Gfortran branch. If you would like to use Gfortran to compile the GEOS-Chem v11-01 public-release code, we recommend that you check out the code as follows:

 git clone https://bitbucket.org/gcst/geos-chem Code.v11-01
 cd Code.v11-01
 git checkout v11-01-ISORROPIA-Gfortran

--Bob Yantosca (talk) 18:30, 26 May 2017 (UTC)

Optimization and/or parallelization issues in ISORROPIA II

During the benchmarking process for GEOS-Chem v9-01-02, we discovered that ISORROPIA does not produce exactly the same output on single processor compared to multi-processor:

Daven Henze wrote:

There is an unidentified OpenMP error which contributes to differences shown [in the v9-01-02a 1-month benchmark], as confirmed by independent tests of consecutive evaluations of a single executable compiled with and without parallelization. The impacts of this OpenMP error are typically less than a fraction of a percent, and negligible in terms of absolute value, with the exception of NIT and NH3. This is an issue to investigate further.

Bob Yantosca wrote:

Here is a quick update on the issue in ISORROPIA.
  1. If you compile with –O0 (all optimization turned off), you get identical results when using ISORROPIA as single vs. multi processor.
  2. When you compile with –O1, then you get different results when using ISORROPIA as single vs. multi processor. This inclines me to believe that the ISORROPIA code is sensitive to optimization.
  3. It appears that issue is occurring at the computation of the GNH3, GHNO3, and GHCL variables in function FUNCG5A within file isoropiaIIcode.f. These variables are computed as:
     GNH3      = MAX(CHI4 - PSI4, TINY)              ! Gas NH3
     GHNO3     = MAX(CHI5 - PSI5, TINY)              ! Gas HNO3
     GHCL      = MAX(CHI6 - PSI6, TINY)              ! Gas HCl
Here the CHI* and PSI* are common-block variables and TINY is a small number (1d-20, but I used 1d-28 for testing). The CHI* and PSI* values are held threadprivate for the OpenMP parallelization. It looks as if there is some kind of roundoff error at the limit of precision that is creeping in here that could be caused by the optimization. For example, at grid box (37,21,11) every once and a while we get this roundoff error:
   with single processor
   ### CHI5 :   4.326986976485327E-008
   ### PSI5 :   4.312588684524357E-008
   ### GHNO3:   1.439829196097068E-010

   with multi processor
   ### CHI5 :   4.326986976485327E-008
   ### PSI5 :   4.312588684524356E-008
   ### GHNO3:   1.439829196097134E-010
with similar roundoff errors (not shown here) happening in the computation of GNH3.
GHNO3 gets saved back into the tracer array STT(I,J,L,IDTHNO3) and GNH3 gets saved back into STT(I,J,L,IDTNH3). So that is why these we see these very small (and seemingly random) differences in HNO3, NH3, and NIT, whereas the other tracers are unaffected.
I don’t have a good sense of why precisely this is occurring, other than to say that ISORROPIA uses a lot of internal common blocks and other obsolete coding practices. The ISORROPIA code is extremely confusing to follow. My guess is that either something is not being held private for the parallelization, or that variables in common blocks are not being initialized properly.

Solution

Shannon Capps wrote:

Please see this attached report that documents my replication of the issue with both the GEOS-Chem and ANISORROPIA versions of ISORROPIA II in a box model....I would propose adding a flag (-fp-model source) for the compilation of ISORROPIA only. Add this to the file ISOROPIA/Makefile as follows:
   isoropiaIIcode.o:  isoropiaIIcode.F  isrpia.inc
           $(F90) $(R8) -fp-model source -c $<
Compiling the GEOS-Chem implementation of ISORROPIA with this flag causes the box model to produce results that "diff" finds identical in all combinations (e.g., m0 v. s0, m3 v. s0, etc.).

The -fp-source option of the Intel Fortran Compiler tells the compiler to only do "safe" optimizations (i.e. nothing that would affect the numerical precision of the output). For some reason, ISORROPIA II is very sensitive to the compiler optimization.

Bob Yantosca has done some additional 3-day test runs to assess the impact of using the -fp-model source option. The results are:

Run # Machine # CPUs Wall time [mm:ss] Scalability Mean OH [1e5 molec/cm3]
1 terra-03-c.as.harvard.edu 4 18:10 3.5935 13.0384629976933
2 terra-04-c.as.harvard.edu 4 09:58 3.8123 13.0384629976933
3 terra-03-c.as.harvard.edu 4 16:38 3.7546 13.0384629976933
4 terra-04-c.as.harvard.edu 4 11:00 3.5090 13.0384629976933
5 terra-11-c.as.harvard.edu 8 09:29 7.3364 13.0384629976933
6 terra-06-c.as.harvard.edu 8 06:18 7.1421 13.0384629976933
7 terra-11-c.as.harvard.edu 8 09:10 7.5530 13.0384629976933
8 terra-01-c.as.harvard.edu 8 05:52 7.4549 13.0384629976933

All of the output files from runs 1-8 were 100% binary identical to each other. Therefore, using -fp-model source seems to eliminate the random numerical noise generated by ISORROPIA.

Persistence of problem

Update 8/26/11: In GEOS-Chem v9-01-02 we shall apply the -fp-model source option to all source code files, not just ISORROPIA. This will help us to ensure the precision of the results as we make widespread changes for ESMF compatibility.

Update 9/1/11: Oops, this may not be solved after all. Runs longer than 3 days still manifest the problem. See this PDF presentation for more details.

--Bob Y. 15:32, 23 November 2011 (EST)

Minor fixes in isorropiaIIcode.f

Havala Pye wrote

There are 2 minor bug fixes for ISORROPIA II. There are two locations in file ISOROPIA/isorropiaIIcode.f where the variable MOLALR(4) should be MOLALR(9). These fixes affect how aerosol phase water is calculated. Since GEOS-Chem typically considers sodium and chloride, these subcases are likely not being used, but the code should be updated in case these subcases are used in the future.

This fix will be implemented in GEOS-Chem v9-01-02.

--Bob Y. 11:10, 1 March 2011 (EST)

Patches for v8-03-01

The following issues in ISORROPIA II have now been fixed as of August 2010:

  1. Bug fix in ISORROPIA for offline aerosol
  2. Patch to fix declaration of IONIC in ISORROPIA II
  3. Bug fix for ND42 diagnostic and ISORROPIA II

--Bob Y. 14:44, 3 August 2010 (EDT)

Outstanding issues

Parallelization issue for nested grid simulations

For those of you who are running one of the GEOS-Chem nested grid simulations, you may encounter a segmentation fault error. This error is caused by GEOS-Chem running up against the memory limits of the machine. One way to solve this problem is to compile ISORROPIA II for single processor only:

Sungshik Patrick Kim wrote:

I tracked down the memory bottlenecks in my [China nested-grid simulation]....The crashes then were caused by memory allocation in ISORROPIA (with all of the private variables that had to be copied), so after turning off the parallelization for the ISORROPIA driver routine, I can actually spin up the model!

To compile ISORROPIA II for single processor, while compiling everything else for multi-processor, you can add an extra comment character ! to the !$OMP commands in GeosCore/isoropiaII_mod.f, such as:

!!$OMP PARALLEL DO
!!$OMP+DEFAULT( SHARED )
!!$OMP+PRIVATE( I,    J,      L,       N,      WI,   WT,  GAS,  TEMPI )
!!$OMP+PRIVATE( RHI,  VOL,    TSO4,    TNH3,   TNA,  TCL, ANO3, GNO3  )
!!$OMP+PRIVATE( TCA,  TMG,    TK,      CNTRL,  SCASI                  )
!!$OMP+PRIVATE( TNO3, AERLIQ, AERSLD,  OTHER,  TNH4, TNIT             )
!!$OMP+PRIVATE( HPLUSTEMP,    NUM_SAV, DEN_SAV, HNO3_DEN              )
!!$OMP+SCHEDULE( DYNAMIC )

...

!!$OMP END PARALLEL DO

--Bob Y. 10:23, 8 April 2011 (EDT)

Aerosol pH

This issue is somewhat resolved.

Havala O. T. Pye wrote:

At the GEOS-Chem User's meeting (2009) Becky Alexander and her student Eric noted that ISORROPIA II sometimes returns a negative H+ concentration. In my tests, the negative values occurred in N. India and Western Russia area during limited times of the year. A quick fix has been implemented. A condition I tracked down was returning a negative [H+] when solving the HSO4 = H + SO4 equilibrium (CALCHS4 in isorropiaIIcode.f) when the actual answer should have been ~1e-27 mol/m3. (I plugged the same input values into a spreadsheet and it returned the same negative answer so it seems like a numerical precision issue.) For my case study, it came down to subtracting two numbers of very similar value (first 8 digits the same) which resulted in a negative. I put in a fix to reset H+ to 1d-30 in CALCHS4 if it goes negative.

--Bob Y. 10:29, 29 January 2010 (EST)