KPP standalone interface from Obin Sturm

CHANGELOG.md

@@ -17,6 +17,7 @@ The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/),
 - Added computation of water concentration to use in photolysis for application of UV absorption by water in Cloud-J v8
 - Added ACO3, ACR, ACRO2, ALK4N{1,2,O}2, ALK4P, ALK7, APAN, APINN, APINO2, APINP, AROCMCHO, AROMCO3, AROMPN, BPINN, BPINO2, BPINON, BPINOO2, BPINOOH, BPINP, BUTN, BUTO2, C4H6, C96N, C96O2, C9602H, EBZ, GCO3, HACTA, LIMAL, LIMKB, LIMKET, LIMKO2, LIMN, LIMNB, LIMO2H, LIMO3, LIMO3H, LIMPAN, MEKCO3, MEKPN, MYRCO, PHAN, PIN, PINAL, PINO3, PINONIC, PINPAN, R7N{1,2}, R7O2, R7P, RNO3, STYR, TLFUO2, TLFUONE, TMB, ZRO2 to `species_database.yml` following Travis et  al. 2024.
 - Added TSOIL1 field to `State_Met` for use in HEMCO soil NOx extension. This should only be read in when the `UseSoilTemperature` option is true in HEMCO config.
+- Added KPP standalone interface
 
 ### Changed
 - Copy values from `State_Chm%KPP_AbsTol` to `ATOL` and `State_Chm%KPP_RelTol` to `RTOL` for fullchem and Hg simulations

GeosCore/CMakeLists.txt

