Skip to content
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

H2 modelling capability #1871

Open
wants to merge 2 commits into
base: main
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
98 changes: 90 additions & 8 deletions GeosCore/drydep_mod.F
Original file line number Diff line number Diff line change
Expand Up @@ -193,6 +193,13 @@ MODULE DRYDEP_MOD
! 16 Mar 2017 - R. Yantosca - Remove references to obsolete vars in Input_Opt
! 24 Aug 2017 - M. Sulprizio- Remove support for GCAP, GEOS-4, GEOS-5 and MERRA
! 01 Dec 2018 - H.P. Lin - Move DEPSAV to State_Chm%DryDepSav
! 26 Jun 2020 - M.P Perez-Pena - Added special case for the dry
! deposition of H2.
! 21 Sep 2020 - M.P Perez-Pena - Added call to dry deposition velocity
! fields H2_DEPVEL from the work developed
! by Yashiro et al. (2011) with the two layer soil
! model to account for soil uptake of H2.
! See special case for H2.
!EOP
!------------------------------------------------------------------------------
!BOC
Expand Down Expand Up @@ -250,6 +257,7 @@ MODULE DRYDEP_MOD
INTEGER :: id_ACET, id_ALD2
INTEGER :: id_NK1
INTEGER :: id_HNO3, id_PAN, id_ISOPND
INTEGER :: id_H2

! Arrays for Baldocchi drydep polynomial coefficients
REAL(fp), TARGET :: DRYCOEFF(NPOLY )
Expand Down Expand Up @@ -287,6 +295,7 @@ MODULE DRYDEP_MOD
! Allocatable arrays
REAL(f8), ALLOCATABLE :: DMID (: )
REAL(f8), ALLOCATABLE :: SALT_V (: )
REAL(f8), ALLOCATABLE :: H2_DEPVEL(:, :) ! H2 dd vel [cm/s]

