Skip to content

Commit

Permalink
Merge PR #2062 (Updates from GMAO for GEOS)
Browse files Browse the repository at this point in the history
This merge brings PR #2062 (Updates from GMAO for GEOS, by @lizziel) into
the GEOS-Chem "no-diff-to-benchmark" development stream.

This adds necessary updates in order to run GEOS-Chem in the NASA GEOS ESM,
and should not affect the GEOS-Chem Classic or GCHP implementations of
GEOS-Chem.

Signed-off-by: Bob Yantosca <[email protected]>
  • Loading branch information
yantosca committed Jan 5, 2024
2 parents def47a0 + 0c89ca1 commit aad13ea
Show file tree
Hide file tree
Showing 18 changed files with 445 additions and 118 deletions.
39 changes: 38 additions & 1 deletion GeosCore/fullchem_mod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -285,6 +285,8 @@ SUBROUTINE Do_FullChem( Input_Opt, State_Chm, State_Diag, &
IF (State_Diag%Archive_KppAutoReducerNVAR) &
State_Diag%KppAutoReducerNVAR = 0.0_f4
IF (State_Diag%Archive_KppcNONZERO) State_Diag%KppcNONZERO = 0.0_f4
IF (State_Diag%Archive_KppNegatives) State_Diag%KppNegatives = 0.0_f4
IF (State_Diag%Archive_KppNegatives0) State_Diag%KppNegatives0 = 0.0_f4
ENDIF

! Also zero satellite diagnostic archival arrays
Expand Down Expand Up @@ -999,20 +1001,42 @@ SUBROUTINE Do_FullChem( Input_Opt, State_Chm, State_Diag, &
ENDIF

#if defined( MODEL_GEOS ) || defined( MODEL_WRF ) || defined( MODEL_CESM )
! Keep track of error boxes
! Keep track of number of error boxes
IF ( State_Diag%Archive_KppError ) THEN
State_Diag%KppError(I,J,L) = State_Diag%KppError(I,J,L) + 1.0
ENDIF
#endif
ENDIF

#if defined( MODEL_GEOS )
! Mark integration as erroneous if negative concentrations so that it will be
! repeated below (cakelle2, 2023/10/26)
IF ( IERR >= 0 .AND. Input_Opt%KppCheckNegatives >= 0 ) THEN
IF ( ( Input_Opt%KppCheckNegatives==0 .AND. State_Met%InStratMeso(I,J,L) ) .OR. &
( L > ( State_Grid%NZ - Input_Opt%KppCheckNegatives) ) ) THEN
IF ( ANY(C < 0.0_dp) ) THEN
IERR = -999
! Include negative concentration boxes within error box diagnostic
IF ( State_Diag%Archive_KppError ) THEN
State_Diag%KppError(I,J,L) = State_Diag%KppError(I,J,L) + 1.0
ENDIF
ENDIF
ENDIF
ENDIF
#endif

!=====================================================================
! HISTORY: Archive KPP solver diagnostics
!
! !TODO: Abstract this into a separate routine
!=====================================================================
IF ( State_Diag%Archive_KppDiags ) THEN

! Check for negative concentrations after first integration
IF ( State_Diag%Archive_KppNegatives0 ) THEN
State_Diag%KppNegatives0(I,J,L) = REAL( COUNT( C < 0.0_dp ), KIND=4 )
ENDIF

! # of integrator calls
IF ( State_Diag%Archive_KppIntCounts ) THEN
State_Diag%KppIntCounts(I,J,L) = ISTATUS(1)
Expand Down Expand Up @@ -1070,6 +1094,12 @@ SUBROUTINE Do_FullChem( Input_Opt, State_Chm, State_Diag, &
RCNTRL(3) = 0.0_dp
C = C_before_integrate

#if defined( MODEL_GEOS )
! In GEOS also inflate the error tolerances (cakelle2, 2023/10/26)
ATOL = 1.0e-2_dp * Input_Opt%KppTolScale
RTOL = 1.0e-2_dp * Input_Opt%KppTolScale
#endif

! Disable auto-reduce solver for the second iteration for safety
IF ( Input_Opt%Use_AutoReduce ) THEN
RCNTRL(12) = -1.0_dp ! without using ICNTRL
Expand Down Expand Up @@ -1235,6 +1265,13 @@ SUBROUTINE Do_FullChem( Input_Opt, State_Chm, State_Diag, &
! Skip if this is not a GEOS-Chem species
IF ( SpcID <= 0 ) CYCLE

! Scan for negatives
IF ( State_Diag%Archive_KppNegatives ) THEN
IF ( C(N) < 0.0_dp ) THEN
State_Diag%KppNegatives(I,J,L) = State_Diag%KppNegatives(I,J,L) + 1.0_f4
ENDIF
ENDIF

! Set negative concentrations to zero
C(N) = MAX( C(N), 0.0_dp )

Expand Down
5 changes: 4 additions & 1 deletion Headers/input_opt_mod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -411,7 +411,10 @@ MODULE Input_Opt_Mod
LOGICAL :: FJX_EXTRAL_ERR = .TRUE.
! Toggle for het rates. If true, turns off three Cl producing het reactions
! in the stratosphere. In MODEL_GEOS, this flag is set in GEOSCHEMchem_GridComp.rc
LOGICAL :: TurnOffHetRates = .FALSE.
LOGICAL :: TurnOffHetRates = .TRUE.
INTEGER :: KppCheckNegatives = -1 ! Check for negatives after KPP integration
REAL(fp) :: KppTolScale = 1.0_fp ! Tolerance scale factor for 2nd KPP integration
LOGICAL :: applyQtend = .FALSE. ! Apply water vapor tendency
#else
LOGICAL :: AlwaysSetH2O
LOGICAL :: TurnOffHetRates
Expand Down
84 changes: 83 additions & 1 deletion Headers/state_diag_mod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -669,6 +669,12 @@ MODULE State_Diag_Mod
REAL(f4), POINTER :: KppSmDecomps(:,:,:)
LOGICAL :: Archive_KppSmDecomps

REAL(f4), POINTER :: KppNegatives(:,:,:)
LOGICAL :: Archive_KppNegatives

REAL(f4), POINTER :: KppNegatives0(:,:,:)
LOGICAL :: Archive_KppNegatives0

!%%%%% KPP auto-reduce solver diagnostics %%%%%
REAL(f4), POINTER :: KppAutoReducerNVAR(:,:,:)
LOGICAL :: Archive_KppAutoReducerNVAR
Expand Down Expand Up @@ -1895,6 +1901,12 @@ SUBROUTINE Zero_State_Diag( State_Diag, RC )
State_Diag%KppSmDecomps => NULL()
State_Diag%Archive_KppSmDecomps = .FALSE.

State_Diag%KppNegatives => NULL()
State_Diag%Archive_KppNegatives = .FALSE.

State_Diag%KppNegatives0 => NULL()
State_Diag%Archive_KppNegatives0 = .FALSE.

State_Diag%KppAutoReducerNVAR => NULL()
State_Diag%Archive_KppAutoReducerNVAR = .FALSE.

Expand Down Expand Up @@ -5996,6 +6008,50 @@ SUBROUTINE Init_State_Diag( Input_Opt, State_Chm, State_Grid, &
RETURN
ENDIF

!-------------------------------------------------------------------
! Number of negative concentrations after KPP integration
!-------------------------------------------------------------------
diagID = 'KppNegatives'
CALL Init_and_Register( &
Input_Opt = Input_Opt, &
State_Chm = State_Chm, &
State_Diag = State_Diag, &
State_Grid = State_Grid, &
DiagList = Diag_List, &
TaggedDiagList = TaggedDiag_List, &
Ptr2Data = State_Diag%KppNegatives, &
archiveData = State_Diag%Archive_KppNegatives, &
diagId = diagId, &
RC = RC )

IF ( RC /= GC_SUCCESS ) THEN
errMsg = TRIM( errMsg_ir ) // TRIM( diagId )
CALL GC_Error( errMsg, RC, thisLoc )
RETURN
ENDIF

!-------------------------------------------------------------------
! Number of negative concentrations after first KPP integration try
!-------------------------------------------------------------------
diagID = 'KppNegatives0'
CALL Init_and_Register( &
Input_Opt = Input_Opt, &
State_Chm = State_Chm, &
State_Diag = State_Diag, &
State_Grid = State_Grid, &
DiagList = Diag_List, &
TaggedDiagList = TaggedDiag_List, &
Ptr2Data = State_Diag%KppNegatives0, &
archiveData = State_Diag%Archive_KppNegatives0, &
diagId = diagId, &
RC = RC )

IF ( RC /= GC_SUCCESS ) THEN
errMsg = TRIM( errMsg_ir ) // TRIM( diagId )
CALL GC_Error( errMsg, RC, thisLoc )
RETURN
ENDIF

!-------------------------------------------------------------------
! AR only -- Number of species in reduced mechanism (NVAR - NRMV)
!-------------------------------------------------------------------
Expand Down Expand Up @@ -6119,7 +6175,7 @@ SUBROUTINE Init_State_Diag( Input_Opt, State_Chm, State_Grid, &
! being requested as diagnostic output when the corresponding
! array has not been allocated.
!-------------------------------------------------------------------
DO N = 1, 35
DO N = 1, 41
! Select the diagnostic ID
SELECT CASE( N )
CASE( 1 )
Expand Down Expand Up @@ -6200,6 +6256,10 @@ SUBROUTINE Init_State_Diag( Input_Opt, State_Chm, State_Grid, &
diagID = 'KppAutoReduceThres'
CASE( 39 )
diagID = 'RxnConst'
CASE( 40 )
diagID = 'KppNegatives'
CASE( 41 )
diagID = 'KppNegatives0'
END SELECT

! Exit if any of the above are in the diagnostic list
Expand Down Expand Up @@ -10460,6 +10520,8 @@ SUBROUTINE Init_State_Diag( Input_Opt, State_Chm, State_Grid, &
State_Diag%Archive_KppLuDecomps .or. &
State_Diag%Archive_KppSubsts .or. &
State_Diag%Archive_KppSmDecomps .or. &
State_Diag%Archive_KppNegatives .or. &
State_Diag%Archive_KppNegatives0 .or. &
State_Diag%Archive_KppAutoReducerNVAR .or. &
State_Diag%Archive_KppAutoReduceThres .or. &
State_Diag%Archive_KppcNONZERO .or. &
Expand Down Expand Up @@ -11893,6 +11955,16 @@ SUBROUTINE Cleanup_State_Diag( State_Diag, RC )
RC = RC )
IF ( RC /= GC_SUCCESS ) RETURN

CALL Finalize( diagId = 'KppNegatives', &
Ptr2Data = State_Diag%KppNegatives, &
RC = RC )
IF ( RC /= GC_SUCCESS ) RETURN

CALL Finalize( diagId = 'KppNegatives0', &
Ptr2Data = State_Diag%KppNegatives0, &
RC = RC )

IF ( RC /= GC_SUCCESS ) RETURN
CALL Finalize( diagId = 'AirMassColumnFull', &
Ptr2Data = State_Diag%AirMassColumnFull, &
RC = RC )
Expand Down Expand Up @@ -13562,6 +13634,16 @@ SUBROUTINE Get_Metadata_State_Diag( am_I_Root, metadataID, Found, &
IF ( isUnits ) Units = 'count'
IF ( isRank ) Rank = 3

ELSE IF ( TRIM( Name_AllCaps ) == 'KPPNEGATIVES' ) THEN
IF ( isDesc ) Desc = 'Number of negative concentrations after KPP integration'
IF ( isUnits ) Units = 'count'
IF ( isRank ) Rank = 3

ELSE IF ( TRIM( Name_AllCaps ) == 'KPPNEGATIVES0' ) THEN
IF ( isDesc ) Desc = 'Number of negative concentrations after first KPP integration attempt'
IF ( isUnits ) Units = 'count'
IF ( isRank ) Rank = 3

ELSE IF ( TRIM( Name_AllCaps ) == 'KPPAUTOREDUCERNVAR' ) THEN
IF ( isDesc ) Desc = 'Number of species in auto-reduced mechanism'
IF ( isUnits ) Units = 'count'
Expand Down
Loading

0 comments on commit aad13ea

Please sign in to comment.