@@ -82,6 +82,7 @@ add_library(GeosCore
   vdiff_mod.F90
   wetscav_mod.F90
   YuIMN_Code.F90
+  kpp_standalone_interface.F90
 
   # Files only included for special cases
   $<$<BOOL:${MODEL_CLASSIC}>:flexgrid_read_mod.F90 get_met_mod.F90 set_boundary_conditions_mod.F90>

GeosCore/fullchem_mod.F90

@@ -140,6 +140,7 @@ SUBROUTINE Do_FullChem( Input_Opt,  State_Chm, State_Diag,                 &
     USE UCX_MOD,                  ONLY : SO4_PHOTFRAC
     USE UCX_MOD,                  ONLY : UCX_NOX
     USE UCX_MOD,                  ONLY : UCX_H2SO4PHOT
+    USE KPP_Standalone_Interface
 #ifdef TOMAS
     USE TOMAS_MOD,                ONLY : H2SO4_RATE
     USE TOMAS_MOD,                ONLY : PSO4AQ_RATE
@@ -178,7 +179,7 @@ SUBROUTINE Do_FullChem( Input_Opt,  State_Chm, State_Diag,                 &
     INTEGER                :: errorCount, previous_units
     REAL(fp)               :: SO4_FRAC,   T,         TIN
     REAL(fp)               :: TOUT,       SR,        LWC
-
+    REAL(dp)               :: KPPH_before_integrate
     ! Strings
     CHARACTER(LEN=255)     :: errMsg,     thisLoc
 
@@ -200,6 +201,7 @@ SUBROUTINE Do_FullChem( Input_Opt,  State_Chm, State_Diag,                 &
     REAL(dp)               :: RCNTRL (20)
     REAL(dp)               :: RSTATE (20)
     REAL(dp)               :: C_before_integrate(NSPEC)
+    REAL(dp)               :: local_RCONST(NREACT)
 
     ! For tagged CO saving
     REAL(fp)               :: LCH4, PCO_TOT, PCO_CH4, PCO_NMVOC
@@ -235,6 +237,8 @@ SUBROUTINE Do_FullChem( Input_Opt,  State_Chm, State_Diag,                 &
     ! (assuming Rosenbrock solver).  Define this locally in order to break
     ! a compile-time dependency.  -- Bob Yantosca (05 May 2022)
     INTEGER,     PARAMETER :: Nhnew = 3
+    ! Add Nhexit, the last timestep length -- Obin Sturm (30 April 2024)
+    INTEGER,     PARAMETER :: Nhexit = 2
 
     ! Suppress printing out KPP error messages after this many errors occur
     INTEGER,     PARAMETER :: INTEGRATE_FAIL_TOGGLE = 20
@@ -438,6 +442,13 @@ SUBROUTINE Do_FullChem( Input_Opt,  State_Chm, State_Diag,                 &
        mapData => NULL()
     ENDIF
 
+    !=======================================================================
+    ! Should we print the full chemical state for any grid cell on this CPU?
+    ! for use with the KPP Standalone
+    ! (psturm, 03/22/24)
+    !=======================================================================
+    CALL Check_Domain( RC )
+
     !========================================================================
     ! Set up integration convergence conditions and timesteps
     ! (cf. M. J. Evans)
@@ -514,6 +525,7 @@ SUBROUTINE Do_FullChem( Input_Opt,  State_Chm, State_Diag,                 &
     !$OMP DEFAULT( SHARED                                                   )&
     !$OMP PRIVATE( I,        J,        L,       N                           )&
     !$OMP PRIVATE( ICNTRL,   C_before_integrate                             )&
+    !$OMP PRIVATE( KPPH_before_integrate,       local_RCONST                )&
     !$OMP PRIVATE( SO4_FRAC, IERR,     RCNTRL,  ISTATUS,   RSTATE           )&
     !$OMP PRIVATE( SpcID,    KppID,    F,       P,         Vloc             )&
     !$OMP PRIVATE( Aout,     Thread,   RC,      S,         LCH4             )&
@@ -569,6 +581,12 @@ SUBROUTINE Do_FullChem( Input_Opt,  State_Chm, State_Diag,                 &
        ! atmosphere if keepActive option is enabled. (hplin, 2/9/22)
        CALL fullchem_AR_SetKeepActive( option=.TRUE. )
 
+       ! Check if the current grid cell in this loop should have its
+       ! full chemical state printed (concentrations, rates, constants)
+       ! for use with the KPP Standalone
+       ! (psturm, 03/22/24)
+       CALL Check_ActiveCell( I, J, L, State_Grid )
+
        ! Start measuring KPP-related routine timing for this grid box
        IF ( State_Diag%Archive_KppTime ) THEN
           call cpu_time(TimeStart)
@@ -990,6 +1008,11 @@ SUBROUTINE Do_FullChem( Input_Opt,  State_Chm, State_Diag,                 &
        ! let us reset concentrations before calling "Integrate" a 2nd time.
        C_before_integrate = C
 
+       ! Do the same for the KPP initial timestep
+       ! Save local rate constants too
+       KPPH_before_integrate = State_Chm%KPPHvalue(I,J,L)
+       local_RCONST          = RCONST
+
        ! Call the Rosenbrock integrator
        ! (with optional auto-reduce functionality)
        CALL Integrate( TIN,    TOUT,    ICNTRL,                              &
@@ -1260,6 +1283,26 @@ SUBROUTINE Do_FullChem( Input_Opt,  State_Chm, State_Diag,                 &
           State_Diag%KppTime(I,J,L) = TimeEnd - TimeStart
        ENDIF
 
+       ! Write chemical state to file for the kpp standalone interface
+       ! No external logic needed, this subroutine exits early if the
+       ! chemical state should not be printed (psturm, 03/23/24)
+       CALL Write_Samples( I, J, L,      C_before_integrate,     &
+                           local_RCONST, KPPH_before_integrate,  &
+                           RSTATE(Nhexit),                       &
+                           State_Grid,   State_Chm,  State_Met,  &
+                           Input_Opt,    ISTATUS(3), RC )
+
+       ! test the force write option on the root node
+       ! example use case: printing chemical state under conditions
+       ! without knowing where those conditions will happen
+       ! IF ( Input_Opt%amIRoot .AND. L == 1 )                     &
+       ! CALL Write_Samples( I, J, L,      C_before_integrate,     &
+       !                     local_RCONST, KPPH_before_integrate,  &
+       !                     RSTATE(Nhexit),                       &
+       !                     State_Grid,   State_Chm,  State_Met,  &
+       !                     Input_Opt,    ISTATUS(3), RC,         &
+       !                     FORCE_WRITE = .TRUE., CELL_NAME = 'root')
+
        !=====================================================================
        ! Check we have no negative values and copy the concentrations
        ! calculated from the C array back into State_Chm%Species%Conc
@@ -2667,6 +2710,7 @@ SUBROUTINE Init_FullChem( Input_Opt, State_Chm, State_Diag, RC )
     USE State_Chm_Mod,            ONLY : ChmState
     USE State_Chm_Mod,            ONLY : Ind_
     USE State_Diag_Mod,           ONLY : DgnState
+    USE KPP_Standalone_Interface, ONLY : Config_KPP_Standalone
 !
 ! !INPUT PARAMETERS:
 !
@@ -2959,6 +3003,16 @@ SUBROUTINE Init_FullChem( Input_Opt, State_Chm, State_Diag, RC )
        ENDIF
     ENDIF
 
+    !--------------------------------------------------------------------
+    ! Initialize grid cells for input to KPP Standalone (Obin Sturm)
+    !--------------------------------------------------------------------
+    CALL Config_KPP_Standalone( Input_Opt, RC )
+    IF ( RC /= GC_SUCCESS ) THEN
+       ErrMsg = 'Error encountered in "KPP_Standalone"!'
+       CALL GC_Error( ErrMsg, RC, ThisLoc )
+       RETURN
+    ENDIF
+
   END SUBROUTINE Init_FullChem
 !EOC
 !------------------------------------------------------------------------------
@@ -2978,6 +3032,7 @@ SUBROUTINE Cleanup_FullChem( RC )
 ! !USES:
 !
     USE ErrCode_Mod
+    USE KPP_Standalone_Interface,   ONLY : Cleanup_KPP_Standalone
 !
 ! !OUTPUT PARAMETERS:
 !
@@ -3027,6 +3082,10 @@ SUBROUTINE Cleanup_FullChem( RC )
        IF ( RC /= GC_SUCCESS ) RETURN
     ENDIF
 
+    ! Deallocate variables from kpp standalone module
+    ! psturm, 03/22/2024
+    CALL Cleanup_KPP_Standalone( RC )
+
   END SUBROUTINE Cleanup_FullChem
 !EOC
 END MODULE FullChem_Mod