Closed yantosca closed 1 year ago
Thanks Bob for working on this! Looking forward to confirming we get the same results.
Results of a test using examples/saprc2006.kpp
with KPP 3.0.0:
! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!
! Log File with Human-Readable Information
!
! Generated by KPP-3.0.0 symbolic chemistry Kinetics PreProcessor
! (https:/github.com/KineticPreProcessor/KPP
! KPP is distributed under GPL, the general public licence
! (http://www.gnu.org/copyleft/gpl.html)
! (C) 1995-1997, V. Damian & A. Sandu, CGRER, Univ. Iowa
! (C) 1997-2022, A. Sandu, Michigan Tech, Virginia Tech
! With important contributions from:
! M. Damian, Villanova University, Philadelphia, PA, USA
! R. Sander, Max-Planck Institute for Chemistry, Mainz, Germany
! M. Long, Renaissance Fiber, LLC, North Carolina, USA
! H. Lin, Harvard University, Cambridge, MA, USA
! R. Yantosca, Harvard University, Cambridge, MA, USA
!
! File : saprc2006.log
! Equation file : saprc2006.kpp
! Output root filename : saprc2006
!
! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
### Options -------------------------------------------
#DECLARE - OFF
#DOUBLE - ON
#DRIVER - general
#DUMMYINDEX - OFF
#EQNTAGS - OFF
#FUNCTION - AGGREGATE
#HESSIAN - ON
#INTEGRATOR - rosenbrock
#JACOBIAN - SPARSE_LU_ROW
#LANGUAGE - Fortran90
#MEX - ON
#MINVERSION - none
#REORDER - ON
#STOCHASTIC - OFF
#STOICMAT - ON
#UPPERCASEF90 - OFF
### Parameters ----------------------------------------
! NSPEC - Number of chemical species
INTEGER, PARAMETER :: NSPEC = 93
! NVAR - Number of Variable species
INTEGER, PARAMETER :: NVAR = 88
! NVARACT - Number of Active species
INTEGER, PARAMETER :: NVARACT = 81
! NFIX - Number of Fixed species
INTEGER, PARAMETER :: NFIX = 5
! NREACT - Number of reactions
INTEGER, PARAMETER :: NREACT = 235
! NVARST - Starting of variables in conc. vect.
INTEGER, PARAMETER :: NVARST = 1
! NFIXST - Starting of fixed in conc. vect.
INTEGER, PARAMETER :: NFIXST = 89
### Species -------------------------------------------
Variable species
1 = H2SO4 (r) 31 = N2O5 (r) 61 = ISOPRENE (r)
2 = BC (r) 32 = HONO (r) 62 = R2O2 (r)
3 = OC (r) 33 = ALK3 (r) 63 = TERP (r)
4 = SSF (r) 34 = TBU_O (r) 64 = METHACRO (r)
5 = SSC (r) 35 = ALK5 (r) 65 = OLE1 (r)
6 = PM10 (r) 36 = ARO2 (r) 66 = ISOPROD (r)
7 = PM25 (r) 37 = HNO4 (r) 67 = OLE2 (r)
8 = DMS (r) 38 = COOH (r) 68 = MVK (r)
9 = DST1 (r) 39 = HOCOO (r) 69 = CCHO (r)
10 = DST2 (r) 40 = BZNO2_O (r) 70 = HCHO (r)
11 = DST3 (r) 41 = MEOH (r) 71 = RNO3 (r)
12 = CO2 (r) 42 = ALK4 (r) 72 = O3P (r)
13 = HCOOH (n) 43 = ARO1 (r) 73 = RCHO (r)
14 = CCO_OH (n) 44 = DCB3 (r) 74 = MEK (r)
15 = RCO_OH (n) 45 = DCB2 (r) 75 = PROD2 (r)
16 = CCO_OOH (n) 46 = CRES (r) 76 = O3 (r)
17 = RCO_OOH (n) 47 = C2H2 (r) 77 = BZCO_O2 (r)
18 = XN (n) 48 = DCB1 (r) 78 = HO2 (r)
19 = XC (n) 49 = BALD (r) 79 = OH (r)
20 = SO2 (r) 50 = ROOH (r) 80 = MA_RCO3 (r)
21 = O1D (r) 51 = NPHE (r) 81 = NO (r)
22 = C2H6 (r) 52 = PHEN (r) 82 = NO3 (r)
23 = BACL (r) 53 = CO (r) 83 = C_O2 (r)
24 = PAN (r) 54 = MGLY (r) 84 = CCO_O2 (r)
25 = PAN2 (r) 55 = ACET (r) 85 = NO2 (r)
26 = PBZN (r) 56 = HNO3 (r) 86 = RO2_N (r)
27 = MA_PAN (r) 57 = ETHENE (r) 87 = RO2_R (r)
28 = H2O2 (r) 58 = C3H6 (r) 88 = RCO_O2 (r)
29 = C3H8 (r) 59 = GLY (r)
30 = ETOH (r) 60 = BZ_O (r)
Fixed species
1 = AIR (r) 3 = H2O (r) 5 = CH4 (r)
2 = O2 (r) 4 = H2 (r)
### Subroutines ---------------------------------------
SUBROUTINE Fun ( V, F, RCT, Vdot, Aout )
SUBROUTINE Fun_SPLIT ( V, F, RCT, Vdot, P_VAR, D_VAR, Aout )
SUBROUTINE CalcStoichNum ( StoichNum )
SUBROUTINE Jac_SP ( V, F, RCT, JVS )
SUBROUTINE Jac_SP_Vec ( JVS, UV, JUV )
SUBROUTINE JacTR_SP_Vec ( JVS, UV, JTUV )
SUBROUTINE KppSolve ( JVS, X )
SUBROUTINE KppSolveTR ( JVS, X, XX )
SUBROUTINE Hessian ( V, F, RCT, HESS )
SUBROUTINE HessTR_Vec ( HESS, U1, U2, HTU )
SUBROUTINE Hess_Vec ( HESS, U1, U2, HU )
SUBROUTINE Initialize ( )
SUBROUTINE Shuffle_user2kpp ( V_USER, V )
SUBROUTINE Shuffle_kpp2user ( V, V_USER )
SUBROUTINE Update_RCONST ( )
SUBROUTINE Update_PHOTO ( )
SUBROUTINE GetMass ( CL, Mass )
SUBROUTINE ReactantProd ( V, F, ARP )
SUBROUTINE JacReactantProd ( V, F, JVRP )
0.0%. T=0.000E+00 SO2= 0.5000E-01; O3= 0.0000E+00; NO= 0.1000E+00; NO2= 0.5000E-01;
2.1%. T=0.360E+04 SO2= 0.5000E-01; O3= 0.0000E+00; NO= 0.9983E-01; NO2= 0.5017E-01;
4.2%. T=0.720E+04 SO2= 0.5000E-01; O3= 0.0000E+00; NO= 0.9965E-01; NO2= 0.5035E-01;
6.2%. T=0.108E+05 SO2= 0.5000E-01; O3= 0.0000E+00; NO= 0.9948E-01; NO2= 0.5052E-01;
8.3%. T=0.144E+05 SO2= 0.5000E-01; O3= 0.0000E+00; NO= 0.9931E-01; NO2= 0.5069E-01;
10.4%. T=0.180E+05 SO2= 0.1141E-04; O3= 0.2865E-01; NO= 0.1259E-02; NO2= 0.5076E-02;
12.5%. T=0.216E+05 SO2=-0.2139E-46; O3= 0.1227E+00; NO= 0.9243E-02; NO2= 0.2088E-01;
14.6%. T=0.252E+05 SO2= 0.1280E-64; O3= 0.1646E+00; NO= 0.1297E-01; NO2= 0.2708E-01;
16.7%. T=0.288E+05 SO2=-0.3515E-81; O3= 0.1798E+00; NO= 0.1448E-01; NO2= 0.2916E-01;
18.8%. T=0.324E+05 SO2= 0.9231E-96; O3= 0.1850E+00; NO= 0.1504E-01; NO2= 0.2982E-01;
20.8%. T=0.360E+05 SO2= 0.2502-109; O3= 0.1867E+00; NO= 0.1523E-01; NO2= 0.3004E-01;
22.9%. T=0.396E+05 SO2= 0.5826-123; O3= 0.1872E+00; NO= 0.1529E-01; NO2= 0.3011E-01;
25.0%. T=0.432E+05 SO2= 0.1413-136; O3= 0.1874E+00; NO= 0.1530E-01; NO2= 0.3013E-01;
27.1%. T=0.468E+05 SO2= 0.3480-150; O3= 0.1874E+00; NO= 0.1530E-01; NO2= 0.3014E-01;
29.2%. T=0.504E+05 SO2= 0.8434-164; O3= 0.1873E+00; NO= 0.1529E-01; NO2= 0.3015E-01;
31.2%. T=0.540E+05 SO2= 0.2279-177; O3= 0.1871E+00; NO= 0.1520E-01; NO2= 0.3018E-01;
33.3%. T=0.576E+05 SO2= 0.1092-190; O3= 0.1862E+00; NO= 0.1496E-01; NO2= 0.3027E-01;
35.4%. T=0.612E+05 SO2= 0.2530-203; O3= 0.1839E+00; NO= 0.1439E-01; NO2= 0.3047E-01;
37.5%. T=0.648E+05 SO2=-0.6294-216; O3= 0.1760E+00; NO= 0.1303E-01; NO2= 0.3098E-01;
39.6%. T=0.684E+05 SO2= 0.3874-236; O3= 0.1337E+00; NO= 0.8058E-02; NO2= 0.3259E-01;
41.7%. T=0.720E+05 SO2= 0.5316-240; O3= 0.8543E-01; NO= 0.4916E-06; NO2= 0.2316E-01;
43.8%. T=0.756E+05 SO2= 0.5316-240; O3= 0.8029E-01; NO= 0.3615E-06; NO2= 0.1522E-01;
45.8%. T=0.792E+05 SO2= 0.5316-240; O3= 0.7693E-01; NO= 0.2658E-06; NO2= 0.1093E-01;
47.9%. T=0.828E+05 SO2= 0.5316-240; O3= 0.7454E-01; NO= 0.2023E-06; NO2= 0.8250E-02;
50.0%. T=0.864E+05 SO2= 0.5316-240; O3= 0.7276E-01; NO= 0.1574E-06; NO2= 0.6397E-02;
52.1%. T=0.900E+05 SO2= 0.5316-240; O3= 0.7141E-01; NO= 0.1241E-06; NO2= 0.5035E-02;
54.2%. T=0.936E+05 SO2= 0.5316-240; O3= 0.7035E-01; NO= 0.9856E-07; NO2= 0.3998E-02;
56.2%. T=0.972E+05 SO2= 0.5316-240; O3= 0.6952E-01; NO= 0.7872E-07; NO2= 0.3194E-02;
58.3%. T=0.101E+06 SO2= 0.5316-240; O3= 0.6886E-01; NO= 0.6312E-07; NO2= 0.2561E-02;
60.4%. T=0.104E+06 SO2= 0.2917-243; O3= 0.1270E+00; NO= 0.1039E-01; NO2= 0.3564E-01;
62.5%. T=0.108E+06 SO2= 0.3614-264; O3= 0.1760E+00; NO= 0.1306E-01; NO2= 0.3100E-01;
64.6%. T=0.112E+06 SO2= 0.1631-277; O3= 0.1839E+00; NO= 0.1440E-01; NO2= 0.3047E-01;
66.7%. T=0.115E+06 SO2=-0.3450-291; O3= 0.1862E+00; NO= 0.1497E-01; NO2= 0.3027E-01;
68.8%. T=0.119E+06 SO2=-0.7546-304; O3= 0.1871E+00; NO= 0.1520E-01; NO2= 0.3018E-01;
70.8%. T=0.122E+06 SO2= 0.6224-316; O3= 0.1874E+00; NO= 0.1529E-01; NO2= 0.3015E-01;
72.9%. T=0.126E+06 SO2=-0.0000E+00; O3= 0.1874E+00; NO= 0.1530E-01; NO2= 0.3014E-01;
75.0%. T=0.130E+06 SO2=-0.0000E+00; O3= 0.1874E+00; NO= 0.1531E-01; NO2= 0.3014E-01;
77.1%. T=0.133E+06 SO2=-0.0000E+00; O3= 0.1874E+00; NO= 0.1530E-01; NO2= 0.3014E-01;
79.2%. T=0.137E+06 SO2=-0.0000E+00; O3= 0.1873E+00; NO= 0.1529E-01; NO2= 0.3015E-01;
81.2%. T=0.140E+06 SO2=-0.0000E+00; O3= 0.1871E+00; NO= 0.1520E-01; NO2= 0.3018E-01;
83.3%. T=0.144E+06 SO2=-0.0000E+00; O3= 0.1862E+00; NO= 0.1496E-01; NO2= 0.3027E-01;
85.4%. T=0.148E+06 SO2=-0.0000E+00; O3= 0.1839E+00; NO= 0.1439E-01; NO2= 0.3047E-01;
87.5%. T=0.151E+06 SO2=-0.0000E+00; O3= 0.1760E+00; NO= 0.1303E-01; NO2= 0.3098E-01;
89.6%. T=0.155E+06 SO2=-0.0000E+00; O3= 0.1337E+00; NO= 0.8058E-02; NO2= 0.3259E-01;
91.7%. T=0.158E+06 SO2=-0.0000E+00; O3= 0.8543E-01; NO= 0.4916E-06; NO2= 0.2316E-01;
93.8%. T=0.162E+06 SO2=-0.0000E+00; O3= 0.8029E-01; NO= 0.3615E-06; NO2= 0.1522E-01;
95.8%. T=0.166E+06 SO2=-0.0000E+00; O3= 0.7693E-01; NO= 0.2658E-06; NO2= 0.1093E-01;
97.9%. T=0.169E+06 SO2=-0.0000E+00; O3= 0.7454E-01; NO= 0.2023E-06; NO2= 0.8249E-02;
100.0%. T=0.173E+06 SO2=-0.0000E+00; O3= 0.7276E-01; NO= 0.1574E-06; NO2= 0.6396E-02;
Results of a test using examples/saprc2006.kpp with KPP 2.2.3 (our starting version):
! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!
! Map File with Human-Readable Information
!
! Generated by KPP-2.2.3 symbolic chemistry Kinetics PreProcessor
! (http://www.cs.vt.edu/~asandu/Software/KPP)
! KPP is distributed under GPL, the general public licence
! (http://www.gnu.org/copyleft/gpl.html)
! (C) 1995-1997, V. Damian & A. Sandu, CGRER, Univ. Iowa
! (C) 1997-2005, A. Sandu, Michigan Tech, Virginia Tech
! With important contributions from:
! M. Damian, Villanova University, USA
! R. Sander, Max-Planck Institute for Chemistry, Mainz, Germany
!
! File : saprc2006.map
! Time : Wed Nov 9 10:14:44 2022
! Working directory : /local/ryantosca/GC/rundirs/epa-kpp/kpp_test/v2.2.3
! Equation file : saprc2006.kpp
! Output root filename : saprc2006
!
! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
### Options -------------------------------------------
FUNCTION - AGGREGATE
JACOBIAN - SPARSE W/ ACCOUNT FOR LU DECOMPOSITION FILL-IN
DOUBLE - ON
REORDER - ON
### Parameters ----------------------------------------
! NSPEC - Number of chemical species
INTEGER, PARAMETER :: NSPEC = 93
! NVAR - Number of Variable species
INTEGER, PARAMETER :: NVAR = 88
! NVARACT - Number of Active species
INTEGER, PARAMETER :: NVARACT = 81
! NFIX - Number of Fixed species
INTEGER, PARAMETER :: NFIX = 5
! NREACT - Number of reactions
INTEGER, PARAMETER :: NREACT = 235
! NVARST - Starting of variables in conc. vect.
INTEGER, PARAMETER :: NVARST = 1
! NFIXST - Starting of fixed in conc. vect.
INTEGER, PARAMETER :: NFIXST = 89
### Species -------------------------------------------
Variable species
1 = H2SO4 (r) 31 = N2O5 (r) 61 = ISOPRENE (r)
2 = BC (r) 32 = HONO (r) 62 = R2O2 (r)
3 = OC (r) 33 = ALK3 (r) 63 = TERP (r)
4 = SSF (r) 34 = TBU_O (r) 64 = METHACRO (r)
5 = SSC (r) 35 = ALK5 (r) 65 = OLE1 (r)
6 = PM10 (r) 36 = ARO2 (r) 66 = ISOPROD (r)
7 = PM25 (r) 37 = HNO4 (r) 67 = OLE2 (r)
8 = DMS (r) 38 = COOH (r) 68 = MVK (r)
9 = DST1 (r) 39 = HOCOO (r) 69 = CCHO (r)
10 = DST2 (r) 40 = BZNO2_O (r) 70 = HCHO (r)
11 = DST3 (r) 41 = MEOH (r) 71 = RNO3 (r)
12 = CO2 (r) 42 = ALK4 (r) 72 = O3P (r)
13 = HCOOH (n) 43 = ARO1 (r) 73 = RCHO (r)
14 = CCO_OH (n) 44 = DCB3 (r) 74 = MEK (r)
15 = RCO_OH (n) 45 = DCB2 (r) 75 = PROD2 (r)
16 = CCO_OOH (n) 46 = CRES (r) 76 = O3 (r)
17 = RCO_OOH (n) 47 = C2H2 (r) 77 = BZCO_O2 (r)
18 = XN (n) 48 = DCB1 (r) 78 = HO2 (r)
19 = XC (n) 49 = BALD (r) 79 = OH (r)
20 = SO2 (r) 50 = ROOH (r) 80 = MA_RCO3 (r)
21 = O1D (r) 51 = NPHE (r) 81 = NO (r)
22 = C2H6 (r) 52 = PHEN (r) 82 = NO3 (r)
23 = BACL (r) 53 = CO (r) 83 = C_O2 (r)
24 = PAN (r) 54 = MGLY (r) 84 = CCO_O2 (r)
25 = PAN2 (r) 55 = ACET (r) 85 = NO2 (r)
26 = PBZN (r) 56 = HNO3 (r) 86 = RO2_N (r)
27 = MA_PAN (r) 57 = ETHENE (r) 87 = RO2_R (r)
28 = H2O2 (r) 58 = C3H6 (r) 88 = RCO_O2 (r)
29 = C3H8 (r) 59 = GLY (r)
30 = ETOH (r) 60 = BZ_O (r)
Fixed species
1 = AIR (r) 3 = H2O (r) 5 = CH4 (r)
2 = O2 (r) 4 = H2 (r)
### Subroutines ---------------------------------------
SUBROUTINE Fun ( V, F, RCT, Vdot )
SUBROUTINE Jac_SP ( V, F, RCT, JVS )
SUBROUTINE Jac_SP_Vec ( JVS, UV, JUV )
SUBROUTINE JacTR_SP_Vec ( JVS, UV, JTUV )
SUBROUTINE KppSolve ( JVS, X )
SUBROUTINE KppSolveTR ( JVS, X, XX )
SUBROUTINE Hessian ( V, F, RCT, HESS )
SUBROUTINE HessTR_Vec ( HESS, U1, U2, HTU )
SUBROUTINE Hess_Vec ( HESS, U1, U2, HU )
SUBROUTINE Initialize ( )
SUBROUTINE Shuffle_user2kpp ( V_USER, V )
SUBROUTINE Shuffle_kpp2user ( V, V_USER )
SUBROUTINE Update_RCONST ( )
SUBROUTINE Update_PHOTO ( )
SUBROUTINE GetMass ( CL, Mass )
SUBROUTINE ReactantProd ( V, F, ARP )
SUBROUTINE JacReactantProd ( V, F, JVRP )
0.0%. T=0.000E+00 SO2= 0.5000E-01; O3= 0.0000E+00; NO= 0.1000E+00; NO2= 0.5000E-01;
2.1%. T=0.360E+04 SO2= 0.5000E-01; O3= 0.0000E+00; NO= 0.9983E-01; NO2= 0.5017E-01;
4.2%. T=0.720E+04 SO2= 0.5000E-01; O3= 0.0000E+00; NO= 0.9965E-01; NO2= 0.5035E-01;
6.2%. T=0.108E+05 SO2= 0.5000E-01; O3= 0.0000E+00; NO= 0.9948E-01; NO2= 0.5052E-01;
8.3%. T=0.144E+05 SO2= 0.5000E-01; O3= 0.0000E+00; NO= 0.9931E-01; NO2= 0.5069E-01;
10.4%. T=0.180E+05 SO2= 0.1141E-04; O3= 0.2865E-01; NO= 0.1259E-02; NO2= 0.5076E-02;
12.5%. T=0.216E+05 SO2=-0.2139E-46; O3= 0.1227E+00; NO= 0.9243E-02; NO2= 0.2088E-01;
14.6%. T=0.252E+05 SO2= 0.1280E-64; O3= 0.1646E+00; NO= 0.1297E-01; NO2= 0.2708E-01;
16.7%. T=0.288E+05 SO2=-0.3515E-81; O3= 0.1798E+00; NO= 0.1448E-01; NO2= 0.2916E-01;
18.8%. T=0.324E+05 SO2= 0.9231E-96; O3= 0.1850E+00; NO= 0.1504E-01; NO2= 0.2982E-01;
20.8%. T=0.360E+05 SO2= 0.2502-109; O3= 0.1867E+00; NO= 0.1523E-01; NO2= 0.3004E-01;
22.9%. T=0.396E+05 SO2= 0.5826-123; O3= 0.1872E+00; NO= 0.1529E-01; NO2= 0.3011E-01;
25.0%. T=0.432E+05 SO2= 0.1413-136; O3= 0.1874E+00; NO= 0.1530E-01; NO2= 0.3013E-01;
27.1%. T=0.468E+05 SO2= 0.3480-150; O3= 0.1874E+00; NO= 0.1530E-01; NO2= 0.3014E-01;
29.2%. T=0.504E+05 SO2= 0.8434-164; O3= 0.1873E+00; NO= 0.1529E-01; NO2= 0.3015E-01;
31.2%. T=0.540E+05 SO2= 0.2279-177; O3= 0.1871E+00; NO= 0.1520E-01; NO2= 0.3018E-01;
33.3%. T=0.576E+05 SO2= 0.1092-190; O3= 0.1862E+00; NO= 0.1496E-01; NO2= 0.3027E-01;
35.4%. T=0.612E+05 SO2= 0.2530-203; O3= 0.1839E+00; NO= 0.1439E-01; NO2= 0.3047E-01;
37.5%. T=0.648E+05 SO2=-0.6294-216; O3= 0.1760E+00; NO= 0.1303E-01; NO2= 0.3098E-01;
39.6%. T=0.684E+05 SO2= 0.3874-236; O3= 0.1337E+00; NO= 0.8058E-02; NO2= 0.3259E-01;
41.7%. T=0.720E+05 SO2= 0.5316-240; O3= 0.8543E-01; NO= 0.4916E-06; NO2= 0.2316E-01;
43.8%. T=0.756E+05 SO2= 0.5316-240; O3= 0.8029E-01; NO= 0.3615E-06; NO2= 0.1522E-01;
45.8%. T=0.792E+05 SO2= 0.5316-240; O3= 0.7693E-01; NO= 0.2658E-06; NO2= 0.1093E-01;
47.9%. T=0.828E+05 SO2= 0.5316-240; O3= 0.7454E-01; NO= 0.2023E-06; NO2= 0.8250E-02;
50.0%. T=0.864E+05 SO2= 0.5316-240; O3= 0.7276E-01; NO= 0.1574E-06; NO2= 0.6397E-02;
52.1%. T=0.900E+05 SO2= 0.5316-240; O3= 0.7141E-01; NO= 0.1241E-06; NO2= 0.5035E-02;
54.2%. T=0.936E+05 SO2= 0.5316-240; O3= 0.7035E-01; NO= 0.9856E-07; NO2= 0.3998E-02;
56.2%. T=0.972E+05 SO2= 0.5316-240; O3= 0.6952E-01; NO= 0.7872E-07; NO2= 0.3194E-02;
58.3%. T=0.101E+06 SO2= 0.5316-240; O3= 0.6886E-01; NO= 0.6312E-07; NO2= 0.2561E-02;
60.4%. T=0.104E+06 SO2= 0.2917-243; O3= 0.1270E+00; NO= 0.1039E-01; NO2= 0.3564E-01;
62.5%. T=0.108E+06 SO2= 0.3614-264; O3= 0.1760E+00; NO= 0.1306E-01; NO2= 0.3100E-01;
64.6%. T=0.112E+06 SO2= 0.1631-277; O3= 0.1839E+00; NO= 0.1440E-01; NO2= 0.3047E-01;
66.7%. T=0.115E+06 SO2=-0.3450-291; O3= 0.1862E+00; NO= 0.1497E-01; NO2= 0.3027E-01;
68.8%. T=0.119E+06 SO2=-0.7546-304; O3= 0.1871E+00; NO= 0.1520E-01; NO2= 0.3018E-01;
70.8%. T=0.122E+06 SO2= 0.6224-316; O3= 0.1874E+00; NO= 0.1529E-01; NO2= 0.3015E-01;
72.9%. T=0.126E+06 SO2=-0.0000E+00; O3= 0.1874E+00; NO= 0.1530E-01; NO2= 0.3014E-01;
75.0%. T=0.130E+06 SO2=-0.0000E+00; O3= 0.1874E+00; NO= 0.1531E-01; NO2= 0.3014E-01;
77.1%. T=0.133E+06 SO2=-0.0000E+00; O3= 0.1874E+00; NO= 0.1530E-01; NO2= 0.3014E-01;
79.2%. T=0.137E+06 SO2=-0.0000E+00; O3= 0.1873E+00; NO= 0.1529E-01; NO2= 0.3015E-01;
81.2%. T=0.140E+06 SO2=-0.0000E+00; O3= 0.1871E+00; NO= 0.1520E-01; NO2= 0.3018E-01;
83.3%. T=0.144E+06 SO2=-0.0000E+00; O3= 0.1862E+00; NO= 0.1496E-01; NO2= 0.3027E-01;
85.4%. T=0.148E+06 SO2=-0.0000E+00; O3= 0.1839E+00; NO= 0.1439E-01; NO2= 0.3047E-01;
87.5%. T=0.151E+06 SO2=-0.0000E+00; O3= 0.1760E+00; NO= 0.1303E-01; NO2= 0.3098E-01;
89.6%. T=0.155E+06 SO2=-0.0000E+00; O3= 0.1337E+00; NO= 0.8058E-02; NO2= 0.3259E-01;
91.7%. T=0.158E+06 SO2=-0.0000E+00; O3= 0.8543E-01; NO= 0.4916E-06; NO2= 0.2316E-01;
93.8%. T=0.162E+06 SO2=-0.0000E+00; O3= 0.8029E-01; NO= 0.3615E-06; NO2= 0.1522E-01;
95.8%. T=0.166E+06 SO2=-0.0000E+00; O3= 0.7693E-01; NO= 0.2658E-06; NO2= 0.1093E-01;
97.9%. T=0.169E+06 SO2=-0.0000E+00; O3= 0.7454E-01; NO= 0.2023E-06; NO2= 0.8249E-02;
100.0%. T=0.173E+06 SO2=-0.0000E+00; O3= 0.7276E-01; NO= 0.1574E-06; NO2= 0.6396E-02;
The concentrations at the end of the KPP 2.2.3 and KPP 3.0.0 test simulations are 100% identical (as they also are at other timesteps during the run):
KPP 2.2.3: 100.0%. T=0.173E+06 SO2=-0.0000E+00; O3= 0.7276E-01; NO= 0.1574E-06; NO2= 0.6396E-02;
KPP 3.0.0: 100.0%. T=0.173E+06 SO2=-0.0000E+00; O3= 0.7276E-01; NO= 0.1574E-06; NO2= 0.6396E-02;
Running into a little snag: I can't find the KPP 2.1 code, the links are broken at: https://people.cs.vt.edu/asandu/Software/Kpp/.
@RolfSander, would you happen to have a KPP 2.1 version on your end?
I would assume that KPP 2.2.3 would produce the same values as KPP 2.1. KPP 2.2.3 only contained "a number of small fixes throughout" as listed here: https://people.cs.vt.edu/asandu/Software/Kpp/
I sent you the code via email.
Thanks @RolfSander, I got the code.
One thing, the KPP 2.1 rosenbrock integrator seems different than the one in KPP 2.2.3. The 2.1 rosenbrock still has all the BLAS/LAPACK calls (WCOPY, WAXPY, etc) whereas the 2.2.3 one does not. I could try to use the 2.2.3 integrator in 2.1 and see if we get the same results.
Also the KPP 2.1 version has:
COMPILER = GNU
#COMPILER = LAHEY
#COMPILER = INTEL
#COMPILER = PGF
#COMPILER = HPUX
FC_GNU = g95
FOPT_GNU = -cpp -O -pg -fbounds-check
which I've replaced with:
FC_GNU = gfortran
FOPT_GNU = -cpp -O -g
to match the other versions.
I was not able to get the same results with KPP 2.1 as KPP 2.2.3:
! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!
! Map File with Human-Readable Information
!
! Generated by KPP-2.1 symbolic chemistry Kinetics PreProcessor
! (http://www.cs.vt.edu/~asandu/Software/KPP)
! KPP is distributed under GPL, the general public licence
! (http://www.gnu.org/copyleft/gpl.html)
! (C) 1995-1997, V. Damian & A. Sandu, CGRER, Univ. Iowa
! (C) 1997-2005, A. Sandu, Michigan Tech, Virginia Tech
! With important contributions from:
! M. Damian, Villanova University, USA
! R. Sander, Max-Planck Institute for Chemistry, Mainz, Germany
!
! File : saprc2006.map
! Time : Wed Nov 9 14:48:38 2022
! Working directory : /local/ryantosca/GC/rundirs/epa-kpp/kpp_test/v2.1
! Equation file : saprc2006.kpp
! Output root filename : saprc2006
!
! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
### Options -------------------------------------------
FUNCTION - AGGREGATE
JACOBIAN - SPARSE W/ ACCOUNT FOR LU DECOMPOSITION FILL-IN
DOUBLE - ON
REORDER - ON
### Parameters ----------------------------------------
! NSPEC - Number of chemical species
INTEGER, PARAMETER :: NSPEC = 93
! NVAR - Number of Variable species
INTEGER, PARAMETER :: NVAR = 88
! NVARACT - Number of Active species
INTEGER, PARAMETER :: NVARACT = 81
! NFIX - Number of Fixed species
INTEGER, PARAMETER :: NFIX = 5
! NREACT - Number of reactions
INTEGER, PARAMETER :: NREACT = 235
! NVARST - Starting of variables in conc. vect.
INTEGER, PARAMETER :: NVARST = 1
! NFIXST - Starting of fixed in conc. vect.
INTEGER, PARAMETER :: NFIXST = 89
### Species -------------------------------------------
Variable species
1 = H2SO4 (r) 31 = N2O5 (r) 61 = ISOPRENE (r)
2 = BC (r) 32 = HONO (r) 62 = R2O2 (r)
3 = OC (r) 33 = ALK3 (r) 63 = TERP (r)
4 = SSF (r) 34 = TBU_O (r) 64 = METHACRO (r)
5 = SSC (r) 35 = ALK5 (r) 65 = OLE1 (r)
6 = PM10 (r) 36 = ARO2 (r) 66 = ISOPROD (r)
7 = PM25 (r) 37 = HNO4 (r) 67 = OLE2 (r)
8 = DMS (r) 38 = COOH (r) 68 = MVK (r)
9 = DST1 (r) 39 = HOCOO (r) 69 = CCHO (r)
10 = DST2 (r) 40 = BZNO2_O (r) 70 = HCHO (r)
11 = DST3 (r) 41 = MEOH (r) 71 = RNO3 (r)
12 = CO2 (r) 42 = ALK4 (r) 72 = O3P (r)
13 = HCOOH (n) 43 = ARO1 (r) 73 = RCHO (r)
14 = CCO_OH (n) 44 = DCB3 (r) 74 = MEK (r)
15 = RCO_OH (n) 45 = DCB2 (r) 75 = PROD2 (r)
16 = CCO_OOH (n) 46 = CRES (r) 76 = O3 (r)
17 = RCO_OOH (n) 47 = C2H2 (r) 77 = BZCO_O2 (r)
18 = XN (n) 48 = DCB1 (r) 78 = HO2 (r)
19 = XC (n) 49 = BALD (r) 79 = OH (r)
20 = SO2 (r) 50 = ROOH (r) 80 = MA_RCO3 (r)
21 = O1D (r) 51 = NPHE (r) 81 = NO (r)
22 = C2H6 (r) 52 = PHEN (r) 82 = NO3 (r)
23 = BACL (r) 53 = CO (r) 83 = C_O2 (r)
24 = PAN (r) 54 = MGLY (r) 84 = CCO_O2 (r)
25 = PAN2 (r) 55 = ACET (r) 85 = NO2 (r)
26 = PBZN (r) 56 = HNO3 (r) 86 = RO2_N (r)
27 = MA_PAN (r) 57 = ETHENE (r) 87 = RO2_R (r)
28 = H2O2 (r) 58 = C3H6 (r) 88 = RCO_O2 (r)
29 = C3H8 (r) 59 = GLY (r)
30 = ETOH (r) 60 = BZ_O (r)
Fixed species
1 = AIR (r) 3 = H2O (r) 5 = CH4 (r)
2 = O2 (r) 4 = H2 (r)
### Subroutines ---------------------------------------
SUBROUTINE Fun ( V, F, RCT, Vdot )
SUBROUTINE Jac_SP ( V, F, RCT, JVS )
SUBROUTINE Jac_SP_Vec ( JVS, UV, JUV )
SUBROUTINE JacTR_SP_Vec ( JVS, UV, JTUV )
SUBROUTINE KppSolve ( JVS, X )
SUBROUTINE KppSolveTR ( JVS, X, XX )
SUBROUTINE Hessian ( V, F, RCT, HESS )
SUBROUTINE HessTR_Vec ( HESS, U1, U2, HTU )
SUBROUTINE Hess_Vec ( HESS, U1, U2, HU )
SUBROUTINE Initialize ( )
SUBROUTINE Shuffle_user2kpp ( V_USER, V )
SUBROUTINE Shuffle_kpp2user ( V, V_USER )
SUBROUTINE Update_RCONST ( )
SUBROUTINE Update_PHOTO ( )
SUBROUTINE GetMass ( CL, Mass )
SUBROUTINE ReactantProd ( V, F, ARP )
SUBROUTINE JacReactantProd ( V, F, JVRP )
0.0%. T=0.000E+00 SO2= 0.5000E-01; O3= 0.0000E+00; NO= 0.1000E+00; NO2= 0.5000E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
2.1%. T=0.360E+04 SO2= 0.5000E-01; O3= 0.0000E+00; NO= 0.9983E-01; NO2= 0.5017E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
4.2%. T=0.720E+04 SO2= 0.5000E-01; O3= 0.0000E+00; NO= 0.9965E-01; NO2= 0.5035E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
6.2%. T=0.108E+05 SO2= 0.5000E-01; O3= 0.0000E+00; NO= 0.9948E-01; NO2= 0.5052E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
8.3%. T=0.144E+05 SO2= 0.5000E-01; O3= 0.0000E+00; NO= 0.9931E-01; NO2= 0.5069E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
10.4%. T=0.180E+05 SO2= 0.1236E-01; O3= 0.2350E-01; NO= 0.2659E-03; NO2= 0.2441E-02; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
12.5%. T=0.216E+05 SO2= 0.4272E-38; O3= 0.1346E+00; NO= 0.1336E-01; NO2= 0.2297E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
14.6%. T=0.252E+05 SO2=-0.4856E-58; O3= 0.1958E+00; NO= 0.2187E-01; NO2= 0.3225E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
16.7%. T=0.288E+05 SO2=-0.2390E-74; O3= 0.2182E+00; NO= 0.2620E-01; NO2= 0.3527E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
18.8%. T=0.324E+05 SO2=-0.2656E-90; O3= 0.2256E+00; NO= 0.2807E-01; NO2= 0.3614E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
20.8%. T=0.360E+05 SO2=-0.1588-105; O3= 0.2278E+00; NO= 0.2873E-01; NO2= 0.3637E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
22.9%. T=0.396E+05 SO2= 0.2011-120; O3= 0.2283E+00; NO= 0.2888E-01; NO2= 0.3643E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
25.0%. T=0.432E+05 SO2=-0.2110-134; O3= 0.2284E+00; NO= 0.2890E-01; NO2= 0.3645E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
27.1%. T=0.468E+05 SO2= 0.2253-148; O3= 0.2284E+00; NO= 0.2889E-01; NO2= 0.3645E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
29.2%. T=0.504E+05 SO2=-0.2445-162; O3= 0.2282E+00; NO= 0.2877E-01; NO2= 0.3645E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
31.2%. T=0.540E+05 SO2= 0.2861-176; O3= 0.2274E+00; NO= 0.2822E-01; NO2= 0.3645E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
33.3%. T=0.576E+05 SO2=-0.3693-190; O3= 0.2249E+00; NO= 0.2674E-01; NO2= 0.3643E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
35.4%. T=0.612E+05 SO2=-0.2009-203; O3= 0.2182E+00; NO= 0.2356E-01; NO2= 0.3631E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
37.5%. T=0.648E+05 SO2=-0.1326-214; O3= 0.2003E+00; NO= 0.1760E-01; NO2= 0.3572E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
39.6%. T=0.684E+05 SO2= 0.6355-232; O3= 0.1426E+00; NO= 0.6969E-02; NO2= 0.3300E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
41.7%. T=0.720E+05 SO2= 0.5765-236; O3= 0.9593E-01; NO= 0.4479E-06; NO2= 0.2105E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
43.8%. T=0.756E+05 SO2= 0.5764-236; O3= 0.9078E-01; NO= 0.3142E-06; NO2= 0.1320E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
45.8%. T=0.792E+05 SO2= 0.5764-236; O3= 0.8753E-01; NO= 0.2232E-06; NO2= 0.9166E-02; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
47.9%. T=0.828E+05 SO2= 0.5764-236; O3= 0.8528E-01; NO= 0.1654E-06; NO2= 0.6735E-02; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
50.0%. T=0.864E+05 SO2= 0.5764-236; O3= 0.8363E-01; NO= 0.1257E-06; NO2= 0.5103E-02; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
52.1%. T=0.900E+05 SO2= 0.5764-236; O3= 0.8240E-01; NO= 0.9706E-07; NO2= 0.3935E-02; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
54.2%. T=0.936E+05 SO2= 0.5764-236; O3= 0.8145E-01; NO= 0.7565E-07; NO2= 0.3066E-02; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
56.2%. T=0.972E+05 SO2= 0.5764-236; O3= 0.8072E-01; NO= 0.5936E-07; NO2= 0.2406E-02; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
58.3%. T=0.101E+06 SO2= 0.5764-236; O3= 0.8015E-01; NO= 0.4683E-07; NO2= 0.1899E-02; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
60.4%. T=0.104E+06 SO2= 0.1286-239; O3= 0.1364E+00; NO= 0.7741E-02; NO2= 0.3422E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
62.5%. T=0.108E+06 SO2= 0.7143-255; O3= 0.2000E+00; NO= 0.1763E-01; NO2= 0.3573E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
64.6%. T=0.112E+06 SO2= 0.3343-265; O3= 0.2182E+00; NO= 0.2356E-01; NO2= 0.3631E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
66.7%. T=0.115E+06 SO2=-0.2184-276; O3= 0.2249E+00; NO= 0.2674E-01; NO2= 0.3643E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
68.8%. T=0.119E+06 SO2= 0.5424-289; O3= 0.2274E+00; NO= 0.2822E-01; NO2= 0.3645E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
70.8%. T=0.122E+06 SO2=-0.1804-301; O3= 0.2282E+00; NO= 0.2877E-01; NO2= 0.3645E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
72.9%. T=0.126E+06 SO2= 0.6232-313; O3= 0.2284E+00; NO= 0.2889E-01; NO2= 0.3645E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
75.0%. T=0.130E+06 SO2= 0.4941-323; O3= 0.2284E+00; NO= 0.2890E-01; NO2= 0.3645E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
77.1%. T=0.133E+06 SO2=-0.0000E+00; O3= 0.2284E+00; NO= 0.2889E-01; NO2= 0.3645E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
79.2%. T=0.137E+06 SO2=-0.0000E+00; O3= 0.2282E+00; NO= 0.2877E-01; NO2= 0.3645E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
81.2%. T=0.140E+06 SO2=-0.0000E+00; O3= 0.2274E+00; NO= 0.2822E-01; NO2= 0.3645E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
83.3%. T=0.144E+06 SO2=-0.0000E+00; O3= 0.2249E+00; NO= 0.2674E-01; NO2= 0.3643E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
85.4%. T=0.148E+06 SO2=-0.0000E+00; O3= 0.2182E+00; NO= 0.2356E-01; NO2= 0.3631E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
87.5%. T=0.151E+06 SO2=-0.0000E+00; O3= 0.2003E+00; NO= 0.1760E-01; NO2= 0.3572E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
89.6%. T=0.155E+06 SO2=-0.0000E+00; O3= 0.1426E+00; NO= 0.6969E-02; NO2= 0.3300E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
91.7%. T=0.158E+06 SO2=-0.0000E+00; O3= 0.9593E-01; NO= 0.4479E-06; NO2= 0.2105E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
93.8%. T=0.162E+06 SO2=-0.0000E+00; O3= 0.9078E-01; NO= 0.3142E-06; NO2= 0.1320E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
95.8%. T=0.166E+06 SO2=-0.0000E+00; O3= 0.8753E-01; NO= 0.2232E-06; NO2= 0.9166E-02; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
97.9%. T=0.169E+06 SO2=-0.0000E+00; O3= 0.8528E-01; NO= 0.1654E-06; NO2= 0.6735E-02; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
100.0%. T=0.173E+06 SO2=-0.0000E+00; O3= 0.8363E-01; NO= 0.1257E-06; NO2= 0.5103E-02; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
even though I made sure that the saprcnov.eqn
, 'saprcnov.def, and
saprcnov.spc` files used by KPP 2.1 were the same as KPP 2.2.3.
Long story short, there seems to have been some signfiicant modifications between KPP 2.1 and KPP 2.2.3 that are difficult to back out.
This is interesting. The differences are not small. Have you done a code diff? Do you need help doing this?
Hi @msl3v, I've attached the diff for all the files but it's very significant between 2.1 and 2.2.3. @yantosca and I are comparing the generated code to see what's causing the differences. If you could help out it would be great as well!
Hello @RolfSander, it looks like we are seeing some significant differences between 2.1 and 2.2.3 output. Perhaps there were a lot of changes between these two versions? From 2.2.3 onwards (which is the last release on the original KPP website) we have zero differences, as verified by @yantosca.
There seems to be a significant divergence starting at 10.4%. T=0.180E+05
.
Was wondering if this is due to Update_Sun. In the KPP 2.1 code there might be a PI missing:
SUBROUTINE Update_SUN()
!USE KPP_ROOT_Parameters
!USE KPP_ROOT_Global
IMPLICIT NONE
KPP_REAL SunRise, SunSet
KPP_REAL Thour, Tlocal, Ttmp
SunRise = 4.5_dp
SunSet = 19.5_dp
Thour = TIME/3600.0_dp
Tlocal = Thour - (INT(Thour)/24)*24
IF ((Tlocal>=SunRise).AND.(Tlocal<=SunSet)) THEN
Ttmp = (2.0*Tlocal-SunRise-SunSet)/(SunSet-SunRise)
IF (Ttmp.GT.0) THEN
Ttmp = Ttmp*Ttmp
ELSE
Ttmp = -Ttmp*Ttmp
END IF
SUN = ( 1.0_dp + COS(PI*Ttmp) )/2.0_dp
ELSE
SUN = 0.0_dp
END IF
END SUBROUTINE Update_SUN
whereas in the KPP 2.2.3 code the PI is defined:
SUBROUTINE Update_SUN()
USE KPP_ROOT_Parameters
USE KPP_ROOT_Global
IMPLICIT NONE
KPP_REAL :: SunRise, SunSet
KPP_REAL :: Thour, Tlocal, Ttmp
! PI - Value of pi
KPP_REAL, PARAMETER :: PI = 3.14159265358979d0
SunRise = 4.5_dp
SunSet = 19.5_dp
Thour = TIME/3600.0_dp
Tlocal = Thour - (INT(Thour)/24)*24
IF ((Tlocal>=SunRise).AND.(Tlocal<=SunSet)) THEN
Ttmp = (2.0*Tlocal-SunRise-SunSet)/(SunSet-SunRise)
IF (Ttmp.GT.0) THEN
Ttmp = Ttmp*Ttmp
ELSE
Ttmp = -Ttmp*Ttmp
END IF
SUN = ( 1.0_dp + COS(PI*Ttmp) )/2.0_dp
ELSE
SUN = 0.0_dp
END IF
END SUBROUTINE Update_SUN
Adding the PI doesn't change the end value:
100.0%. T=0.173E+06 SO2=-0.0000E+00; O3= 0.8363E-01; NO= 0.1257E-06; NO2= 0.5103E-02; = 0.0000E+00
It appears that the culprit is in the KPP-2.1 version of Rates.F90
. The issue is this line:
RCONST(38) = (EP3(3.08e-34,-2800.0,2.59e-38,-3180.0))
which is generated correctly in KPP-2.2 and above, but is generated as
RCONST(38) = (EP3(3.08e-34,-2800.0,2.59e-54,-3180.0))
The 2.59e-54
also is underflowing the real type. I'm not sure why it does not generate in KPP-2.1 correctly. @yantosca, I assume you used the same input files (.eqn?) for saprcnov in KPP-2.1 and KPP-2.2.3? Maybe it is just an input problem.
The diff between KPP-2.1 and KPP-2.2.3 is big! Maybe we should simply point out in the paper that our project started with 2.2.3 and still produces the same results.
We could warn users who still use 2.1 that they will see some differences caused by changed made until 2.2.3.
@jimmielin @RolfSander: I see the issue, in the saprcnov.eqn in KPP 2.1, we have these lines:
> //{38} HO2 + HO2 + H2O = H2O2 : EP3(3.08e-34,-2800.0,2.59e-54,-3180.0);
> {38} HO2 + HO2 + H2O = H2O2 : EP3(3.08e-34,-2800.0,2.59e-38,-3180.0);
Let me try to rebuild the mechanism again and make sure that it's 2.59e-54.
I've rebuilt the test mechansim making sure that the rxn rates in KPP v2.1 and v2.2.3 are identical:
If we diff the saprc2006_Rates.f90 files in v2.1 and 2.2.3 we no loner get any differences in RCONST:
5c5
< ! Generated by KPP-2.1 symbolic chemistry Kinetics PreProcessor
---
> ! Generated by KPP-2.2.3 symbolic chemistry Kinetics PreProcessor
16,17c16,17
< ! Time : Thu Nov 10 10:01:10 2022
< ! Working directory : /local/ryantosca/GC/rundirs/epa-kpp/kpp_test/v2.1
---
> ! Time : Wed Nov 9 10:14:44 2022
> ! Working directory : /local/ryantosca/GC/rundirs/epa-kpp/kpp_test/v2.2.3
58c58
< REAL(dp) K0,K2,K3
---
> REAL(kind=dp) K0,K2,K3
68c68
< REAL(dp) K1, K2
---
> REAL(kind=dp) K1, K2
76c76
< REAL(dp) K0, K1
---
> REAL(kind=dp) K0, K1
87c87
< ELEMENTAL REAL(dp) FUNCTION k_3rd(temp,cair,k0_300K,n,kinf_300K,m,fc)
---
> ELEMENTAL REAL(kind=dp) FUNCTION k_3rd(temp,cair,k0_300K,n,kinf_300K,m,fc)
91,98c91,98
< REAL(dp), INTENT(IN) :: temp ! temperature [K]
< REAL(dp), INTENT(IN) :: cair ! air concentration [molecules/cm3]
< REAL, INTENT(IN) :: k0_300K ! low pressure limit at 300 K
< REAL, INTENT(IN) :: n ! exponent for low pressure limit
< REAL, INTENT(IN) :: kinf_300K ! high pressure limit at 300 K
< REAL, INTENT(IN) :: m ! exponent for high pressure limit
< REAL, INTENT(IN) :: fc ! broadening factor (usually fc=0.6)
< REAL :: zt_help, k0_T, kinf_T, k_ratio
---
> REAL(kind=dp), INTENT(IN) :: temp ! temperature [K]
> REAL(kind=dp), INTENT(IN) :: cair ! air concentration [molecules/cm3]
> REAL, INTENT(IN) :: k0_300K ! low pressure limit at 300 K
> REAL, INTENT(IN) :: n ! exponent for low pressure limit
> REAL, INTENT(IN) :: kinf_300K ! high pressure limit at 300 K
> REAL, INTENT(IN) :: m ! exponent for high pressure limit
> REAL, INTENT(IN) :: fc ! broadening factor (usually fc=0.6)
> REAL(kind=dp) :: zt_help, k0_T, kinf_T, k_ratio
110c110
< ELEMENTAL REAL(dp) FUNCTION k_arr (k_298,tdep,temp)
---
> ELEMENTAL REAL(kind=dp) FUNCTION k_arr (k_298,tdep,temp)
115c115
< REAL(dp), INTENT(IN) :: temp ! temperature
---
> REAL(kind=dp), INTENT(IN) :: temp ! temperature
148,150c148,152
< REAL(kind=dp) SunRise, SunSet
< REAL(kind=dp) Thour, Tlocal, Ttmp
<
---
> REAL(kind=dp) :: SunRise, SunSet
> REAL(kind=dp) :: Thour, Tlocal, Ttmp
> ! PI - Value of pi
> REAL(kind=dp), PARAMETER :: PI = 3.14159265358979d0
>
and the results are now identical to the v2.2.3 and v3.0.0 runs:
! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!
! Map File with Human-Readable Information
!
! Generated by KPP-2.1 symbolic chemistry Kinetics PreProcessor
! (http://www.cs.vt.edu/~asandu/Software/KPP)
! KPP is distributed under GPL, the general public licence
! (http://www.gnu.org/copyleft/gpl.html)
! (C) 1995-1997, V. Damian & A. Sandu, CGRER, Univ. Iowa
! (C) 1997-2005, A. Sandu, Michigan Tech, Virginia Tech
! With important contributions from:
! M. Damian, Villanova University, USA
! R. Sander, Max-Planck Institute for Chemistry, Mainz, Germany
!
! File : saprc2006.map
! Time : Thu Nov 10 10:01:10 2022
! Working directory : /local/ryantosca/GC/rundirs/epa-kpp/kpp_test/v2.1
! Equation file : saprc2006.kpp
! Output root filename : saprc2006
!
! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
### Options -------------------------------------------
FUNCTION - AGGREGATE
JACOBIAN - SPARSE W/ ACCOUNT FOR LU DECOMPOSITION FILL-IN
DOUBLE - ON
REORDER - ON
### Parameters ----------------------------------------
! NSPEC - Number of chemical species
INTEGER, PARAMETER :: NSPEC = 93
! NVAR - Number of Variable species
INTEGER, PARAMETER :: NVAR = 88
! NVARACT - Number of Active species
INTEGER, PARAMETER :: NVARACT = 81
! NFIX - Number of Fixed species
INTEGER, PARAMETER :: NFIX = 5
! NREACT - Number of reactions
INTEGER, PARAMETER :: NREACT = 235
! NVARST - Starting of variables in conc. vect.
INTEGER, PARAMETER :: NVARST = 1
! NFIXST - Starting of fixed in conc. vect.
INTEGER, PARAMETER :: NFIXST = 89
### Species -------------------------------------------
Variable species
1 = H2SO4 (r) 31 = N2O5 (r) 61 = ISOPRENE (r)
2 = BC (r) 32 = HONO (r) 62 = R2O2 (r)
3 = OC (r) 33 = ALK3 (r) 63 = TERP (r)
4 = SSF (r) 34 = TBU_O (r) 64 = METHACRO (r)
5 = SSC (r) 35 = ALK5 (r) 65 = OLE1 (r)
6 = PM10 (r) 36 = ARO2 (r) 66 = ISOPROD (r)
7 = PM25 (r) 37 = HNO4 (r) 67 = OLE2 (r)
8 = DMS (r) 38 = COOH (r) 68 = MVK (r)
9 = DST1 (r) 39 = HOCOO (r) 69 = CCHO (r)
10 = DST2 (r) 40 = BZNO2_O (r) 70 = HCHO (r)
11 = DST3 (r) 41 = MEOH (r) 71 = RNO3 (r)
12 = CO2 (r) 42 = ALK4 (r) 72 = O3P (r)
13 = HCOOH (n) 43 = ARO1 (r) 73 = RCHO (r)
14 = CCO_OH (n) 44 = DCB3 (r) 74 = MEK (r)
15 = RCO_OH (n) 45 = DCB2 (r) 75 = PROD2 (r)
16 = CCO_OOH (n) 46 = CRES (r) 76 = O3 (r)
17 = RCO_OOH (n) 47 = C2H2 (r) 77 = BZCO_O2 (r)
18 = XN (n) 48 = DCB1 (r) 78 = HO2 (r)
19 = XC (n) 49 = BALD (r) 79 = OH (r)
20 = SO2 (r) 50 = ROOH (r) 80 = MA_RCO3 (r)
21 = O1D (r) 51 = NPHE (r) 81 = NO (r)
22 = C2H6 (r) 52 = PHEN (r) 82 = NO3 (r)
23 = BACL (r) 53 = CO (r) 83 = C_O2 (r)
24 = PAN (r) 54 = MGLY (r) 84 = CCO_O2 (r)
25 = PAN2 (r) 55 = ACET (r) 85 = NO2 (r)
26 = PBZN (r) 56 = HNO3 (r) 86 = RO2_N (r)
27 = MA_PAN (r) 57 = ETHENE (r) 87 = RO2_R (r)
28 = H2O2 (r) 58 = C3H6 (r) 88 = RCO_O2 (r)
29 = C3H8 (r) 59 = GLY (r)
30 = ETOH (r) 60 = BZ_O (r)
Fixed species
1 = AIR (r) 3 = H2O (r) 5 = CH4 (r)
2 = O2 (r) 4 = H2 (r)
### Subroutines ---------------------------------------
SUBROUTINE Fun ( V, F, RCT, Vdot )
SUBROUTINE Jac_SP ( V, F, RCT, JVS )
SUBROUTINE Jac_SP_Vec ( JVS, UV, JUV )
SUBROUTINE JacTR_SP_Vec ( JVS, UV, JTUV )
SUBROUTINE KppSolve ( JVS, X )
SUBROUTINE KppSolveTR ( JVS, X, XX )
SUBROUTINE Hessian ( V, F, RCT, HESS )
SUBROUTINE HessTR_Vec ( HESS, U1, U2, HTU )
SUBROUTINE Hess_Vec ( HESS, U1, U2, HU )
SUBROUTINE Initialize ( )
SUBROUTINE Shuffle_user2kpp ( V_USER, V )
SUBROUTINE Shuffle_kpp2user ( V, V_USER )
SUBROUTINE Update_RCONST ( )
SUBROUTINE Update_PHOTO ( )
SUBROUTINE GetMass ( CL, Mass )
SUBROUTINE ReactantProd ( V, F, ARP )
SUBROUTINE JacReactantProd ( V, F, JVRP )
0.0%. T=0.000E+00 SO2= 0.5000E-01; O3= 0.0000E+00; NO= 0.1000E+00; NO2= 0.5000E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
2.1%. T=0.360E+04 SO2= 0.5000E-01; O3= 0.0000E+00; NO= 0.9983E-01; NO2= 0.5017E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
4.2%. T=0.720E+04 SO2= 0.5000E-01; O3= 0.0000E+00; NO= 0.9965E-01; NO2= 0.5035E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
6.2%. T=0.108E+05 SO2= 0.5000E-01; O3= 0.0000E+00; NO= 0.9948E-01; NO2= 0.5052E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
8.3%. T=0.144E+05 SO2= 0.5000E-01; O3= 0.0000E+00; NO= 0.9931E-01; NO2= 0.5069E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
10.4%. T=0.180E+05 SO2= 0.1141E-04; O3= 0.2865E-01; NO= 0.1259E-02; NO2= 0.5076E-02; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
12.5%. T=0.216E+05 SO2=-0.2053E-46; O3= 0.1227E+00; NO= 0.9243E-02; NO2= 0.2088E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
14.6%. T=0.252E+05 SO2= 0.1228E-64; O3= 0.1646E+00; NO= 0.1297E-01; NO2= 0.2708E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
16.7%. T=0.288E+05 SO2=-0.3373E-81; O3= 0.1798E+00; NO= 0.1448E-01; NO2= 0.2916E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
18.8%. T=0.324E+05 SO2= 0.8857E-96; O3= 0.1850E+00; NO= 0.1504E-01; NO2= 0.2982E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
20.8%. T=0.360E+05 SO2= 0.2400-109; O3= 0.1867E+00; NO= 0.1523E-01; NO2= 0.3004E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
22.9%. T=0.396E+05 SO2= 0.5590-123; O3= 0.1872E+00; NO= 0.1529E-01; NO2= 0.3011E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
25.0%. T=0.432E+05 SO2= 0.1356-136; O3= 0.1874E+00; NO= 0.1530E-01; NO2= 0.3013E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
27.1%. T=0.468E+05 SO2= 0.3339-150; O3= 0.1874E+00; NO= 0.1530E-01; NO2= 0.3014E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
29.2%. T=0.504E+05 SO2= 0.8092-164; O3= 0.1873E+00; NO= 0.1529E-01; NO2= 0.3015E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
31.2%. T=0.540E+05 SO2= 0.2186-177; O3= 0.1871E+00; NO= 0.1520E-01; NO2= 0.3018E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
33.3%. T=0.576E+05 SO2= 0.1048-190; O3= 0.1862E+00; NO= 0.1496E-01; NO2= 0.3027E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
35.4%. T=0.612E+05 SO2= 0.2428-203; O3= 0.1839E+00; NO= 0.1439E-01; NO2= 0.3047E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
37.5%. T=0.648E+05 SO2=-0.6039-216; O3= 0.1760E+00; NO= 0.1303E-01; NO2= 0.3098E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
39.6%. T=0.684E+05 SO2= 0.3717-236; O3= 0.1337E+00; NO= 0.8058E-02; NO2= 0.3259E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
41.7%. T=0.720E+05 SO2= 0.5100-240; O3= 0.8543E-01; NO= 0.4916E-06; NO2= 0.2316E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
43.8%. T=0.756E+05 SO2= 0.5100-240; O3= 0.8029E-01; NO= 0.3615E-06; NO2= 0.1522E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
45.8%. T=0.792E+05 SO2= 0.5100-240; O3= 0.7693E-01; NO= 0.2658E-06; NO2= 0.1093E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
47.9%. T=0.828E+05 SO2= 0.5100-240; O3= 0.7454E-01; NO= 0.2023E-06; NO2= 0.8250E-02; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
50.0%. T=0.864E+05 SO2= 0.5100-240; O3= 0.7276E-01; NO= 0.1574E-06; NO2= 0.6397E-02; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
52.1%. T=0.900E+05 SO2= 0.5100-240; O3= 0.7141E-01; NO= 0.1241E-06; NO2= 0.5035E-02; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
54.2%. T=0.936E+05 SO2= 0.5100-240; O3= 0.7035E-01; NO= 0.9856E-07; NO2= 0.3998E-02; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
56.2%. T=0.972E+05 SO2= 0.5100-240; O3= 0.6952E-01; NO= 0.7872E-07; NO2= 0.3194E-02; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
58.3%. T=0.101E+06 SO2= 0.5100-240; O3= 0.6886E-01; NO= 0.6312E-07; NO2= 0.2561E-02; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
60.4%. T=0.104E+06 SO2= 0.2799-243; O3= 0.1270E+00; NO= 0.1039E-01; NO2= 0.3564E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
62.5%. T=0.108E+06 SO2= 0.3467-264; O3= 0.1760E+00; NO= 0.1306E-01; NO2= 0.3100E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
64.6%. T=0.112E+06 SO2= 0.1565-277; O3= 0.1839E+00; NO= 0.1440E-01; NO2= 0.3047E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
66.7%. T=0.115E+06 SO2=-0.3310-291; O3= 0.1862E+00; NO= 0.1497E-01; NO2= 0.3027E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
68.8%. T=0.119E+06 SO2=-0.7239-304; O3= 0.1871E+00; NO= 0.1520E-01; NO2= 0.3018E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
70.8%. T=0.122E+06 SO2= 0.5971-316; O3= 0.1874E+00; NO= 0.1529E-01; NO2= 0.3015E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
72.9%. T=0.126E+06 SO2=-0.0000E+00; O3= 0.1874E+00; NO= 0.1530E-01; NO2= 0.3014E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
75.0%. T=0.130E+06 SO2=-0.0000E+00; O3= 0.1874E+00; NO= 0.1531E-01; NO2= 0.3014E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
77.1%. T=0.133E+06 SO2=-0.0000E+00; O3= 0.1874E+00; NO= 0.1530E-01; NO2= 0.3014E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
79.2%. T=0.137E+06 SO2=-0.0000E+00; O3= 0.1873E+00; NO= 0.1529E-01; NO2= 0.3015E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
81.2%. T=0.140E+06 SO2=-0.0000E+00; O3= 0.1871E+00; NO= 0.1520E-01; NO2= 0.3018E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
83.3%. T=0.144E+06 SO2=-0.0000E+00; O3= 0.1862E+00; NO= 0.1496E-01; NO2= 0.3027E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
85.4%. T=0.148E+06 SO2=-0.0000E+00; O3= 0.1839E+00; NO= 0.1439E-01; NO2= 0.3047E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
87.5%. T=0.151E+06 SO2=-0.0000E+00; O3= 0.1760E+00; NO= 0.1303E-01; NO2= 0.3098E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
89.6%. T=0.155E+06 SO2=-0.0000E+00; O3= 0.1337E+00; NO= 0.8058E-02; NO2= 0.3259E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
91.7%. T=0.158E+06 SO2=-0.0000E+00; O3= 0.8543E-01; NO= 0.4916E-06; NO2= 0.2316E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
93.8%. T=0.162E+06 SO2=-0.0000E+00; O3= 0.8029E-01; NO= 0.3615E-06; NO2= 0.1522E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
95.8%. T=0.166E+06 SO2=-0.0000E+00; O3= 0.7693E-01; NO= 0.2658E-06; NO2= 0.1093E-01; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
97.9%. T=0.169E+06 SO2=-0.0000E+00; O3= 0.7454E-01; NO= 0.2023E-06; NO2= 0.8249E-02; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
100.0%. T=0.173E+06 SO2=-0.0000E+00; O3= 0.7276E-01; NO= 0.1574E-06; NO2= 0.6396E-02; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
The concentrations at the last timestep are now:
KPP 3.0.0 100.0%. T=0.173E+06 SO2=-0.0000E+00; O3= 0.7276E-01; NO= 0.1574E-06; NO2= 0.6396E-02;
KPP 2.2.3: 100.0%. T=0.173E+06 SO2=-0.0000E+00; O3= 0.7276E-01; NO= 0.1574E-06; NO2= 0.6396E-02;
KPP 2.1: 100.0%. T=0.173E+06 SO2=-0.0000E+00; O3= 0.7276E-01; NO= 0.1574E-06; NO2= 0.6396E-02; \0\0\0\0\0\0\0\0\0\0\0\0= 0.0000E+00;
NOTE: KPP 2.1 has an extra column of output that was removed from later versions.
So the answer is yes, we do get identical results. In this test, I forced KPP 2.1 to use the same integrator as KPP 2.2.3 (as the 2.1 integrator had all those BLAS/LAPACK calls), if that is a caveat.
Thanks for catching that @jimmielin. Sorry I didn't pick up on that!
I think we are good @jimmielin @RolfSander @msl3v.
Thanks for spotting this! I checked the origin for this number. It seems to come from Tab. 2f in doi 10.1029/97JD00849. Thus, "-54" is the correct exponent.
Thanks @yantosca and @RolfSander! Looks like we are good to go as it was an input issue, so KPP 3.0.0 with same inputs does have identical results all the way back to KPP 2.1. Great success!
Superb debugging.
I'll close out this issue now as it is resolved.
I have added this issue as a place where we can post results from KPP validation simulations with the sample mechanisms. This will allow us to determine if KPP 3.0.0 produces identical results w/r/t older versions.
Tagging @jimmielin @RolfSander @msl3v