!=================================================================
! MODULE ROUTINES -- follow below the "CONTAINS" statement
Expand Down Expand Up @@ -316,6 +325,7 @@ SUBROUTINE DO_DRYDEP( am_I_Root, Input_Opt, State_Chm,
USE ErrCode_Mod
USE HCO_ERROR_MOD
USE HCO_INTERFACE_MOD, ONLY : HcoState
USE HCO_Calc_Mod, ONLY : HCO_EvalFld
USE HCO_DIAGN_MOD, ONLY : Diagn_Update
USE Input_Opt_Mod, ONLY : OptInput
USE Species_Mod, ONLY : Species
Expand Down Expand Up @@ -430,14 +440,18 @@ SUBROUTINE DO_DRYDEP( am_I_Root, Input_Opt, State_Chm,
!
! Scalars
INTEGER :: I, J, L, D, N, NDVZ
REAL(f8) :: DVZ, THIK
REAL(f8) :: DVZ, THIK, DVH2
CHARACTER(LEN=255) :: ErrMsg, ThisLoc

CHARACTER(LEN=255) :: LOC='DRY DEP (drydep_mod.F)'
REAL(fp) :: YMID

! Arrays
LOGICAL :: LSNOW (State_Grid%NX,State_Grid%NY) ! Flag for snow/ice on the sfc
LOGICAL :: DESERT(State_Grid%NX,State_Grid%NY) ! Flag for denoting a desert area
REAL(f8) :: CZ1 (State_Grid%NX,State_Grid%NY) ! Midpt ht of 1st model lev[m]
REAL(f8) :: TC0 (State_Grid%NX,State_Grid%NY) ! Grid box sfc temperature [K]
REAL(f8) :: FLAND (State_Grid%NX,State_Grid%NY) ! Fraction of land in the grid
REAL(f8) :: FSNOW (State_Grid%NX,State_Grid%NY) ! Fraction of snow in the grid
REAL(f8) :: ZH (State_Grid%NX,State_Grid%NY) ! PBL height [m]
REAL(f8) :: OBK (State_Grid%NX,State_Grid%NY) ! Monin-Obhukov Length [m]
REAL(f8) :: CFRAC (State_Grid%NX,State_Grid%NY) ! Column cloud frac [unitless]
Expand Down Expand Up @@ -485,7 +499,8 @@ SUBROUTINE DO_DRYDEP( am_I_Root, Input_Opt, State_Chm,
! Added sfc pressure as PRESSU and 10m windspeed as W10
! (jaegle 5/11/11, mpayer 1/10/12)
CALL METERO( State_Grid, State_Met, CZ1, TC0, OBK, CFRAC,
& RADIAT, AZO, USTAR, ZH, LSNOW,
& RADIAT, AZO, USTAR, ZH, LSNOW, DESERT,
& FLAND, FSNOW,
& RHB, PRESSU, W10, SUNCOS_MID )

! Call DEPVEL to compute dry deposition velocities [m/s]
Expand All @@ -503,18 +518,32 @@ SUBROUTINE DO_DRYDEP( am_I_Root, Input_Opt, State_Chm,
RETURN
ENDIF

! Call the dry deposition velocities for H2
CALL HCO_EvalFld( am_I_Root, HcoState, 'H2_DEPVEL',
& H2_DEPVEL, RC )

IF ( RC /= GC_SUCCESS ) THEN
ErrMsg = 'Could not find offline H2DV in HEMCO data list!'
CALL GC_Error( ErrMsg, RC, LOC )
RETURN
ENDIF

!=================================================================
! Compute dry deposition frequencies; archive diagnostics
!=================================================================
!$OMP PARALLEL DO
!$OMP+DEFAULT( SHARED )
!$OMP+PRIVATE( I, J, THIK, D, N, NDVZ, DVZ, SpcInfo )
!$OMP+PRIVATE( YMID, DVH2 )
DO J = 1, State_Grid%NY
DO I = 1, State_Grid%NX

! THIK = thickness of surface layer [m]
THIK = State_Met%BXHEIGHT(I,J,1)

! Latitude of grid box
YMID = State_Grid%YMid(I,J)

! Now we calculate drydep throughout the entire PBL.
! Make sure that the PBL depth is greater than or equal
! to the thickness of the 1st layer (rjp, bmy, 7/21/03)
Expand Down Expand Up @@ -624,6 +653,23 @@ SUBROUTINE DO_DRYDEP( am_I_Root, Input_Opt, State_Chm,
ENDIF
ENDIF

!-----------------------------------------------------------
! Special treatment for H2
! Sep. 21 2020
!-----------------------------------------------------------
! Added special case to read the dry deposition of H2 using
! HEMCO. HEMCO re-grids the information of the spatial
! distributed dry deposition velocity fields provided by
! Yashiro et al. (2011) (result from the implementation of a
! two layer soil model). The fields (netCDF files 2D data
! with a resolution of T42 (~2.8deg, 64x128 grid) are
! specified in the NON-EMISSION-DATA section in the
! HEMCO_Config.rc file.

IF ( N == id_H2 ) THEN
DVZ = H2_DEPVEL(I, J)
ENDIF

!-----------------------------------------------------------
! Compute drydep frequency and update diagnostics
!-----------------------------------------------------------
Expand Down Expand Up @@ -682,7 +728,8 @@ END SUBROUTINE DO_DRYDEP
! !INTERFACE:
!
SUBROUTINE METERO( State_Grid, State_Met, CZ1, TC0, OBK, CFRAC,
& RADIAT, AZO, USTR, ZH, LSNOW,
& RADIAT, AZO, USTR, ZH, LSNOW, DESERT,
& FLAND, FSNOW,
& RHB, PRESSU, W10, SUNCOS_MID )
!
! !USES:
Expand All @@ -698,8 +745,11 @@ SUBROUTINE METERO( State_Grid, State_Met, CZ1, TC0, OBK, CFRAC,
! !OUTPUT PARAMETERS:
!
LOGICAL, INTENT(OUT) :: LSNOW (State_Grid%NX,State_Grid%NY) ! Flag for denoting snow/ice
LOGICAL, INTENT(OUT) :: DESERT(State_Grid%NX,State_Grid%NY) ! Flag for denoting a desert area
REAL(f8), INTENT(OUT) :: CZ1 (State_Grid%NX,State_Grid%NY) ! Midpt ht of 1st model lev [m]
REAL(f8), INTENT(OUT) :: TC0 (State_Grid%NX,State_Grid%NY) ! Grid box sfc temp [K]
REAL(f8), INTENT(OUT) :: FLAND (State_Grid%NX,State_Grid%NY) ! Fraction of land in the grid
REAL(f8), INTENT(OUT) :: FSNOW (State_Grid%NX,State_Grid%NY) ! Fraction of snow in the grid
REAL(f8), INTENT(OUT) :: OBK (State_Grid%NX,State_Grid%NY) ! Monin-Obhukov length [m]
REAL(f8), INTENT(OUT) :: CFRAC (State_Grid%NX,State_Grid%NY) ! Column cloud fraction [unitless]
REAL(f8), INTENT(OUT) :: RADIAT(State_Grid%NX,State_Grid%NY) ! Solar radiation @ ground [W/m2]
Expand Down Expand Up @@ -753,6 +803,10 @@ SUBROUTINE METERO( State_Grid, State_Met, CZ1, TC0, OBK, CFRAC,
! 25 Jul 2014 - R. Yantosca - Now remove reference to function SFCWINDSQR
! 17 May 2016 - M. Sulprizio- Remove IJLOOP and change dimension of output
! arrays from (MAXIJ) to (IIPAR,JJPAR)
! 27 Feb 2020 - M.P Perez-Pena - Include a separate definition for desert
! areas based on the GWETTOP variable
! 25 Jun 2020 - M.P Perez-Pena - Added fraction of land and fraction of
! snow
!EOP
!------------------------------------------------------------------------------
!BOC
Expand Down Expand Up @@ -801,6 +855,15 @@ SUBROUTINE METERO( State_Grid, State_Met, CZ1, TC0, OBK, CFRAC,
! Set logical LSNOW if snow and sea ice (ALBEDO > 0.4)
LSNOW(I,J) = ( State_Met%ALBD(I,J) > 0.4 )

! Set logical DESERT if low soil top wetness ( GWETTOP < 0.35 )
DESERT(I,J) = ( State_Met%GWETTOP(I,J) < 0.35 )

! Fraction of land
FLAND(I,J) = State_Met%FRLAND(I,J)

! Fraction of snow
FSNOW(I,J) = State_Met%FRSNO(I,J)

! Monin-Obhukov length [m]
OBK(I,J) = GET_OBK( I, J, State_Met )

Expand Down Expand Up @@ -851,12 +914,14 @@ END SUBROUTINE METERO
!\\
! !INTERFACE:
!
SUBROUTINE DEPVEL( am_I_Root, Input_Opt, State_Chm, State_Diag,

SUBROUTINE DEPVEL( am_I_Root, Input_Opt, State_Chm, State_Diag,
& State_Grid, State_Met, RADIAT, TEMP, SUNCOS,
& F0, HSTAR, XMW, AIROSOL,
& USTAR, CZ1, OBK, CFRAC,
& ZH, LSNOW, DVEL, ZO,
& F0, HSTAR, XMW, AIROSOL,
& USTAR, CZ1, OBK, CFRAC,
& ZH, LSNOW, DVEL, ZO,
& RHB, PRESSU, W10, RC )

!
! !USES:
!
Expand Down Expand Up @@ -4024,6 +4089,8 @@ SUBROUTINE INIT_DRYDEP( am_I_Root, Input_Opt, State_Chm,
!
USE ErrCode_Mod
USE Input_Opt_Mod, ONLY : OptInput
USE HCO_INTERFACE_MOD, ONLY : HcoState
USE HCO_Calc_Mod, ONLY : HCO_EvalFld
USE Species_Mod, ONLY : Species
USE State_Chm_Mod, ONLY : ChmState
USE State_Chm_Mod, ONLY : Ind_
Expand Down Expand Up @@ -4118,6 +4185,7 @@ SUBROUTINE INIT_DRYDEP( am_I_Root, Input_Opt, State_Chm,
! PAN, and ISOPN (from K. Travis, J. Fisher)
! 05 Oct 2017 - R. Yantosca - Now accept State_Diag as an argument
! 07 Nov 2017 - R. Yantosca - Now return error condition to calling program
! 21 Sep 2020 - M. Perez-Pena- Add call to H2_DEPVEL
!EOP
!------------------------------------------------------------------------------
!BOC
Expand All @@ -4126,10 +4194,12 @@ SUBROUTINE INIT_DRYDEP( am_I_Root, Input_Opt, State_Chm,
!
LOGICAL :: LDRYD
LOGICAL :: IS_Hg
INTEGER :: AS, INDICES(4), YYYYMMDD, HHMMSS
INTEGER :: N

! Strings
CHARACTER(LEN=255) :: Msg, ErrMsg, ThisLoc
CHARACTER(LEN=255) :: LOC='DRY DEP (drydep_mod.F)'

! Objects
TYPE(Species), POINTER :: SpcInfo
Expand All @@ -4145,6 +4215,7 @@ SUBROUTINE INIT_DRYDEP( am_I_Root, Input_Opt, State_Chm,
NUMDEP = 0
id_ACET = 0
id_ALD2 = 0
id_H2 = 0
id_HNO3 = IND_('HNO3' )
id_PAN = IND_('PAN' )
id_ISOPND = IND_('ISOPND')
Expand Down Expand Up @@ -4232,6 +4303,12 @@ SUBROUTINE INIT_DRYDEP( am_I_Root, Input_Opt, State_Chm,
IF ( RC /= GC_SUCCESS ) RETURN
XMW = 0e+0_f8

! Dry deposition velocity for H2 spatial fields
ALLOCATE( H2_DEPVEL(State_Grid%NX,State_Grid%NY),
& STAT=AS)
IF ( AS /= 0 ) CALL ALLOC_ERR( 'H2_DEPVEL' )
H2_DEPVEL = 0e+0_fp

ENDIF

!=================================================================
Expand Down Expand Up @@ -4308,6 +4385,10 @@ SUBROUTINE INIT_DRYDEP( am_I_Root, Input_Opt, State_Chm,
! Flag the species ID of ALD2 for use above.
id_ALD2 = SpcInfo%ModelId

CASE( 'H2' )
! Flag the species ID of H2 for use above.
id_H2 = SpcInfo%ModelId

CASE( 'NITs', 'NITS' )
! DEPNAME for NITs has to be in all caps, for
! backwards compatibility with older code.
Expand Down Expand Up @@ -4500,6 +4581,7 @@ SUBROUTINE CLEANUP_DRYDEP
IF ( ALLOCATED( NTRAIND ) ) DEALLOCATE( NTRAIND )
IF ( ALLOCATED( SALT_V ) ) DEALLOCATE( SALT_V )
IF ( ALLOCATED( XMW ) ) DEALLOCATE( XMW )
IF ( ALLOCATED( H2_DEPVEL ) ) DEALLOCATE( H2_DEPVEL )

END SUBROUTINE CLEANUP_DRYDEP
!EOC
Expand Down
26 changes: 11 additions & 15 deletions GeosCore/ucx_mod.F
Original file line number Diff line number Diff line change
Expand Up @@ -147,7 +147,7 @@ MODULE UCX_MOD

! OCS and H2 surface boundary conditions
REAL(fp),PARAMETER :: OCS_MR = 5.00e-10_fp ! OCS is fixed to 500 pptv
REAL(fp),PARAMETER :: H2_MR = 5.00e-7_fp ! H2 is fixed to 500 ppbv
!REAL(fp),PARAMETER :: H2_MR = 5.00e-7_fp ! H2 is fixed to 500 ppbv
!
! PRIVATE TYPES:
!
Expand Down Expand Up @@ -1437,9 +1437,9 @@ SUBROUTINE READ_SFC( am_I_Root, Input_Opt )
CASE ('OCS')
! Simple boundary condition
GRID_EMIT(:,:,N) = OCS_MR
CASE ('H2')
!CASE ('H2')
! Simple boundary condition
GRID_EMIT(:,:,N) = H2_MR
! GRID_EMIT(:,:,N) = H2_MR
CASE ('CH4')
! Now obtain monthly mean surface CH4 concentrations
! in main.F (mps, 1/18/18)
Expand Down Expand Up @@ -1492,9 +1492,9 @@ SUBROUTINE READ_SFC( am_I_Root, Input_Opt )
CASE ('OCS')
! Simple boundary condition
GRID_EMIT(:,:,N) = OCS_MR
CASE ('H2')
!CASE ('H2')
! Simple boundary condition
GRID_EMIT(:,:,N) = H2_MR
! GRID_EMIT(:,:,N) = H2_MR
CASE ('CH4')
! Now obtain monthly mean surface CH4 concentrations
! in main.F (mps, 1/18/18)
Expand Down Expand Up @@ -5007,10 +5007,6 @@ SUBROUTINE SFCMR_INIT( am_I_Root, Input_Opt )
LOCAL_NAME = TRIM('OCS')
CALL SFCMR_READ( am_I_Root, Input_Opt,
& LOCAL_NAME, ID, iMR=OCS_MR )
CASE ('H2')
LOCAL_NAME = TRIM('H2')
CALL SFCMR_READ( am_I_Root, Input_Opt,
& LOCAL_NAME, ID, iMR=H2_MR )
CASE DEFAULT
LOCAL_NAME = TRIM(TRAC_EMIT(I))
CALL SFCMR_READ( am_I_Root, Input_Opt, LOCAL_NAME, ID )
Expand Down Expand Up @@ -5639,7 +5635,7 @@ SUBROUTINE INIT_UCX( am_I_Root, Input_Opt, State_Chm, State_Diag,
id_ClOO = Ind_('ClOO' )
id_H1211 = Ind_('H1211' )
id_H1301 = Ind_('H1301' )
id_H2 = Ind_('H2', 'A') ! Need advected species index for H2
id_H2 = Ind_('H2' ) ! Need advected species index for H2
id_H2402 = Ind_('H2402' )
id_H2O = Ind_('H2O' )
id_HBr = Ind_('HBr' )
Expand Down Expand Up @@ -5897,11 +5893,11 @@ SUBROUTINE INIT_UCX( am_I_Root, Input_Opt, State_Chm, State_Diag,
TRAC_EMIT(N_EMIT) = 'CH3Br'
ENDIF

IF (id_H2.gt.0) THEN
N_EMIT = N_EMIT + 1
IDT_EMIT( N_EMIT) = id_H2
TRAC_EMIT(N_EMIT) = 'H2'
ENDIF
!IF (id_H2.gt.0) THEN
! N_EMIT = N_EMIT + 1
! IDT_EMIT( N_EMIT) = id_H2
! TRAC_EMIT(N_EMIT) = 'H2'
!ENDIF

IF (N_EMIT.gt.0) THEN
ALLOCATE( GRID_EMIT( State_Grid%NX, State_Grid%NY, N_EMIT ),
Expand Down
2 changes: 2 additions & 0 deletions HEMCO/Core/hco_config_mod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -939,6 +939,8 @@ SUBROUTINE Config_ReadCont( am_I_Root, HcoConfig, IU_HCO, &
Dta%MustFind = .TRUE.
ELSEIF ( TRIM(TmCycle) == "EC" ) THEN
Dta%CycleFlag = HCO_CFLAG_EXACT
Dta%MustFind = .TRUE.
Dta%UseSimYear= .TRUE.
ELSEIF ( TRIM(TmCycle) == "ECF" ) THEN
Dta%CycleFlag = HCO_CFLAG_EXACT
Dta%MustFind = .TRUE.
Expand Down
14 changes: 14 additions & 0 deletions HEMCO/Extensions/hcox_gfed_include_gfed3.H
Original file line number Diff line number Diff line change
Expand Up @@ -52,6 +52,9 @@
! - Updated the emission factors for NMVOCs and
! NOx from extratropical forests, savannas and
! agricultural fires from Akagi et al. (2011,ACP)
! 05 Feb 2020 - M.P Perez-Pena - Added H2
! 24 Mar 2022 - M.P Perez-Pena - Added biome-specific emission factors for
! H2 biomass burning emissions
!EOP
!------------------------------------------------------------------------------
!BOC
Expand Down Expand Up @@ -316,4 +319,15 @@ GFED3_EMFAC(29,3)=GFED3_EMFAC(2,3)
GFED3_EMFAC(29,4)=GFED3_EMFAC(2,4)
GFED3_EMFAC(29,5)=GFED3_EMFAC(2,5)
GFED3_EMFAC(29,6)=GFED3_EMFAC(2,6)

! H2
! H2 Added from Akagi (2011) and Andrae (2019)
GFED3_SPEC_NAME(30)="H2"
GFED3_EMFAC(30,1)=2.60E-03_hp ! AGR
GFED3_EMFAC(30,2)=3.10E-03_hp ! TROP. DEFOREST
GFED3_EMFAC(30,3)=2.10E-03_hp ! TEMPERATE
GFED3_EMFAC(30,4)=1.20E-03_hp ! PEAT
GFED3_EMFAC(30,5)=1.70E-03_hp ! SAV
GFED3_EMFAC(30,6)=1.60E-03_hp ! BOR

!EOC
13 changes: 13 additions & 0 deletions HEMCO/Extensions/hcox_gfed_include_gfed4.H
Original file line number Diff line number Diff line change
Expand Up @@ -45,6 +45,9 @@
! 29 Mar 2018 - K. Travis - Updated to latest GFED4 emission factors from
! https://www.geo.vu.nl/~gwerf/GFED/GFED4/ancill/
! - Added HAC, GLYX, Hg0
! 05 Feb 2020 - M.P Perez-Pena - Added H2
! 24 Mar 2022 - M.P Perez-Pena - Added biome-specific emission factors for
! H2 biomass burning emissions
!EOP
!------------------------------------------------------------------------------
!BOC
Expand Down Expand Up @@ -354,4 +357,14 @@ Inst%GFED4_EMFAC(34,3)=9.95E-04_hp
Inst%GFED4_EMFAC(34,4)=9.95E-04_hp
Inst%GFED4_EMFAC(34,5)=9.95E-04_hp
Inst%GFED4_EMFAC(34,6)=9.95E-04_hp

! H2 Added from Akagi (2011) and Andrae (2019)
GFED4_SPEC_NAME(35)="H2"
Inst%GFED4_EMFAC(35,1)=1.70E-03_hp ! SAVA
Inst%GFED4_EMFAC(35,2)=1.60E-03_hp ! BORF
Inst%GFED4_EMFAC(35,3)=2.10E-03_hp ! TEMP
Inst%GFED4_EMFAC(35,4)=3.10E-03_hp ! DEFO
Inst%GFED4_EMFAC(35,5)=1.20E-03_hp ! PEAT
Inst%GFED4_EMFAC(35,6)=2.60E-03_hp ! AGRI

!EOC