LCOV - code coverage report
Current view: | top level - hemco/HEMCO/src/Extensions - hcox_soilnox_mod.F90 (source / functions) |  |  | Hit | Total | Coverage
Test: | coverage.info |  | Lines: | 0 | 616 | 0.0 %
Date: | 2024-12-17 22:39:59 |  | Functions: | 0 | 22 | 0.0 %

              Line data    Source code
           1             : !------------------------------------------------------------------------------
           2             : !                   Harmonized Emissions Component (HEMCO)                    !
           3             : !------------------------------------------------------------------------------
           4             : !BOP
           5             : !
           6             : ! !MODULE: hcox_soilnox_mod.F90
           7             : !
           8             : ! !DESCRIPTION: Module HCOX\_SoilNOx\_Mod contains routines to compute soil
           9             : !  NOx emissions.  We follow the implementation in GEOS-Chem by Hudman
          10             : !  et al 2012.
          11             : !\\
          12             : !\\
          13             : ! !INTERFACE:
          14             : !
          15             : MODULE HCOX_SoilNOx_Mod
          16             : !
          17             : ! !USES:
          18             : !
          19             :   USE HCO_ERROR_Mod
          20             :   USE HCO_CHARTOOLS_MOD
          21             :   USE HCO_DIAGN_Mod
          22             :   USE HCOX_TOOLS_MOD
          23             :   USE HCOX_State_Mod,     ONLY : Ext_State
          24             :   USE HCO_STATE_Mod,      ONLY : HCO_State
          25             :
          26             :   IMPLICIT NONE
          27             :   PRIVATE
          28             : !
          29             : ! !PUBLIC MEMBER FUNCTIONS:
          30             : !
          31             :   PUBLIC :: HCOX_SoilNOx_Run
          32             :   PUBLIC :: HCOX_SoilNOx_Init
          33             :   PUBLIC :: HCOX_SoilNOx_Final
          34             : !
          35             : ! !REMARKS:
          36             : ! This is a HEMCO extension module that uses many of the HEMCO core
          37             : ! utilities.
          38             : !                                                                             .
          39             : !  Original codes from:
          40             : !    HARVARD ATMOSPHERIC CHEMISTRY MODELING GROUP
          41             : !    MODULE FOR SOIL NOX EMISSIONS
          42             : !    by Yuhang Wang, Gerry Gardner, and Prof. Daniel Jacob
          43             : !  Updated model code:
          44             : !    by  Rynda Hudman, Neil Moore, Randall Martin, and Bram Maasakkers
          45             : !                                                                             .
          46             : !  The soil NOx code has been updated from the original implementation
          47             : !  of Yienger & Levy [1995] from  Wang et al., [1998] as summarized below.
          48             : !                                                                             .
          49             : !  Old:
          50             : !  ENOx   = f( T, biome, w/d)  x Pulse(precip) x canopy uptake + FERT
          51             : !                                                                             .
          52             : !  New:
          53             : !  ENOx   = f( T, biome, WFPS, Fert)  x Pulse(dryspell) x canopy uptake
          54             : !                                                                             .
          55             : !  1 - Update moisture treatment: soil moisture as a continuous variable
          56             : !  using WFPS rather than discrete wet/dry states and purely exponential
          57             : !  T impact (impact = -1. Tg N/yr)
          58             : !                                                                             .
          59             : !  2 - Update to Fertilizer:  new fertilizer maps including chemical and
          60             : !  manure fertilizer from Potter et al., [2010] distributed using MODIS EVI
          61             : !  seasonality, online-N deposition as a fertilizer source, and N-fertilizer
          62             : !  source subject to T, WFPS, and pulsing like other N (impact = +1.3 Tg N/yr)
          63             : !                                                                             .
          64             : !  3- Update Pulsing Scheme: Yan et al., [2005] (shorter, stronger pulses)
          65             : !  (impact = +1. Tg N/yr). Also added restart file containing dry spell
          66             : !  information to properly account for dry spell length in continuing runs.
          67             : !                                                                             .
          68             : !  References:
          69             : !  ============================================================================
          70             : !  (1 ) Wang, Y., D.J. Jacob, and J.A. Logan,  Global simulation of
          71             : !        tropospheric O3-NOx-hydrocarbon chemistry, 1. Model formulation,
          72             : !        J. Geophys. Res., 103/D9, 10, 713-10,726, 1998.
          73             : !  (2 ) Yienger, J.J, and H. Levy, Empirical model of global soil-biogenic
          74             : !        NOx emissions, J. Geophys. Res., 100, D6, 11,447-11464, June 20, 1995.
          75             : !  (3 ) Yan, X., T. Ohara, and H. Akimoto, Statistical modeling of global
          76             : !        soil NOx emissions, Global Biogeochem. Cycles, 19, GB3019,
          77             : !        doi:10.1029/2004GB002276, 2005.
          78             : !  (4 ) Potter, P., Ramankutty, N., Bennett, E., and Donner, S.:
          79             : !        Characterizing the Spatial Patterns of Global Fertilizer Application
          80             : !        and Manure Production, Earth Interactions, 14, 1-22,
          81             : !        10.1175/2009EI288.1, 2010.
          82             : !  (5 ) Moore, N.E., Improving global bottom-up biogenic soil NOx inventories,
          83             : !        Master's Thesis, Dalhousie University, 2007.
          84             : !  (6 ) Hudman, R.C., N.E. Moore, A.K. Mebust, R.V. Martin, A.R. Russell,
          85             : !        L.C. Valin, and R.C Cohen, Steps toward a mechanistic model of global
          86             : !        soil nitric oxide emissions: implementation and space
          87             : !        based-constraints, Atmos. Chem. Phys., 12, 7779-7795,
          88             : !        doi:10.5194/acp-12-7779-2012, 2012.
          89             : !
          90             : ! !REVISION HISTORY:
          91             : !  See https://github.com/geoschem/hemco for complete history
          92             : !EOP
          93             : !------------------------------------------------------------------------------
          94             : !BOC
          95             : !
          96             : ! !MODULE VARIABLES:
          97             : !
          98             :   ! Derived type to hold MODIS land types
          99             :   TYPE MODL
         100             :      REAL(hp), POINTER             :: VAL              (:,:)
         101             :   ENDTYPE MODL
         102             :
         103             :   TYPE :: MyInst
         104             :      INTEGER                        :: Instance
         105             :      INTEGER                        :: ExtNr          ! Extension number
         106             :      INTEGER                        :: IDTNO          ! NO tracer ID
         107             :      LOGICAL                        :: LFERTILIZERNOX ! Use fertilizer NOx?
         108             :      REAL(hp)                       :: FERT_SCALE     ! fertilizer scale factor
         109             :
         110             :      ! Dry period length (from restart)
         111             :      REAL(sp), ALLOCATABLE          :: DRYPERIOD(:,:  )
         112             :
         113             :      ! Pulse factors (from restart)
         114             :      REAL(sp), ALLOCATABLE          :: PFACTOR(:,:  )
         115             :      REAL(sp), ALLOCATABLE          :: GWET_PREV(:,:  )
         116             :
         117             :      ! Deposition reservoir (from restart)
         118             :      REAL(sp), ALLOCATABLE          :: DEP_RESERVOIR(:,:  )
         119             :
         120             :      ! NOx in the canopy
         121             :      REAL(hp), ALLOCATABLE          :: CANOPYNOX(:,:,:)
         122             :
         123             :      ! MODIS landtype
         124             :      TYPE(MODL), POINTER            :: LANDTYPE(:) => NULL()
         125             :
         126             :      ! Soil fertilizer (kg/m3)
         127             :      REAL(hp), POINTER              :: SOILFERT(:,:) => NULL()
         128             :
         129             :      ! Fraction of arid and non-arid land
         130             :      REAL(hp), POINTER              :: CLIMARID(:,:) => NULL()
         131             :      REAL(hp), POINTER              :: CLIMNARID(:,:) => NULL()
         132             :
         133             :      ! DRYCOEFF (if read from settings in configuration file)
         134             :      REAL(hp), POINTER              :: DRYCOEFF(:)
         135             :
         136             :      ! Overall scale factor to be applied to total soil NOx emissions. Must
         137             :      ! be defined in the HEMCO configuration file as extension attribute
         138             :      ! 'Scaling_NO'
         139             :      REAL(sp),          ALLOCATABLE :: SpcScalVal(:)
         140             :      CHARACTER(LEN=61), ALLOCATABLE :: SpcScalFldNme(:)
         141             :
         142             :      ! Diagnostics
         143             :      REAL(sp), POINTER              :: FertNO_Diag(:,:)
         144             :
         145             :      TYPE(MyInst), POINTER          :: NextInst => NULL()
         146             :   END TYPE MyInst
         147             :
         148             :   ! Pointer to all instances
         149             :   TYPE(MyInst), POINTER             :: AllInst => NULL()
         150             : !
         151             : ! !DEFINED PARAMETERS:
         152             : !
         153             :   ! # of MODIS/Koppen biome types
         154             :   INTEGER, PARAMETER               :: NBIOM = 24
         155             :
         156             :   ! Max. # of allowed drycoeff vars
         157             :   INTEGER, PARAMETER               :: MaxDryCoeff = 50
         158             :
         159             :   ! Canopy wind extinction coefficients
         160             :   ! (cf. Yienger & Levy [1995], Sec 5), now a function of the
         161             :   ! MODIS/KOPPEN biometype (J.D. Maasakkers)
         162             :   REAL(hp),  PARAMETER, PRIVATE :: SOILEXC(NBIOM) =                        (/&
         163             :         0.10_hp, 0.50_hp, 0.10_hp, 0.10_hp, 0.10_hp, 0.10_hp, 0.10_hp,       &
         164             :         0.10_hp, 0.10_hp, 1.00_hp, 1.00_hp, 1.00_hp, 1.00_hp, 2.00_hp,       &
         165             :         4.00_hp, 4.00_hp, 4.00_hp, 4.00_hp, 4.00_hp, 4.00_hp, 4.00_hp,       &
         166             :         2.00_hp, 0.10_hp, 2.00_hp                                          /)
         167             :
         168             :   ! Steinkamp and Lawrence, 2011 A values, wet biome coefficients
         169             :   ! for each of the 24 soil biomes [ng N/m2/s].
         170             :   REAL(hp),  PARAMETER, PRIVATE  :: A_BIOME(NBIOM) =                      (/&
         171             :         0.00_hp, 0.00_hp, 0.00_hp, 0.00_hp, 0.00_hp, 0.06_hp, 0.09_hp,      &
         172             :         0.09_hp, 0.01_hp, 0.84_hp, 0.84_hp, 0.24_hp, 0.42_hp, 0.62_hp,      &
         173             :         0.03_hp, 0.36_hp, 0.36_hp, 0.35_hp, 1.66_hp, 0.08_hp, 0.44_hp,      &
         174             :         0.57_hp, 0.57_hp, 0.57_hp                                         /)
         175             :
         176             :   ! "A" coefficients for converting surface temp to soil temp
         177             :   ! for each of the 24 soil biomes
         178             :   REAL(hp),  PARAMETER, PRIVATE :: SOILTA(NBIOM)  =                       (/&
         179             :         0.00_hp, 0.92_hp, 0.00_hp, 0.66_hp, 0.66_hp, 0.66_hp, 0.66_hp,      &
         180             :         0.66_hp, 0.66_hp, 0.66_hp, 0.66_hp, 0.66_hp, 0.66_hp, 0.66_hp,      &
         181             :         0.84_hp, 0.84_hp, 0.84_hp, 0.84_hp, 0.84_hp, 0.84_hp, 0.84_hp,      &
         182             :         1.03_hp, 1.03_hp, 1.03_hp                                         /)
         183             :
         184             :   ! "B" coefficients for converting surface temp to soil temp
         185             :   ! for each of the 24 soil biomes
         186             :   REAL(hp),  PARAMETER, PRIVATE :: SOILTB(NBIOM)  =                       (/&
         187             :         0.00_hp, 4.40_hp, 0.00_hp, 8.80_hp, 8.80_hp, 8.80_hp, 8.80_hp,      &
         188             :         8.80_hp, 8.80_hp, 8.80_hp, 8.80_hp, 8.80_hp, 8.80_hp, 8.80_hp,      &
         189             :         3.60_hp, 3.60_hp, 3.60_hp, 3.60_hp, 3.60_hp, 3.60_hp, 3.60_hp,      &
         190             :         2.90_hp, 2.90_hp, 2.90_hp                                         /)
         191             :
         192             :   ! MODIS/Koppen resistance values
         193             :   INTEGER, PARAMETER, PRIVATE :: SNIMODIS(NBIOM) =                        (/&
         194             :            1,    2,    3,    4,    5,    6,    7,    8,    9,   10,         &
         195             :           11,   12,   13,   14,   15,   16,   17,   18,   19,   20,         &
         196             :           21,   22,   23,   24                                            /)
         197             :
         198             :   INTEGER, PARAMETER, PRIVATE :: SNIRI(NBIOM) =                           (/&
         199             :         9999,  200, 9999, 9999, 9999, 9999,  200,  200,  200,  200,         &
         200             :          200,  200,  200,  200,  200,  200,  200,  400,  400,  200,         &
         201             :          200,  200, 9999,  200                                            /)
         202             :
         203             :   INTEGER, PARAMETER, PRIVATE :: SNIRLU(NBIOM) =                          (/&
         204             :         9999, 9000, 9999, 9999, 9999, 9999, 9000, 9000, 9000, 9000,         &
         205             :         9000, 9000, 9000, 9000, 9000, 1000, 9000, 9000, 9000, 9000,         &
         206             :         1000, 9000, 9999, 9000                                            /)
         207             :
         208             :   INTEGER, PARAMETER, PRIVATE :: SNIRAC(NBIOM) =                          (/&
         209             :            0,  300,    0,    0,    0,    0,  100,  100,  100,  100,         &
         210             :          100,  100,  100,  100, 2000, 2000, 2000, 2000, 2000, 2000,         &
         211             :         2000,  200,  100,  200                                            /)
         212             :
         213             :   INTEGER, PARAMETER, PRIVATE :: SNIRGSS(NBIOM) =                         (/&
         214             :            0,    0,  100, 1000,  100, 1000,  350,  350,  350,  350,         &
         215             :          350,  350,  350,  350,  500,  200,  500,  500,  500,  500,         &
         216             :          200,  150,  400,  150                                            /)
         217             :
         218             :   INTEGER, PARAMETER, PRIVATE :: SNIRGSO(NBIOM) =                         (/&
         219             :         2000, 1000, 3500,  400, 3500,  400,  200,  200,  200,  200,         &
         220             :          200,  200,  200,  200,  200,  200,  200,  200,  200,  200,         &
         221             :          200,  150,  300,  150                                            /)
         222             :
         223             :   INTEGER, PARAMETER, PRIVATE :: SNIRCLS(NBIOM) =                         (/&
         224             :        9999, 2500, 9999, 9999, 9999, 9999, 2000, 2000, 2000, 2000,          &
         225             :        2000, 2000, 2000, 2000, 2000, 9999, 2000, 2000, 2000, 2000,          &
         226             :        9999, 2000, 9999, 2000                                             /)
         227             :
         228             :   INTEGER, PARAMETER, PRIVATE :: SNIRCLO(NBIOM) =                         (/&
         229             :        9999, 1000, 1000, 9999, 1000, 9999, 1000, 1000, 1000, 1000,          &
         230             :        1000, 1000, 1000, 1000, 1000, 9999, 1000, 1000, 1000, 1000,          &
         231             :        9999, 1000, 9999, 1000                                             /)
         232             :
         233             :   INTEGER, PARAMETER, PRIVATE :: SNIVSMAX(NBIOM) =                        (/&
         234             :          10,  100,  100,   10,  100,   10,  100,  100,  100,  100,          &
         235             :         100,  100,  100,  100,  100,  100,  100,  100,  100,  100,          &
         236             :         100,  100,  100,  100                                             /)
         237             :
         238             : !  ! Conversion factor from kg NO to ng N
         239             : !  REAL(hp),  PARAMETER, PRIVATE :: kgNO_to_ngN = 4.666d11 !(14/30 * 1e12)
         240             :
         241             : CONTAINS
         242             : !EOC
         243             : !------------------------------------------------------------------------------
         244             : !                   Harmonized Emissions Component (HEMCO)                    !
         245             : !------------------------------------------------------------------------------
         246             : !BOP
         247             : !
         248             : ! !IROUTINE: HCOX_SoilNOx_Run
         249             : !
         250             : ! !DESCRIPTION: Subroutine HcoX\_SoilNox\_Run is the driver routine to
         251             : ! calculate ship NOx emissions for the current time step. Emissions in
         252             : ! [kg/m2/s] are added to the emissions array of the passed
         253             : !\\
         254             : !\\
         255             : ! !INTERFACE:
         256             : !
         257           0 :   SUBROUTINE HCOX_SoilNOx_Run( ExtState, HcoState, RC )
         258             : !
         259             : ! !USES:
         260             : !
         261             :     USE HCO_Types_Mod,      ONLY : DiagnCont
         262             :     USE HCO_CLOCK_MOD,      ONLY : HcoClock_First
         263             :     USE HCO_CLOCK_MOD,      ONLY : HcoClock_Rewind
         264             :     USE HCO_FLuxArr_Mod,    ONLY : HCO_EmisAdd
         265             :     USE HCO_EmisList_Mod,   ONLY : HCO_GetPtr
         266             :     USE HCO_Calc_Mod,       ONLY : HCO_EvalFld
         267             :     USE HCO_ExtList_Mod,    ONLY : GetExtOpt
         268             :     USE HCO_ExtList_Mod,    ONLY : HCO_GetOpt
         269             :     USE HCO_Restart_Mod,    ONLY : HCO_RestartGet
         270             :     USE HCO_Restart_Mod,    ONLY : HCO_RestartWrite
         271             : !
         272             : ! !INPUT PARAMETERS:
         273             : !
         274             :     TYPE(Ext_State), POINTER        :: ExtState    ! Module options
         275             :
         276             : !
         277             : ! !INPUT/OUTPUT PARAMETERS:
         278             : !
         279             :     TYPE(HCO_State), POINTER        :: HcoState   ! Output obj
         280             :     INTEGER,         INTENT(INOUT)  :: RC
         281             : !
         282             : ! !REVISION HISTORY:
         283             : !  05 Nov 2013 - C. Keller - Initial Version
         284             : !  See https://github.com/geoschem/hemco for complete history
         285             : !EOP
         286             : !------------------------------------------------------------------------------
         287             : !BOC
         288             : !
         289             : ! !LOCAL VARIABLES:
         290             : !
         291             :     INTEGER                  :: I, J, N
         292           0 :     REAL(hp), TARGET         :: FLUX_2D(HcoState%NX,HcoState%NY)
         293           0 :     REAL(hp), TARGET         :: Tmp2D  (HcoState%NX,HcoState%NY)
         294           0 :     REAL(sp)                 :: Def2D  (HcoState%NX,HcoState%NY)
         295             :     REAL(hp)                 :: FERTDIAG, DEP_FERT, SOILFRT
         296             :     REAL*4                   :: TSEMIS
         297             :     REAL(hp)                 :: UNITCONV, IJFLUX
         298           0 :     REAL(dp), ALLOCATABLE    :: VecDp(:)
         299             :     LOGICAL                  :: FIRST
         300             :     LOGICAL                  :: aIR, FOUND
         301             :     CHARACTER(LEN= 31)       :: DiagnName
         302             :     CHARACTER(LEN=255)       :: MSG, DMY, LOC
         303             :     TYPE(MyInst),    POINTER :: Inst
         304             :
         305             :     !=================================================================
         306             :     ! HCOX_SoilNOx_RUN begins here!
         307             :     !=================================================================
         308           0 :     LOC = 'HCOX_SoilNOx_RUN (HCOX_SOILNOX_MOD.F90)'
         309             :
         310             :     ! Enter
         311           0 :     CALL HCO_ENTER( HcoState%Config%Err, LOC, RC )
         312           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         313           0 :         CALL HCO_ERROR( 'ERROR 0', RC, THISLOC=LOC )
         314           0 :         RETURN
         315             :     ENDIF
         316             :
         317             :     ! Return if extension disabled
         318           0 :     IF ( ExtState%SoilNOx < 0 ) RETURN
         319             :
         320             :     ! Nullify
         321           0 :     Inst   => NULL()
         322             :
         323             :     ! Get Instance
         324           0 :     CALL InstGet ( ExtState%SoilNox, Inst, RC )
         325           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         326           0 :        WRITE(MSG,*) 'Cannot find soil NOx instance Nr. ', ExtState%SoilNOx
         327           0 :        CALL HCO_ERROR(MSG,RC)
         328           0 :        RETURN
         329             :     ENDIF
         330             :
         331             :     ! Conversion factor from ng N to kg NO
         332           0 :     UNITCONV = 1.0e-12_hp / 14.0_hp * HcoState%Spc(Inst%IDTNO)%MW_g
         333             :
         334             :     !-----------------------------------------------------------------
         335             :     ! On first call, set pointers to all arrays needed by SoilNOx
         336             :     !-----------------------------------------------------------------
         337           0 :     FIRST = HcoClock_First ( HcoState%Clock, .TRUE. )
         338             :
         339             :     !IF ( FIRST ) THEN
         340           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK1',  Inst%LANDTYPE(1)%VAL,  RC )
         341           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         342           0 :            CALL HCO_ERROR( 'ERROR 1', RC, THISLOC=LOC )
         343           0 :            RETURN
         344             :        ENDIF
         345           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK2',  Inst%LANDTYPE(2)%VAL,  RC )
         346           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         347           0 :            CALL HCO_ERROR( 'ERROR 2', RC, THISLOC=LOC )
         348           0 :            RETURN
         349             :        ENDIF
         350           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK3',  Inst%LANDTYPE(3)%VAL,  RC )
         351           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         352           0 :            CALL HCO_ERROR( 'ERROR 3', RC, THISLOC=LOC )
         353           0 :            RETURN
         354             :        ENDIF
         355           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK4',  Inst%LANDTYPE(4)%VAL,  RC )
         356           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         357           0 :            CALL HCO_ERROR( 'ERROR 4', RC, THISLOC=LOC )
         358           0 :            RETURN
         359             :        ENDIF
         360           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK5',  Inst%LANDTYPE(5)%VAL,  RC )
         361           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         362           0 :            CALL HCO_ERROR( 'ERROR 5', RC, THISLOC=LOC )
         363           0 :            RETURN
         364             :        ENDIF
         365           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK6',  Inst%LANDTYPE(6)%VAL,  RC )
         366           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         367           0 :            CALL HCO_ERROR( 'ERROR 6', RC, THISLOC=LOC )
         368           0 :            RETURN
         369             :        ENDIF
         370           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK7',  Inst%LANDTYPE(7)%VAL,  RC )
         371           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         372           0 :            CALL HCO_ERROR( 'ERROR 7', RC, THISLOC=LOC )
         373           0 :            RETURN
         374             :        ENDIF
         375           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK8',  Inst%LANDTYPE(8)%VAL,  RC )
         376           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         377           0 :            CALL HCO_ERROR( 'ERROR 8', RC, THISLOC=LOC )
         378           0 :            RETURN
         379             :        ENDIF
         380           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK9',  Inst%LANDTYPE(9)%VAL,  RC )
         381           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         382           0 :            CALL HCO_ERROR( 'ERROR 9', RC, THISLOC=LOC )
         383           0 :            RETURN
         384             :        ENDIF
         385           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK10', Inst%LANDTYPE(10)%VAL, RC )
         386           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         387           0 :            CALL HCO_ERROR( 'ERROR 10', RC, THISLOC=LOC )
         388           0 :            RETURN
         389             :        ENDIF
         390           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK11', Inst%LANDTYPE(11)%VAL, RC )
         391           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         392           0 :            CALL HCO_ERROR( 'ERROR 11', RC, THISLOC=LOC )
         393           0 :            RETURN
         394             :        ENDIF
         395           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK12', Inst%LANDTYPE(12)%VAL, RC )
         396           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         397           0 :            CALL HCO_ERROR( 'ERROR 12', RC, THISLOC=LOC )
         398           0 :            RETURN
         399             :        ENDIF
         400           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK13', Inst%LANDTYPE(13)%VAL, RC )
         401           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         402           0 :            CALL HCO_ERROR( 'ERROR 13', RC, THISLOC=LOC )
         403           0 :            RETURN
         404             :        ENDIF
         405           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK14', Inst%LANDTYPE(14)%VAL, RC )
         406           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         407           0 :            CALL HCO_ERROR( 'ERROR 14', RC, THISLOC=LOC )
         408           0 :            RETURN
         409             :        ENDIF
         410           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK15', Inst%LANDTYPE(15)%VAL, RC )
         411           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         412           0 :            CALL HCO_ERROR( 'ERROR 15', RC, THISLOC=LOC )
         413           0 :            RETURN
         414             :        ENDIF
         415           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK16', Inst%LANDTYPE(16)%VAL, RC )
         416           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         417           0 :            CALL HCO_ERROR( 'ERROR 16', RC, THISLOC=LOC )
         418           0 :            RETURN
         419             :        ENDIF
         420           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK17', Inst%LANDTYPE(17)%VAL, RC )
         421           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         422           0 :            CALL HCO_ERROR( 'ERROR 17', RC, THISLOC=LOC )
         423           0 :            RETURN
         424             :        ENDIF
         425           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK18', Inst%LANDTYPE(18)%VAL, RC )
         426           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         427           0 :            CALL HCO_ERROR( 'ERROR 18', RC, THISLOC=LOC )
         428           0 :            RETURN
         429             :        ENDIF
         430           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK19', Inst%LANDTYPE(19)%VAL, RC )
         431           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         432           0 :            CALL HCO_ERROR( 'ERROR 19', RC, THISLOC=LOC )
         433           0 :            RETURN
         434             :        ENDIF
         435           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK20', Inst%LANDTYPE(20)%VAL, RC )
         436           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         437           0 :            CALL HCO_ERROR( 'ERROR 20', RC, THISLOC=LOC )
         438           0 :            RETURN
         439             :        ENDIF
         440           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK21', Inst%LANDTYPE(21)%VAL, RC )
         441           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         442           0 :            CALL HCO_ERROR( 'ERROR 21', RC, THISLOC=LOC )
         443           0 :            RETURN
         444             :        ENDIF
         445           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK22', Inst%LANDTYPE(22)%VAL, RC )
         446           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         447           0 :            CALL HCO_ERROR( 'ERROR 22', RC, THISLOC=LOC )
         448           0 :            RETURN
         449             :        ENDIF
         450           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK23', Inst%LANDTYPE(23)%VAL, RC )
         451           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         452           0 :            CALL HCO_ERROR( 'ERROR 23', RC, THISLOC=LOC )
         453           0 :            RETURN
         454             :        ENDIF
         455           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_LANDK24', Inst%LANDTYPE(24)%VAL, RC )
         456           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         457           0 :            CALL HCO_ERROR( 'ERROR 24', RC, THISLOC=LOC )
         458           0 :            RETURN
         459             :        ENDIF
         460           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_FERT',    Inst%SOILFERT,         RC )
         461           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         462           0 :            CALL HCO_ERROR( 'ERROR 25', RC, THISLOC=LOC )
         463           0 :            RETURN
         464             :        ENDIF
         465           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_ARID',    Inst%CLIMARID,         RC )
         466           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         467           0 :            CALL HCO_ERROR( 'ERROR 26', RC, THISLOC=LOC )
         468           0 :            RETURN
         469             :        ENDIF
         470           0 :        CALL HCO_EvalFld( HcoState, 'SOILNOX_NONARID', Inst%CLIMNARID,        RC )
         471           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         472           0 :            CALL HCO_ERROR( 'ERROR 27', RC, THISLOC=LOC )
         473           0 :            RETURN
         474             :        ENDIF
         475             :
         476           0 :     IF ( FIRST ) THEN
         477             :        ! Check if ExtState variables DRYCOEFF is defined. Otherwise, try to
         478             :        ! read it from settings.
         479           0 :        IF ( .NOT. ASSOCIATED(ExtState%DRYCOEFF) ) THEN
         480             :           CALL GetExtOpt( HcoState%Config, Inst%ExtNr, 'DRYCOEFF', &
         481           0 :                            OptValChar=DMY, FOUND=FOUND, RC=RC )
         482           0 :           IF ( .NOT. FOUND ) THEN
         483           0 :              CALL HCO_ERROR( 'DRYCOEFF not defined', RC )
         484           0 :              RETURN
         485             :           ENDIF
         486           0 :           ALLOCATE(VecDp(MaxDryCoeff))
         487             :           CALL HCO_CharSplit( DMY, HCO_GetOpt(HcoState%Config%ExtList,'Separator'), &
         488           0 :                               HCO_GetOpt(HcoState%Config%ExtList,'Wildcard'), VecDp, N, RC )
         489           0 :           IF ( RC /= HCO_SUCCESS ) THEN
         490           0 :               CALL HCO_ERROR( 'ERROR 28', RC, THISLOC=LOC )
         491           0 :               RETURN
         492             :           ENDIF
         493           0 :           ALLOCATE(Inst%DRYCOEFF(N))
         494           0 :           Inst%DRYCOEFF(1:N) = VecDp(1:N)
         495           0 :           ExtState%DRYCOEFF => Inst%DRYCOEFF
         496           0 :           DEALLOCATE(VecDp)
         497             :        ENDIF
         498             :     ENDIF
         499             :
         500             :     !---------------------------------------------------------------
         501             :     ! Fill restart variables. Restart variables can be specified
         502             :     ! in the HEMCO configuration file. In an ESMF environment, they
         503             :     ! can also be defined as internal state fields. The internal
         504             :     ! state fields take precedence over fields read through the
         505             :     ! HEMCO interface.
         506             :     ! The restart variables must be read on the first time step or
         507             :     ! after the rewinding the clock.
         508             :     !---------------------------------------------------------------
         509             :
         510             : #if !defined ( MODEL_GEOS )
         511           0 :     IF ( FIRST .OR. HcoClock_Rewind( HcoState%Clock, .TRUE. ) ) THEN
         512             : #else
         513             :     IF ( .TRUE. ) THEN
         514             : #endif
         515             :
         516             :        ! DEP_RESERVOIR. Read in kg NO/m3
         517             :        CALL HCO_EvalFld( HcoState, 'DEP_RESERVOIR_DEFAULT', &
         518           0 :                          Tmp2D, RC, FOUND=FOUND )
         519           0 :        IF ( FOUND ) THEN
         520           0 :           Def2D = Tmp2D
         521             :        ELSE
         522           0 :           Def2D = 1.0e-4_sp
         523             :        ENDIF
         524             :        CALL HCO_RestartGet( HcoState, 'DEP_RESERVOIR', &
         525           0 :                             Inst%DEP_RESERVOIR, RC, Def2D=Def2D )
         526           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         527           0 :            CALL HCO_ERROR( 'ERROR 29', RC, THISLOC=LOC )
         528           0 :            RETURN
         529             :        ENDIF
         530             :
         531             :        ! GWET_PREV [unitless]
         532             :        CALL HCO_RestartGet( HcoState, 'GWET_PREV', &
         533           0 :                             Inst%GWET_PREV, RC, FILLED=FOUND )
         534           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         535           0 :            CALL HCO_ERROR( 'ERROR 30', RC, THISLOC=LOC )
         536           0 :            RETURN
         537             :        ENDIF
         538           0 :        IF ( .NOT. FOUND ) THEN
         539           0 :           Inst%GWET_PREV = 0.0_sp
         540           0 :           IF ( HcoState%amIRoot ) THEN
         541           0 :              MSG = 'Cannot find GWET_PREV restart variable - initialized to 0.0!'
         542           0 :              CALL HCO_WARNING(HcoState%Config%Err,MSG,RC)
         543             :           ENDIF
         544             :        ENDIF
         545             :
         546             :        ! PFACTOR [unitless]
         547             :        CALL HCO_RestartGet( HcoState, 'PFACTOR', &
         548           0 :                             Inst%PFACTOR, RC, FILLED=FOUND )
         549           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         550           0 :            CALL HCO_ERROR( 'ERROR 31', RC, THISLOC=LOC )
         551           0 :            RETURN
         552             :        ENDIF
         553           0 :        IF ( .NOT. FOUND ) THEN
         554           0 :           Inst%PFACTOR = 1.0_sp
         555           0 :           IF ( HcoState%amIRoot ) THEN
         556           0 :              MSG = 'Cannot find PFACTOR restart variable - initialized to 1.0!'
         557           0 :              CALL HCO_WARNING(HcoState%Config%Err,MSG,RC)
         558             :           ENDIF
         559             :        ENDIF
         560             :
         561             :        ! DRYPERIOD [unitless]
         562             :        CALL HCO_RestartGet( HcoState, 'DRYPERIOD', &
         563           0 :                             Inst%DRYPERIOD, RC, FILLED=FOUND )
         564           0 :        IF ( RC /= HCO_SUCCESS ) THEN
         565           0 :            CALL HCO_ERROR( 'ERROR 32', RC, THISLOC=LOC )
         566           0 :            RETURN
         567             :        ENDIF
         568           0 :        IF ( .NOT. FOUND ) THEN
         569           0 :           Inst%DRYPERIOD = 0.0_sp
         570           0 :           IF ( HcoState%amIRoot ) THEN
         571           0 :              MSG = 'Cannot find DRYPERIOD restart variable - initialized to 0.0!'
         572           0 :              CALL HCO_WARNING(HcoState%Config%Err,MSG,RC)
         573             :           ENDIF
         574             :        ENDIF
         575             :
         576             :     ENDIF
         577             :
         578             :     !---------------------------------------------------------------
         579             :     ! Now need to call GET_CANOPY_NOX to break ugly dependency between
         580             :     ! drydep and soil NOx emissions. (bmy, 6/22/09)
         581             :     ! Now a function of the new MODIS/Koppen biome map (J.D. Maasakkers)
         582           0 :     CALL GET_CANOPY_NOX( HcoState, ExtState, Inst, RC )
         583           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         584           0 :         CALL HCO_ERROR( 'ERROR 33', RC, THISLOC=LOC )
         585           0 :         RETURN
         586             :     ENDIF
         587             :
         588             :     ! Init
         589           0 :     TSEMIS  = HcoState%TS_EMIS
         590           0 :     FLUX_2D = 0e+0_hp
         591             :
         592             :     ! Loop over each land grid-box, removed loop over landpoints
         593             : !$OMP PARALLEL DO                                                            &
         594             : !$OMP DEFAULT( SHARED )                                                      &
         595             : !$OMP PRIVATE( I, J, Dep_Fert, SoilFrt, FertDiag, IJflux                    )&
         596             : !$OMP COLLAPSE( 2                                                           )
         597           0 :     DO J = 1, HcoState%NY
         598           0 :     DO I = 1, HcoState%NX
         599             :
         600             :        ! Zero private variables for safety's sake
         601           0 :        Dep_Fert = 0.0_hp
         602           0 :        SoilFrt  = 0.0_hp
         603           0 :        FertDiag = 0.0_hp
         604           0 :        IJflux   = 0.0_hp
         605             :
         606             :        ! Do not calculate soil NOx emissions if there is no soil in
         607             :        ! the gridbox
         608           0 :        IF ( Inst%LANDTYPE(1)%VAL(I,J) == 1.0_hp ) CYCLE
         609             :
         610             :        ! Get Deposited Fertilizer DEP_FERT [kg NO/m2]
         611           0 :        CALL Get_Dep_N( I, J, ExtState, HcoState, Inst, Dep_Fert )
         612             :
         613             :        ! Get N fertilizer reservoir associated with chemical and
         614             :        ! manure fertilizer [kg NO/m2]
         615           0 :        IF ( Inst%LFERTILIZERNOX ) THEN
         616           0 :           SOILFRT = Inst%SOILFERT( I, J )
         617             :        ENDIF
         618             :
         619             :        ! Put in constraint if dry period gt 1 yr, keep at 1yr to
         620             :        ! avoid unrealistic pulse
         621           0 :        IF ( Inst%DRYPERIOD(I,J) > 8760.0_sp ) THEN
         622           0 :           Inst%DRYPERIOD(I,J) = 8760.0_sp
         623             :        ENDIF
         624             :
         625             :        ! Return NO emissions from soils [kg NO/m2/s]
         626             :        CALL Soil_NOx_Emission( ExtState,            Inst,                    &
         627             :                                TsEmis,              I,                       &
         628             :                                J,                   SoilFrt,                 &
         629           0 :                                Inst%GWET_PREV(I,J), Inst%DryPeriod(I,J),     &
         630           0 :                                Inst%PFACTOR(I,J),   IJflux,                  &
         631             :                                Dep_Fert,            FertDiag,                &
         632           0 :                                UnitConv,            Inst%CanopyNOx(I,J,:)   )
         633             :
         634             :        ! Write out
         635           0 :        Flux_2D(I,J) = MAX( IJflux, 0.0_hp )
         636             :
         637             :        ! Update diagnostics
         638           0 :        Inst%FertNO_Diag(I,J) = FertDiag
         639             :
         640             :     ENDDO !J
         641             :     ENDDO !I
         642             : !$OMP END PARALLEL DO
         643             :
         644             :     !-----------------------------------------------------------------
         645             :     ! EVENTUALLY ADD SCALE FACTORS
         646             :     !-----------------------------------------------------------------
         647             :
         648             :     ! Eventually apply species specific scale factor
         649           0 :     IF ( Inst%SpcScalVal(1) /= 1.0_sp ) THEN
         650           0 :        FLUX_2D = FLUX_2D * Inst%SpcScalVal(1)
         651             :     ENDIF
         652             :
         653             :     ! Eventually apply spatiotemporal scale factors
         654           0 :     CALL HCOX_SCALE ( HcoState, FLUX_2D, TRIM(Inst%SpcScalFldNme(1)), RC )
         655           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         656           0 :         CALL HCO_ERROR( 'ERROR 34', RC, THISLOC=LOC )
         657           0 :         RETURN
         658             :     ENDIF
         659             :
         660             :     !-----------------------------------------------------------------
         661             :     ! PASS TO HEMCO STATE AND UPDATE DIAGNOSTICS
         662             :     !-----------------------------------------------------------------
         663             :
         664             :     ! Add flux to emission array
         665             :     CALL HCO_EmisAdd( HcoState, FLUX_2D, Inst%IDTNO, &
         666           0 :                       RC,       ExtNr=Inst%ExtNr )
         667           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         668           0 :        CALL HCO_ERROR( 'HCO_EmisAdd error', RC )
         669           0 :        RETURN
         670             :     ENDIF
         671             :
         672             :     ! ----------------------------------------------------------------
         673             :     ! Eventually copy internal values to ESMF internal state object
         674             :     ! ----------------------------------------------------------------
         675             :
         676             :     ! DEP_RESERVOIR [kg/m3]
         677           0 :     CALL HCO_RestartWrite( HcoState, 'DEP_RESERVOIR', Inst%DEP_RESERVOIR, RC )
         678           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         679           0 :         CALL HCO_ERROR( 'ERROR 35', RC, THISLOC=LOC )
         680           0 :         RETURN
         681             :     ENDIF
         682             :
         683             :     ! GWET_PREV [unitless]
         684           0 :     CALL HCO_RestartWrite( HcoState, 'GWET_PREV', Inst%GWET_PREV, RC )
         685           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         686           0 :         CALL HCO_ERROR( 'ERROR 36', RC, THISLOC=LOC )
         687           0 :         RETURN
         688             :     ENDIF
         689             :
         690             :     ! PFACTOR [unitless]
         691           0 :     CALL HCO_RestartWrite( HcoState, 'PFACTOR', Inst%PFACTOR, RC )
         692           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         693           0 :         CALL HCO_ERROR( 'ERROR 37', RC, THISLOC=LOC )
         694           0 :         RETURN
         695             :     ENDIF
         696             :
         697             :     ! DRYPERIOD [unitless]
         698           0 :     CALL HCO_RestartWrite( HcoState, 'DRYPERIOD', Inst%DRYPERIOD, RC )
         699           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         700           0 :         CALL HCO_ERROR( 'ERROR 38', RC, THISLOC=LOC )
         701           0 :         RETURN
         702             :     ENDIF
         703             :
         704             :     ! Leave w/ success
         705           0 :     Inst => NULL()
         706           0 :     CALL HCO_LEAVE( HcoState%Config%Err,RC )
         707             :
         708           0 :   END SUBROUTINE HCOX_SoilNox_Run
         709             : !EOC
         710             : !------------------------------------------------------------------------------
         711             : !                   Harmonized Emissions Component (HEMCO)                    !
         712             : !------------------------------------------------------------------------------
         713             : !BOP
         714             : !
         715             : ! !IROUTINE: HCOX_SoilNOx_Init
         716             : !
         717             : ! !DESCRIPTION: Subroutine HcoX\_SoilNox\_Init initializes the HEMCO
         718             : ! SOILNOX extension.
         719             : !\\
         720             : !\\
         721             : ! !INTERFACE:
         722             : !
         723           0 :   SUBROUTINE HCOX_SoilNOx_Init( HcoState, ExtName, ExtState, RC )
         724             : !
         725             : ! !USES:
         726             : !
         727             :     USE HCO_ExtList_Mod,        ONLY : GetExtNr, GetExtOpt
         728             :     USE HCO_ExtList_Mod,        ONLY : GetExtSpcVal
         729             :     USE HCO_STATE_MOD,          ONLY : HCO_GetExtHcoID
         730             :     USE HCO_Restart_Mod,        ONLY : HCO_RestartDefine
         731             : !
         732             : ! !ARGUMENTS:
         733             : !
         734             :     TYPE(HCO_State),  POINTER        :: HcoState   ! Output obj
         735             :     CHARACTER(LEN=*), INTENT(IN   )  :: ExtName    ! Extension name
         736             :     TYPE(Ext_State),  POINTER        :: ExtState   ! Module options
         737             :     INTEGER,          INTENT(INOUT)  :: RC
         738             : !
         739             : ! !REMARKS:
         740             : !
         741             : ! !REVISION HISTORY:
         742             : !  05 Nov 2013 - C. Keller   - Initial Version
         743             : !  See https://github.com/geoschem/hemco for complete history
         744             : !EOP
         745             : !------------------------------------------------------------------------------
         746             : !BOC
         747             : !
         748             : ! !LOCAL VARIABLES:
         749             : !
         750             :     INTEGER                        :: ExtNr
         751             :     CHARACTER(LEN=255)             :: MSG, LOC
         752           0 :     CHARACTER(LEN=31), ALLOCATABLE :: SpcNames(:)
         753           0 :     INTEGER, ALLOCATABLE           :: HcoIDs(:)
         754             :     INTEGER                        :: nSpc, I, J, II, AS
         755             :     TYPE(MyInst), POINTER          :: Inst
         756             :
         757             :     !=================================================================
         758             :     ! HCOX_SoilNOx_INIT begins here!
         759             :     !=================================================================
         760           0 :     LOC = 'HCOX_SoilNOx_INIT (HCOX_SOILNOX_MOD.F90)'
         761             :
         762             :     ! Extension Nr.
         763           0 :     ExtNr = GetExtNr( HcoState%Config%ExtList, TRIM(ExtName) )
         764           0 :     IF ( ExtNr <= 0 ) RETURN
         765             :
         766             :     ! Enter
         767           0 :     CALL HCO_ENTER( HcoState%Config%Err, LOC, RC )
         768           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         769           0 :         CALL HCO_ERROR( 'ERROR 39', RC, THISLOC=LOC )
         770           0 :         RETURN
         771             :     ENDIF
         772             :
         773             :     ! Create instance
         774           0 :     Inst => NULL()
         775           0 :     CALL InstCreate ( ExtNr, ExtState%SoilNox, Inst, RC )
         776           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         777           0 :        CALL HCO_ERROR ( 'Cannot create soil NOx instance', RC )
         778           0 :        RETURN
         779             :     ENDIF
         780             :
         781             :     ! ----------------------------------------------------------------------
         782             :     ! Get species IDs and settings
         783             :     ! ----------------------------------------------------------------------
         784             :
         785             :     ! Read settings specified in configuration file
         786             :     ! Note: the specified strings have to match those in
         787             :     !       the config. file!
         788             :     CALL GetExtOpt( HcoState%Config, ExtNr, 'Use fertilizer NOx', &
         789           0 :                      OptValBool=Inst%LFERTILIZERNOX, RC=RC )
         790           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         791           0 :         CALL HCO_ERROR( 'ERROR 40', RC, THISLOC=LOC )
         792           0 :         RETURN
         793             :     ENDIF
         794             :
         795             :     ! Get global scale factor
         796           0 :     Inst%FERT_SCALE = HCOX_SoilNOx_GetFertScale()
         797             :
         798             :     ! Get HEMCO species IDs
         799           0 :     CALL HCO_GetExtHcoID( HcoState, ExtNr, HcoIDs, SpcNames, nSpc, RC )
         800           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         801           0 :         CALL HCO_ERROR( 'ERROR 41', RC, THISLOC=LOC )
         802           0 :         RETURN
         803             :     ENDIF
         804           0 :     IF ( nSpc /= 1 ) THEN
         805           0 :        MSG = 'Module soil NOx accepts only one species!'
         806           0 :        CALL HCO_ERROR(MSG, RC )
         807           0 :        RETURN
         808             :     ENDIF
         809           0 :     Inst%IDTNO = HcoIDs(1)
         810             :
         811             :     ! Get species scale factor
         812             :     CALL GetExtSpcVal( HcoState%Config, ExtNr, nSpc, &
         813           0 :                        SpcNames, 'Scaling', 1.0_sp, Inst%SpcScalVal, RC )
         814           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         815           0 :         CALL HCO_ERROR( 'ERROR 42', RC, THISLOC=LOC )
         816           0 :         RETURN
         817             :     ENDIF
         818             :
         819             :     CALL GetExtSpcVal( HcoState%Config, ExtNr, nSpc, &
         820           0 :                        SpcNames, 'ScaleField', HCOX_NOSCALE, Inst%SpcScalFldNme, RC )
         821           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         822           0 :         CALL HCO_ERROR( 'ERROR 43', RC, THISLOC=LOC )
         823           0 :         RETURN
         824             :     ENDIF
         825             :
         826             :     ! Verbose mode
         827           0 :     IF ( HcoState%amIRoot ) THEN
         828             :
         829             :        ! Write the name of the extension regardless of the verbose setting
         830           0 :        msg = 'Using HEMCO extension: SoilNOx (soil NOx emissions)'
         831           0 :        IF ( HCO_IsVerb( HcoState%Config%Err ) ) THEN
         832           0 :           CALL HCO_Msg( HcoState%Config%Err, msg, sep1='-' ) ! with separator
         833             :        ELSE
         834           0 :           CALL HCO_Msg( msg, verb=.TRUE.                   ) ! w/o separator
         835             :        ENDIF
         836             :
         837             :        ! Write all other messages as debug printout only
         838           0 :        WRITE(MSG,*) '   - NOx species            : ', TRIM(SpcNames(1)), Inst%IDTNO
         839           0 :        CALL HCO_MSG(HcoState%Config%Err,MSG)
         840           0 :        WRITE(MSG,*) '   - NOx scale factor       : ', Inst%SpcScalVal(1)
         841           0 :        CALL HCO_MSG(HcoState%Config%Err,MSG)
         842           0 :        WRITE(MSG,*) '   - NOx scale field        : ', TRIM(Inst%SpcScalFldNme(1))
         843           0 :        CALL HCO_MSG(HcoState%Config%Err,MSG)
         844           0 :        WRITE(MSG,*) '   - Use fertilizer NOx     : ', Inst%LFERTILIZERNOX
         845           0 :        CALL HCO_MSG(HcoState%Config%Err,MSG)
         846           0 :        WRITE(MSG,*) '   - Fertilizer scale factor: ', Inst%FERT_SCALE
         847           0 :        CALL HCO_MSG(HcoState%Config%Err,MSG,SEP2='-')
         848             :     ENDIF
         849             :
         850             :     ! ----------------------------------------------------------------------
         851             :     ! Set module variables
         852             :     ! ----------------------------------------------------------------------
         853             :
         854             :     ! horizontal dimensions
         855           0 :     I = HcoState%NX
         856           0 :     J = HcoState%NY
         857             :
         858           0 :     ALLOCATE( Inst%FertNO_Diag( I, J ), STAT=AS )
         859           0 :     IF ( AS /= 0 ) THEN
         860           0 :        CALL HCO_ERROR( 'FertNO_Diag', RC )
         861           0 :        RETURN
         862             :     ENDIF
         863           0 :     Inst%FertNO_Diag = 0.0_sp
         864             :
         865           0 :     ALLOCATE( Inst%DRYPERIOD( I, J ), STAT=AS )
         866           0 :     IF ( AS /= 0 ) THEN
         867           0 :        CALL HCO_ERROR( 'DRYPERIOD', RC )
         868           0 :        RETURN
         869             :     ENDIF
         870           0 :     Inst%DRYPERIOD     = 0.0_sp
         871             :
         872           0 :     ALLOCATE( Inst%PFACTOR( I, J ), STAT=AS )
         873           0 :     IF ( AS /= 0 ) THEN
         874           0 :        CALL HCO_ERROR( 'PFACTOR', RC )
         875           0 :        RETURN
         876             :     ENDIF
         877           0 :     Inst%PFACTOR       = 0.0_sp
         878             :
         879           0 :     ALLOCATE( Inst%GWET_PREV( I, J ), STAT=AS )
         880           0 :     IF ( AS /= 0 ) THEN
         881           0 :        CALL HCO_ERROR( 'GWET_PREV', RC )
         882           0 :        RETURN
         883             :     ENDIF
         884           0 :     Inst%GWET_PREV     = 0.0_sp
         885             :
         886           0 :     ALLOCATE( Inst%DEP_RESERVOIR( I, J ), STAT=AS )
         887           0 :     IF ( AS /= 0 ) THEN
         888           0 :        CALL HCO_ERROR( 'DEP_RESERVOIR', RC )
         889           0 :        RETURN
         890             :     ENDIF
         891           0 :     Inst%DEP_RESERVOIR = 0.0_sp
         892             :
         893           0 :     ALLOCATE( Inst%CANOPYNOX( I, J, NBIOM ), STAT=AS )
         894           0 :     IF ( AS /= 0 ) THEN
         895           0 :        CALL HCO_ERROR( 'CANOPYNOX', RC )
         896           0 :        RETURN
         897             :     ENDIF
         898           0 :     Inst%CANOPYNOX     = 0e+0_hp
         899             :
         900             :     ! Reserve 24 pointers for land fractions for each Koppen category
         901           0 :     ALLOCATE ( Inst%LANDTYPE(NBIOM), STAT=AS )
         902             :     IF ( AS /= 0 ) THEN
         903           0 :        CALL HCO_ERROR( 'LANDTYPE', RC )
         904           0 :        RETURN
         905             :     ENDIF
         906           0 :     DO II = 1,NBIOM
         907           0 :        ALLOCATE( Inst%LANDTYPE(II)%VAL( I, J ), STAT=AS )
         908             :        IF ( AS /= 0 ) THEN
         909           0 :           CALL HCO_ERROR( 'LANDTYPE array', RC )
         910           0 :           RETURN
         911             :        ENDIF
         912           0 :        Inst%LANDTYPE(II)%Val = 0.0_hp
         913             :     ENDDO
         914             :
         915             :     ALLOCATE ( Inst%SOILFERT  ( I, J ), &
         916             :                Inst%CLIMARID  ( I, J ), &
         917           0 :                Inst%CLIMNARID ( I, J ), STAT=AS )
         918           0 :     IF ( AS /= 0 ) THEN
         919           0 :        CALL HCO_ERROR( 'SOILFERT', RC )
         920           0 :        RETURN
         921             :     ENDIF
         922           0 :     Inst%SOILFERT  = 0.0_hp
         923           0 :     Inst%CLIMARID  = 0.0_hp
         924           0 :     Inst%CLIMNARID = 0.0_hp
         925             :
         926             :     ! ----------------------------------------------------------------------
         927             :     ! Set diagnostics
         928             :     ! ----------------------------------------------------------------------
         929             :     CALL Diagn_Create( HcoState  = HcoState,              &
         930             :                        cName     = 'EmisNO_Fert',         &
         931             :                        ExtNr     = ExtNr,                 &
         932             :                        Cat       = -1,                    &
         933             :                        Hier      = -1,                    &
         934             :                        HcoID     = -1,                    &
         935             :                        SpaceDim  = 2,                     &
         936             :                        OutUnit   = 'kg/m2/s',             &
         937             :                        AutoFill  = 0,                     &
         938             :                        Trgt2D    = Inst%FertNO_Diag,      &
         939           0 :                        RC        = RC )
         940           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         941           0 :         CALL HCO_ERROR( 'ERROR 44', RC, THISLOC=LOC )
         942           0 :         RETURN
         943             :     ENDIF
         944             :
         945             :     CALL HCO_RestartDefine( HcoState, 'PFACTOR', &
         946           0 :                             Inst%PFACTOR, '1',  RC )
         947           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         948           0 :         CALL HCO_ERROR( 'ERROR 45', RC, THISLOC=LOC )
         949           0 :         RETURN
         950             :     ENDIF
         951             :
         952             :     CALL HCO_RestartDefine( HcoState, 'DRYPERIOD', &
         953           0 :                             Inst%DRYPERIOD, '1',  RC )
         954           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         955           0 :         CALL HCO_ERROR( 'ERROR 46', RC, THISLOC=LOC )
         956           0 :         RETURN
         957             :     ENDIF
         958             :
         959             :     CALL HCO_RestartDefine( HcoState, 'GWET_PREV', &
         960           0 :                             Inst%GWET_PREV, '1',  RC )
         961           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         962           0 :         CALL HCO_ERROR( 'ERROR 47', RC, THISLOC=LOC )
         963           0 :         RETURN
         964             :     ENDIF
         965             :
         966             :     CALL HCO_RestartDefine( HcoState, '   DEP_RESERVOIR', &
         967           0 :                             Inst%DEP_RESERVOIR, 'kg/m3', RC )
         968           0 :     IF ( RC /= HCO_SUCCESS ) THEN
         969           0 :         CALL HCO_ERROR( 'ERROR 48', RC, THISLOC=LOC )
         970           0 :         RETURN
         971             :     ENDIF
         972             :
         973             :     ! ----------------------------------------------------------------------
         974             :     ! Set HEMCO extensions variables
         975             :     ! ----------------------------------------------------------------------
         976             :
         977             :     ! Activate required met fields
         978           0 :     ExtState%T2M%DoUse       = .TRUE.
         979           0 :     ExtState%GWETTOP%DoUse   = .TRUE.
         980           0 :     ExtState%SUNCOS%DoUse    = .TRUE.
         981           0 :     ExtState%U10M%DoUse      = .TRUE.
         982           0 :     ExtState%V10M%DoUse      = .TRUE.
         983           0 :     ExtState%LAI%DoUse       = .TRUE.
         984           0 :     ExtState%FRLAND%DoUse    = .TRUE.
         985           0 :     ExtState%FRLANDIC%DoUse  = .TRUE.
         986           0 :     ExtState%FROCEAN%DoUse   = .TRUE.
         987           0 :     ExtState%FRSEAICE%DoUse  = .TRUE.
         988           0 :     ExtState%FRLAKE%DoUse    = .TRUE.
         989           0 :     ExtState%SNODP%DoUse     = .TRUE.
         990           0 :     ExtState%RADSWG%DoUse    = .TRUE.
         991           0 :     ExtState%CLDFRC%DoUse    = .TRUE.
         992             :
         993             :     ! Activate required deposition parameter
         994           0 :     ExtState%DRY_TOTN%DoUse  = .TRUE.
         995           0 :     ExtState%WET_TOTN%DoUse  = .TRUE.
         996             :
         997             :     ! Leave w/ success
         998           0 :     Inst => NULL()
         999           0 :     IF ( ALLOCATED(HcoIDs  ) ) DEALLOCATE(HcoIDs  )
        1000           0 :     IF ( ALLOCATED(SpcNames) ) DEALLOCATE(SpcNames)
        1001           0 :     CALL HCO_LEAVE( HcoState%Config%Err,RC )
        1002             :
        1003           0 :   END SUBROUTINE HCOX_SoilNOx_Init
        1004             : !EOC
        1005             : !------------------------------------------------------------------------------
        1006             : !                   Harmonized Emissions Component (HEMCO)                    !
        1007             : !------------------------------------------------------------------------------
        1008             : !BOP
        1009             : !
        1010             : ! !IROUTINE: HCOX_SoilNOx_Final
        1011             : !
        1012             : ! !DESCRIPTION: Subroutine HcoX\_SoilNOx\_Final finalizes the HEMCO
        1013             : ! SOILNOX extension.
        1014             : !\\
        1015             : !\\
        1016             : ! !INTERFACE:
        1017             : !
        1018           0 :   SUBROUTINE HCOX_SoilNOx_Final( HcoState, ExtState, RC )
        1019             : !
        1020             : ! !USES
        1021             : !
        1022             : !
        1023             : ! !INPUT PARAMETERS:
        1024             : !
        1025             :     TYPE(HCO_State), POINTER        :: HcoState      ! HEMCO State obj
        1026             :     TYPE(Ext_State),  POINTER       :: ExtState      ! Extension state
        1027             : !
        1028             : ! !INPUT/OUTPUT PARAMETERS:
        1029             : !
        1030             :     INTEGER,         INTENT(INOUT)  :: RC
        1031             : !
        1032             : ! !REVISION HISTORY:
        1033             : !  05 Nov 2013 - C. Keller - Initial Version
        1034             : !  See https://github.com/geoschem/hemco for complete history
        1035             : !EOP
        1036             : !------------------------------------------------------------------------------
        1037             : !BOC
        1038             : !
        1039             : ! !LOCAL VARIABLES:
        1040             : !
        1041             :
        1042             :     !=================================================================
        1043             :     ! HCOX_SoilNOx_FINAL begins here!
        1044             :     !=================================================================
        1045           0 :     CALL InstRemove ( ExtState )
        1046             :
        1047           0 :   END SUBROUTINE HCOX_SoilNox_Final
        1048             : !EOC
        1049             : !------------------------------------------------------------------------------
        1050             : !                   Harmonized Emissions Component (HEMCO)                    !
        1051             : !------------------------------------------------------------------------------
        1052             : !BOP
        1053             : !
        1054             : ! !IROUTINE: HCOX_SoilNOx_GetFertScale
        1055             : !
        1056             : ! !DESCRIPTION: Function HCOX\_SoilNOx\_GETFERTSCALE returns the scale factor
        1057             : ! applied to fertilizer NOx emissions.
        1058             : !\\
        1059             : !\\
        1060             : ! !INTERFACE:
        1061             : !
        1062           0 :   FUNCTION HCOX_SoilNOx_GetFertScale() RESULT ( FERT_SCALE )
        1063             : !
        1064             : ! !ARGUMENTS:
        1065             : !
        1066             :     REAL(hp) :: FERT_SCALE
        1067             : !
        1068             : ! !REMARKS:
        1069             : !
        1070             : ! !REVISION HISTORY:
        1071             : !  11 Dec 2013 - C. Keller   - Initial version
        1072             : !  See https://github.com/geoschem/hemco for complete history
        1073             : !EOP
        1074             : !------------------------------------------------------------------------------
        1075             : !BOC
        1076             : !
        1077             : ! !LOCAL VARIABLES:
        1078             : !
        1079             :
        1080             :     ! Scale factor so that fertilizer emission = 1.8 Tg N/yr
        1081             :     ! (Stehfest and Bouwman, 2006)
        1082             :     ! before canopy reduction
        1083           0 :     FERT_SCALE = 0.0068_hp
        1084             :     ! Value calculated by running the 2x2.5 model
        1085             :     ! For now, use this value for all resolutions since regular soil NOx
        1086             :     ! emissions change with resolution as well (J.D. Maasakkers)
        1087             :
        1088           0 :   END FUNCTION HCOX_SoilNOx_GetFertScale
        1089             : !EOC
        1090             : !------------------------------------------------------------------------------
        1091             : !                   Harmonized Emissions Component (HEMCO)                    !
        1092             : !------------------------------------------------------------------------------
        1093             : !BOP
        1094             : !
        1095             : ! !IROUTINE: Soil_NOx_Emission
        1096             : !
        1097             : ! !DESCRIPTION: Subroutine Soil\_NOx\_Emission computes the emission of soil
        1098             : !  and fertilizer NOx for the GEOS-Chem model.
        1099             : !\\
        1100             : !\\
        1101             : ! !INTERFACE:
        1102             : !
        1103           0 :   SUBROUTINE Soil_NOx_Emission( ExtState,  Inst,     TS_EMIS,   I,           &
        1104             :                                 J,         SOILFRT,  GWET_PREV, DRYPERIOD,   &
        1105             :                                 PFACTOR,   SOILNOx,  DEPN,      FERTDIAG,    &
        1106           0 :                                 UNITCONV,  R_CANOPY                         )
        1107             : !
        1108             : ! !INPUT PARAMETERS:
        1109             : !
        1110             :     TYPE(Ext_State), POINTER     :: ExtState   ! Module options
        1111             :     TYPE(MyInst),    POINTER     :: Inst       ! Instance object
        1112             :     REAL*4,          INTENT(IN)  :: TS_EMIS    ! Emission timestep [s]
        1113             :     INTEGER,         INTENT(IN)  :: I          ! grid box lon index
        1114             :     INTEGER,         INTENT(IN)  :: J          ! grid box lat index
        1115             :     REAL(hp),        INTENT(IN)  :: DEPN       ! Dry Dep Fert term [kg/m2]
        1116             :     REAL(hp),        INTENT(IN)  :: SOILFRT    ! Fertilizer emissions [kg/m2]
        1117             :     REAL(hp),        INTENT(IN)  :: UNITCONV   ! ng N to kg NO
        1118             :     REAL(hp),        INTENT(IN)  :: R_CANOPY(:)! Resist of canopy to NOx [1/s]
        1119             : !
        1120             : ! !OUTPUT PARAMETERS:
        1121             : !
        1122             :     REAL(hp),        INTENT(OUT) :: SOILNOx    ! Soil NOx emissions [kg/m2/s]
        1123             :     REAL(sp),        INTENT(OUT) :: GWET_PREV  ! Soil Moisture Prev timestep
        1124             :     REAL(sp),        INTENT(OUT) :: DRYPERIOD  ! Dry period length in hours
        1125             :     REAL(sp),        INTENT(OUT) :: PFACTOR    ! Pulsing Factor
        1126             :     REAL(hp),        INTENT(OUT) :: FERTDIAG   ! Fert emissions [kg/m2/s]
        1127             : !
        1128             : ! !REMARKS:
        1129             : !  R_CANOPY is computed in routine GET_CANOPY_NOX of "canopy_nox_mod.f".
        1130             : !  This was originally in the GEOS-Chem dry deposition code, but was split
        1131             : !  off in order to avoid an ugly code dependency between the dry deposition
        1132             : !  and soil NOx codes.
        1133             : !
        1134             : !  As of v9-02, this module uses the MODIS/Koppen biome types instead
        1135             : !  of the Olson land type / biome type, making it different from the original
        1136             : !  dry deposition code (J.D. Maasakkers)
        1137             : !
        1138             : ! !REVISION HISTORY:
        1139             : !  31 Jan 2011 - R. Hudman   - New Model added
        1140             : !  See https://github.com/geoschem/hemco for complete history
        1141             : !EOP
        1142             : !------------------------------------------------------------------------------
        1143             : !BOC
        1144             : !
        1145             : ! !LOCAL VARIABLES:
        1146             : !
        1147             :     INTEGER   :: K
        1148             :     REAL(hp)  :: BASE_TERM, CRF_TERM,  PULSE
        1149             :     REAL(hp)  :: TC,        TEMP_TERM, WINDSQR
        1150             :     REAL(hp)  :: WET_TERM,  A_FERT,    A_BIOM
        1151             :     REAL(hp)  :: LAI,       SUNCOS,    GWET
        1152             :     REAL(hp)  :: ARID,      NARID
        1153             :
        1154             :     !=================================================================
        1155             :     ! Initialize
        1156             :     !=================================================================
        1157             :
        1158             :     ! Initialize
        1159           0 :     SOILNOX        = 0.0_hp
        1160           0 :     FERTDIAG       = 0.0_hp
        1161             :
        1162             :     ! Surface temperature [C]
        1163           0 :     TC             = ExtState%T2M%Arr%Val(I,J) - 273.15_hp
        1164             :
        1165             :     ! Surface wind speed, squared
        1166           0 :     WINDSQR        = ExtState%U10M%Arr%Val(I,J)**2 + &
        1167           0 :                      ExtState%V10M%Arr%Val(I,J)**2
        1168             :
        1169             :     ! Leaf area index
        1170           0 :     LAI = ExtState%LAI%Arr%Val(I,J)
        1171             :
        1172             :     ! Cosine of Solar Zenit Angle
        1173           0 :     SUNCOS = ExtState%SUNCOS%Arr%Val(I,J)
        1174             :
        1175             :     ! Top soil wetness [unitless]
        1176           0 :     GWET = ExtState%GWETTOP%Arr%Val(I,J)
        1177             :
        1178             :     !=================================================================
        1179             :     ! Compute soil NOx emissions
        1180             :     !=================================================================
        1181             :
        1182             :     ! Cumulative multiplication factor (over baseline emissions)
        1183             :     ! that accounts for soil pulsing
        1184           0 :     PULSE = PULSING( GWET, TS_EMIS, GWET_PREV, PFACTOR, DRYPERIOD )
        1185             :
        1186             :     ! ------Loop Over MODIS/Koppen  Landtypes
        1187           0 :     DO K = 1, 24
        1188             :
        1189             :        ! Temperature-dependent term of soil NOx emissions [unitless]
        1190             :        ! Use GWET instead of climo wet/dry
        1191           0 :        TEMP_TERM = SOILTEMP( K, TC, GWET )
        1192             :
        1193             :        ! Soil moisture scaling of soil NOx emissions
        1194           0 :        ARID      = Inst%CLIMARID(I,J)
        1195           0 :        NARID     = Inst%CLIMNARID(I,J)
        1196           0 :        WET_TERM  = SOILWET( GWET , ARID, NARID )
        1197             :
        1198             :        ! Fertilizer emission [kg/m2/s]
        1199           0 :        A_FERT    = FERTADD( SOILFRT, DEPN )
        1200             :
        1201             :        ! Scale fertilizer emissions as specified
        1202             :        ! (scale needed to force fert emiss of 1.8 Tg N/yr w/o canopy uptake)
        1203           0 :        A_FERT    = A_FERT * Inst%FERT_SCALE
        1204             :
        1205             :        ! Canopy reduction factor
        1206           0 :        CRF_TERM  = SOILCRF( K, LAI, R_CANOPY(K), WINDSQR, SUNCOS )
        1207             :
        1208             :        ! Base emission. ng N/m2/s --> kg NO/m2/s
        1209           0 :        A_BIOM    = A_BIOME(K) * UNITCONV
        1210             :
        1211             :        ! SOILNOX includes fertilizer
        1212             :        SOILNOX   = ( SOILNOX                                                 &
        1213             :                  + ( A_BIOM + A_FERT )                                       &
        1214             :                  * ( TEMP_TERM * WET_TERM * PULSE )                          &
        1215           0 :                  * Inst%LANDTYPE(K)%VAL(I,J)                                 &
        1216           0 :                  * ( 1.0_hp - CRF_TERM  )                                   )
        1217             :
        1218             :        ! FERTDIAG, only used for the fertilizer diagnostic
        1219             :        FERTDIAG  = ( FERTDIAG                                                &
        1220             :                  + ( A_FERT )                                                &
        1221             :                  * ( TEMP_TERM * WET_TERM * PULSE )                          &
        1222           0 :                  * Inst%LANDTYPE(K)%VAL(I,J)                                 &
        1223           0 :                  * ( 1.0_hp - CRF_TERM  )                                   )
        1224             :
        1225             :     ENDDO
        1226             :
        1227           0 :   END SUBROUTINE Soil_NOx_Emission
        1228             : !EOC
        1229             : !------------------------------------------------------------------------------
        1230             : !                   Harmonized Emissions Component (HEMCO)                    !
        1231             : !------------------------------------------------------------------------------
        1232             : !BOP
        1233             : !
        1234             : ! !IROUTINE: Get_Canopy_NOx
        1235             : !
        1236             : ! !DESCRIPTION: Subroutine Get\_Canopy\_NOx computes the bulk surface
        1237             : !  resistance of the canopy to NOx.  This computation was originally done
        1238             : !  within legacy routine DEPVEL (in "drydep\_mod.f").  Moving this computation
        1239             : !  to Get\_Canopy\_NOx now allows for a totally clean separation between
        1240             : !  dry deposition routines and emissions routines in GEOS-Chem.
        1241             : !\\
        1242             : !\\
        1243             : ! !INTERFACE:
        1244             : !
        1245           0 :   SUBROUTINE Get_Canopy_NOx( HcoState, ExtState, Inst, RC )
        1246             : !
        1247             : ! !USES:
        1248             : !
        1249             :     USE Drydep_Toolbox_Mod, ONLY : BIOFIT
        1250             :     USE HCO_GeoTools_Mod,   ONLY : HCO_LANDTYPE
        1251             : !
        1252             : ! !ARGUMENTS:
        1253             : !
        1254             :     TYPE(HCO_State), POINTER        :: HcoState
        1255             :     TYPE(Ext_State), POINTER        :: ExtState
        1256             :     TYPE(MyInst),    POINTER        :: Inst
        1257             :     INTEGER,         INTENT(INOUT)  :: RC
        1258             : !
        1259             : ! !REMARKS:
        1260             : !  For backwards compatibility, the bulk surface resistance is stored
        1261             : !  in common block array CANOPYNOX in "commsoil.h".  Leave it like this
        1262             : !  for the time being...we'll clean it up when we fix all of the soil
        1263             : !  NOx routines.
        1264             : !
        1265             : ! !REVISION HISTORY:
        1266             : !  14 Jun 2012 - J.D. Maasakkers - Rewritten as a function of the
        1267             : !                                     MODIS/Koppen biometype
        1268             : !  See https://github.com/geoschem/hemco for complete history
        1269             : !EOP
        1270             : !------------------------------------------------------------------------------
        1271             : !BOC
        1272             : !
        1273             : ! !DEFINED PARAMETERS:
        1274             : !
        1275             :     ! Molecular weight of water [kg]
        1276             :     REAL(hp), PARAMETER :: XMWH2O = 18e-3_hp
        1277             :
        1278             :     ! Surface pressure??? [Pa]
        1279             :     REAL(hp), PARAMETER :: PRESS  = 1.5e+5_hp
        1280             : !
        1281             : ! !LOCAL VARIABLES:
        1282             : !
        1283             :     ! Scalars
        1284             :     INTEGER             :: I,     J,     K,      KK
        1285             :     INTEGER             :: DCSZ
        1286             :     REAL(hp)            :: F0,    HSTAR, XMW
        1287             :     REAL(hp)            :: DTMP1, DTMP2, DTMP3,  DTMP4, GFACT, GFACI
        1288             :     REAL(hp)            :: RT,    RAD0,  RIX,    RIXX,  RDC,   RLUXX
        1289             :     REAL(hp)            :: RGSX,  RCLX,  TEMPK,  TEMPC
        1290             :     REAL(hp)            :: LAI,   SUNCOS, CLDFRC
        1291             :     REAL(hp)            :: BIO_RESULT
        1292             :
        1293             :     ! Arrays
        1294             :     REAL(hp)            :: RI  (NBIOM)
        1295             :     REAL(hp)            :: RLU (NBIOM)
        1296             :     REAL(hp)            :: RAC (NBIOM)
        1297             :     REAL(hp)            :: RGSS(NBIOM)
        1298             :     REAL(hp)            :: RGSO(NBIOM)
        1299             :     REAL(hp)            :: RCLS(NBIOM)
        1300             :     REAL(hp)            :: RCLO(NBIOM)
        1301             :
        1302             :     !=================================================================
        1303             :     ! GET_CANOPY_NOX begins here!
        1304             :     !=================================================================
        1305             :
        1306             :     ! Set physical parameters
        1307           0 :     HSTAR = 0.01e+0_hp              ! Henry's law constant
        1308           0 :     F0    = 0.1e+0_hp               ! Reactivity factor for biological oxidation
        1309           0 :     XMW   = 46e-3_hp               ! Molecular wt of NO2 (kg)
        1310             :
        1311             :     ! Get size of DRYCOEFF (will be passed to routine BIOFIT)
        1312           0 :     DCSZ = SIZE( ExtState%DRYCOEFF )
        1313             :
        1314             :     ! Loop over surface boxes
        1315           0 :     DO J = 1, HcoState%NY
        1316           0 :     DO I = 1, HcoState%NX
        1317             :
        1318             :        ! Surface temperature [K] and [C]
        1319           0 :        TEMPK = ExtState%T2M%Arr%Val(I,J)
        1320           0 :        TEMPC = ExtState%T2M%Arr%Val(I,J) - 273.15_hp
        1321             :
        1322             :        ! Compute bulk surface resistance for gases.
        1323             :        !
        1324             :        !  Adjust external surface resistances for temperature;
        1325             :        !  from Wesely [1989], expression given in text on p. 1296.
        1326           0 :        RT = 1000.0_hp * EXP( -TEMPC - 4.0_hp )
        1327             :
        1328             :        !--------------------------------------------------------------
        1329             :        ! Get surface resistances - loop over biome types K
        1330             :        !
        1331             :        ! The land types within each grid square are defined using the
        1332             :        ! Olson land-type database.  Each of the Olson land types is
        1333             :        ! assigned a corresponding "deposition land type" with
        1334             :        ! characteristic values of surface resistance components.
        1335             :        ! There are 74 Olson land-types but only 11 deposition
        1336             :        ! land-types (i.e., many of the Olson land types share the
        1337             :        ! same deposition characteristics).  Surface resistance
        1338             :        ! components for the "deposition land types" are from Wesely
        1339             :        ! [1989] except for tropical forests [Jacob and Wofsy, 1990]
        1340             :        ! and for tundra [Jacob et al., 1992].  All surface resistance
        1341             :        ! components are normalized to a leaf area index of unity.
        1342             :        !--------------------------------------------------------------
        1343             :        !Loop over all biometypes
        1344           0 :        DO K = 1, NBIOM
        1345             :
        1346             :           ! Skip if not present
        1347           0 :           IF ( Inst%LANDTYPE(K)%VAL(I,J) == 0.0_hp ) CYCLE
        1348             :
        1349             :           ! Set second loop variable to K to allow snow/ice correction
        1350           0 :           KK = K
        1351             :
        1352             :           ! If the surface is snow or ice, then set K=3
        1353           0 :           IF ( ( ExtState%SNODP%Arr%Val(I,J) > 0.2 ) .OR.     &
        1354           0 :                ( HCO_LANDTYPE(ExtState%FRLAND%Arr%Val(I,J),   &
        1355           0 :                               ExtState%FRLANDIC%Arr%Val(I,J), &
        1356           0 :                               ExtState%FROCEAN%Arr%Val(I,J),  &
        1357           0 :                               ExtState%FRSEAICE%Arr%Val(I,J), &
        1358           0 :                               ExtState%FRLAKE%Arr%Val(I,J) ) == 2 ) ) KK = 3
        1359             :
        1360             :           ! USE new MODIS/KOPPEN Biometypes to read data
        1361             :
        1362             :           ! Read the internal resistance RI (minimum stomatal resistance
        1363             :           ! for water vapor, per unit area of leaf) from the IRI array;
        1364             :           ! a '9999' value means no deposition to stomata so we impose a
        1365             :           ! very large value for RI.
        1366           0 :           RI(K) = DBLE( SNIRI(KK) )
        1367           0 :           IF ( RI(K) >= 9999.e+0_hp ) RI(K)= 1.e+12_hp
        1368             :
        1369             :           ! Cuticular resistances IRLU read in from 'drydep.table'
        1370             :           ! are per unit area of leaf; divide them by the leaf area index
        1371             :           ! to get a cuticular resistance for the bulk canopy.  If IRLU is
        1372             :           !'9999' it means there are no cuticular surfaces on which to
        1373             :           ! deposit so we impose a very large value for RLU.
        1374           0 :           IF ( SNIRLU(KK) >= 9999 .OR. &
        1375             :                ExtState%LAI%Arr%Val(I,J) <= 0e+0_hp ) THEN
        1376           0 :              RLU(K)  = 1.e+6_hp
        1377             :           ELSE
        1378           0 :              RLU(K)= DBLE(SNIRLU(KK)) / ExtState%LAI%Arr%Val(I,J) + RT
        1379             :           ENDIF
        1380             :
        1381             :           ! The following are the remaining resistances for the Wesely
        1382             :           ! resistance-in-series model for a surface canopy
        1383             :           ! (see Atmos. Environ. paper, Fig.1).
        1384           0 :           RAC(K)  = MAX( DBLE( SNIRAC(KK)  ),      1e+0_hp )
        1385           0 :           RGSS(K) = MAX( DBLE( SNIRGSS(KK) ) + RT, 1e+0_hp )
        1386           0 :           RGSO(K) = MAX( DBLE( SNIRGSO(KK) ) + RT, 1e+0_hp )
        1387           0 :           RCLS(K) =      DBLE( SNIRCLS(KK) ) + RT
        1388           0 :           RCLO(K) =      DBLE( SNIRCLO(KK) ) + RT
        1389             :
        1390           0 :           IF (  RAC(K) >= 9999.e+0_hp ) RAC(K)  = 1e+12_hp
        1391           0 :           IF ( RGSS(K) >= 9999.e+0_hp ) RGSS(K) = 1e+12_hp
        1392           0 :           IF ( RGSO(K) >= 9999.e+0_hp ) RGSO(K) = 1e+12_hp
        1393           0 :           IF ( RCLS(K) >= 9999.e+0_hp ) RCLS(K) = 1e+12_hp
        1394           0 :           IF ( RCLO(K) >= 9999.e+0_hp ) RCLO(K) = 1e+12_hp
        1395             :
        1396             :           !-------------------------------------------------------------
        1397             :           ! Adjust stomatal resistances for insolation and temperature:
        1398             :           !
        1399             :           ! Temperature adjustment is from Wesely [1989], equation (3).
        1400             :           !
        1401             :           ! Light adjustment by the function BIOFIT is described by Wang
        1402             :           ! [1996].  It combines:
        1403             :           !
        1404             :           ! - Local dependence of stomal resistance on the intensity I
        1405             :           !   of light impinging the leaf; this is expressed as a
        1406             :           !   multiplicative factor I/(I+b) to the stomatal resistance
        1407             :           !   where b = 50 W m-2
        1408             :           !   (equation (7) of Baldocchi et al. [1987])
        1409             :           ! - Radiative transfer of direct and diffuse radiation in the
        1410             :           !   canopy using equations (12)-(16) from Guenther et al.
        1411             :           !   [1995]
        1412             :           ! - Separate accounting of sunlit and shaded leaves using
        1413             :           !   equation (12) of Guenther et al. [1995]
        1414             :           ! - Partitioning of the radiation at the top of the canopy
        1415             :           !   into direct and diffuse components using a
        1416             :           !   parameterization to results from an atmospheric radiative
        1417             :           !   transfer model [Wang, 1996]
        1418             :           !
        1419             :           ! The dependent variables of the function BIOFIT are the leaf
        1420             :           ! area index (XYLAI), the cosine of zenith angle (SUNCOS) and
        1421             :           ! the fractional cloud cover (CFRAC).  The factor GFACI
        1422             :           ! integrates the light dependence over the canopy depth; so
        1423             :           ! be scaled by LAI to yield a bulk canopy value because that's
        1424             :           ! already done in the GFACI formulation.
        1425             :           !-------------------------------------------------------------
        1426             :
        1427             :           ! Radiation @ sfc [W/m2]
        1428           0 :           RAD0 = ExtState%RADSWG%Arr%Val(I,J)
        1429             :
        1430             :           ! Internal resistance
        1431           0 :           RIX  = RI(K)
        1432             :
        1433             :           ! Skip the following block if the resistance RIX is high
        1434           0 :           IF ( RIX < 9999e+0_hp ) THEN
        1435           0 :              GFACT = 100.0e+0_hp
        1436             :
        1437           0 :              IF ( TEMPC > 0.e+0_hp .AND. TEMPC < 40.e+0_hp) THEN
        1438           0 :                 GFACT = 400.e+0_hp / TEMPC / ( 40.0e+0_hp - TEMPC )
        1439             :              ENDIF
        1440             :
        1441           0 :              GFACI = 100.e+0_hp
        1442             :
        1443           0 :              IF ( RAD0 > 0e+0_hp .AND. ExtState%LAI%Arr%Val(I,J) > 0e+0_hp ) THEN
        1444             :
        1445           0 :                 LAI    = ExtState%LAI%Arr%Val(I,J)
        1446           0 :                 SUNCOS = ExtState%SUNCOS%Arr%Val(I,J)
        1447           0 :                 CLDFRC = ExtState%CLDFRC%Arr%Val(I,J)
        1448             :
        1449             :                 BIO_RESULT = BIOFIT( ExtState%DRYCOEFF, LAI, &
        1450           0 :                               SUNCOS, CLDFRC, SIZE(ExtState%DRYCOEFF) )
        1451             :
        1452           0 :                 GFACI= 1e+0_hp / BIO_RESULT
        1453             :              ENDIF
        1454             :
        1455           0 :              RIX = RIX * GFACT * GFACI
        1456             :           ENDIF
        1457             :
        1458             :           ! Compute aerodynamic resistance to lower elements in lower
        1459             :           ! part of the canopy or structure, assuming level terrain -
        1460             :           ! equation (5) of Wesely [1989].
        1461           0 :           RDC = 100.0_hp * ( 1.0_hp + 1000.0_hp / ( RAD0 + 10.0_hp ) )
        1462             :
        1463             :           ! Loop over species; species-dependent corrections to resistances
        1464             :           ! are from equations (6)-(9) of Wesely [1989].
        1465             :           !
        1466             :           ! NOTE: here we only consider NO2 (bmy, 6/22/09)
        1467             :           RIXX   = RIX * DIFFG( TEMPK, PRESS, XMWH2O ) / &
        1468             :                          DIFFG( TEMPK, PRESS, XMW    )   &
        1469           0 :                  + 1.e+0_hp / ( HSTAR/3000.e+0_hp + 100.e+0_hp*F0  )
        1470             :
        1471           0 :           RLUXX  = 1.e+12_hp
        1472             :
        1473           0 :           IF ( RLU(K) < 9999.e+0_hp ) THEN
        1474           0 :              RLUXX = RLU(K) / ( HSTAR / 1.0e+05_hp + F0 )
        1475             :           ENDIF
        1476             :
        1477             :           ! To prevent virtually zero resistance to species with huge HSTAR,
        1478             :           ! such as HNO3, a minimum value of RLUXX needs to be set.
        1479             :           ! The rationality of the existence of such a minimum is
        1480             :           ! demonstrated by the observed relationship between Vd(NOy-NOx)
        1481             :           ! and Ustar in Munger et al.[1996]; Vd(HNO3) never exceeds 2 cm/s
        1482             :           ! in observations. The corresponding minimum resistance is 50 s/m.
        1483             :           ! was introduced by J.Y. Liang on 7/9/95.
        1484           0 :           RGSX = 1e+0_hp / ( HSTAR/1e+5_hp/RGSS(K) + F0/RGSO(K) )
        1485           0 :           RCLX = 1e+0_hp / ( HSTAR/1e+5_hp/RCLS(K) + F0/RCLO(K) )
        1486             :
        1487             :           ! Get the bulk surface resistance of the canopy
        1488             :           ! from the network of resistances in parallel and in series
        1489             :           ! (Fig. 1 of Wesely [1989])
        1490           0 :           DTMP1 = 1.e+0_hp / RIXX
        1491           0 :           DTMP2 = 1.e+0_hp / RLUXX
        1492           0 :           DTMP3 = 1.e+0_hp / ( RAC(K) + RGSX )
        1493           0 :           DTMP4 = 1.e+0_hp / ( RDC      + RCLX )
        1494             :
        1495             :           ! Save the within canopy depvel of NOx, used in calculating
        1496             :           ! the canopy reduction factor for soil emissions [1/s]
        1497           0 :           Inst%CANOPYNOX(I,J,K) = DTMP1 + DTMP2 + DTMP3 + DTMP4
        1498             :
        1499             :        ENDDO !K
        1500             :     ENDDO !I
        1501             :     ENDDO !J
        1502             :
        1503             :     ! Return w/ success
        1504           0 :     RC = HCO_SUCCESS
        1505             :
        1506           0 :   END SUBROUTINE Get_Canopy_NOx
        1507             : !EOC
        1508             : !------------------------------------------------------------------------------
        1509             : !                   Harmonized Emissions Component (HEMCO)                    !
        1510             : !------------------------------------------------------------------------------
        1511             : !BOP
        1512             : !
        1513             : ! !IROUTINE: DiffG
        1514             : !
        1515             : ! !DESCRIPTION: Function DiffG calculates the molecular diffusivity [m2/s] in
        1516             : !  air for a gas X of molecular weight XM [kg] at temperature TK [K] and
        1517             : !  pressure PRESS [Pa].
        1518             : !\\
        1519             : !\\
        1520             : ! !INTERFACE:
        1521             : !
        1522           0 :   FUNCTION DiffG( TK, PRESS, XM ) RESULT( DIFF_G )
        1523             : !
        1524             : ! !INPUT PARAMETERS:
        1525             : !
        1526             :     REAL(hp), INTENT(IN) :: TK      ! Temperature [K]
        1527             :     REAL(hp), INTENT(IN) :: PRESS   ! Pressure [hPa]
        1528             :     REAL(hp), INTENT(IN) :: XM      ! Molecular weight of gas [kg]
        1529             : !
        1530             : ! !RETURN VALUE:
        1531             : !
        1532             :     REAL(hp)             :: DIFF_G  ! Molecular diffusivity [m2/s]
        1533             : !
        1534             : ! !REMARKS:
        1535             : !  We specify the molecular weight of air (XMAIR) and the hard-sphere molecular
        1536             : !  radii of air (RADAIR) and of the diffusing gas (RADX).  The molecular
        1537             : !  radius of air is given in a Table on p. 479 of Levine [1988].  The Table
        1538             : !  also gives radii for some other molecules.  Rather than requesting the user
        1539             : !  to supply a molecular radius we specify here a generic value of 2.E-10 m for
        1540             : !  all molecules, which is good enough in terms of calculating the diffusivity
        1541             : !  as long as molecule is not too big.
        1542             : !
        1543             : ! !REVISION HISTORY:
        1544             : !  22 Jun 2009 - R. Yantosca - Copied from "drydep_mod.f"
        1545             : !  See https://github.com/geoschem/hemco for complete history
        1546             : !EOP
        1547             : !------------------------------------------------------------------------------
        1548             : !BOC
        1549             : !
        1550             : ! !DEFINED PARAMETERS:
        1551             : !
        1552             :     REAL(hp), PARAMETER  :: XMAIR  = 28.8e-3_hp
        1553             :     REAL(hp), PARAMETER  :: RADAIR = 1.2e-10_hp
        1554             :     REAL(hp), PARAMETER  :: PI     = 3.14159265358979323e+0_hp
        1555             :     REAL(hp), PARAMETER  :: RADX   = 1.5e-10_hp
        1556             :     REAL(hp), PARAMETER  :: RGAS   = 8.32e+0_hp
        1557             :     REAL(hp), PARAMETER  :: AVOGAD = 6.022140857e+23_hp
        1558             : !
        1559             : ! !LOCAL VARIABLES:
        1560             : !
        1561             :     REAL(hp)             :: AIRDEN, Z, DIAM, FRPATH, SPEED
        1562             :
        1563             :     !=================================================================
        1564             :     ! DIFFG begins here!
        1565             :     !=================================================================
        1566             :
        1567             :     ! Air density
        1568           0 :     AIRDEN = ( PRESS * AVOGAD ) / ( RGAS * TK )
        1569             :
        1570             :     ! DIAM is the collision diameter for gas X with air.
        1571           0 :     DIAM   = RADX + RADAIR
        1572             :
        1573             :     ! Calculate the mean free path for gas X in air:
        1574             :     ! eq. 8.5 of Seinfeld [1986];
        1575           0 :     Z      = XM  / XMAIR
        1576           0 :     FRPATH = 1e+0_hp /( PI * SQRT( 1e+0_hp + Z ) * AIRDEN*( DIAM**2 ) )
        1577             :
        1578             :     ! Calculate average speed of gas X; eq. 15.47 of Levine [1988]
        1579           0 :     SPEED  = SQRT( 8e+0_hp * RGAS * TK / ( PI * XM ) )
        1580             :
        1581             :     ! Calculate diffusion coefficient of gas X in air;
        1582             :     ! eq. 8.9 of Seinfeld [1986]
        1583           0 :     DIFF_G = ( 3e+0_hp * PI / 32e+0_hp ) * ( 1e+0_hp + Z ) * FRPATH * SPEED
        1584             :
        1585           0 :   END FUNCTION DiffG
        1586             : !EOC
        1587             : !------------------------------------------------------------------------------
        1588             : !                   Harmonized Emissions Component (HEMCO)                    !
        1589             : !------------------------------------------------------------------------------
        1590             : !BOP
        1591             : !
        1592             : ! !IROUTINE: Get_Dep_N
        1593             : !
        1594             : ! !DESCRIPTION: Subroutine GET\_DEP\_N sums dry and wet deposition since prev.
        1595             : ! timestep and calculates contribution to fertilizer N source. Output is in
        1596             : ! kg NO/m2.
        1597             : !\\
        1598             : !\\
        1599             : ! !INTERFACE:
        1600             : !
        1601           0 :   SUBROUTINE Get_Dep_N( I, J, ExtState, HcoState, Inst, DEP_FERT )
        1602             : !
        1603             : ! !INPUT PARAMETERS:
        1604             : !
        1605             :     INTEGER,  INTENT(IN)         :: I
        1606             :     INTEGER,  INTENT(IN)         :: J
        1607             :     TYPE(Ext_State), POINTER     :: ExtState
        1608             :     TYPE(HCO_State), POINTER     :: HcoState
        1609             :     TYPE(MyInst),    POINTER     :: Inst
        1610             : !
        1611             : ! !INPUT/OUTPUT PARAMETERS:
        1612             : !
        1613             :     ! Dep emitted as Fert [kgNO/m2]
        1614             :     REAL(hp) ,  INTENT(INOUT) :: DEP_FERT
        1615             : !
        1616             : ! !REVISION HISTORY:
        1617             : !  See https://github.com/geoschem/hemco for complete history
        1618             : !EOP
        1619             : !------------------------------------------------------------------------------
        1620             : !BOC
        1621             : !
        1622             : ! !DEFINED PARAMETERS:
        1623             : !
        1624             :     REAL(hp),  PARAMETER :: TAU_MONTHS   = 6.0_hp ! Decay rate of dep. N [mon]
        1625             :     REAL(hp),  PARAMETER :: SECPERDAY    = 86400.0_hp
        1626             :     REAL(hp),  PARAMETER :: DAYSPERMONTH = 30.0_hp
        1627             : !
        1628             : ! !LOCAL VARIABLES:
        1629             : !
        1630             :     REAL(hp)             :: DRYN  ! Dry dep. N since prev timestep
        1631             :                                   ! Units ng N/m2/s
        1632             :     REAL(hp)             :: WETN  ! Wet dep. N since prev timestep
        1633             :     REAL(hp)             :: DEPN  ! dep. N since prev timestep
        1634             :
        1635             :     REAL(hp)             :: C1
        1636             :     REAL(hp)             :: C2
        1637             :     REAL(hp)             :: TAU_SEC
        1638             :     REAL(hp)             :: TS_SEC
        1639             :
        1640             :     !Total all N species & convert molec/cm2/s --> kg NO/m2/s
        1641           0 :     DRYN = SOURCE_DRYN( I, J, ExtState, HcoState, Inst )
        1642             :
        1643             :     !Total all N species & convert kg/s --> kg NO/m2/s
        1644           0 :     WETN = SOURCE_WETN( I, J, ExtState, HcoState )
        1645             :
        1646             :     ! Sum wet and dry deposition [kg NO/m2/s]
        1647           0 :     DEPN = DRYN + WETN
        1648             :
        1649             :     ! Emission Timestep in seconds
        1650           0 :     TS_SEC = HcoState%TS_EMIS
        1651             :
        1652             :     !Do mass balance (see Intro to Atm Chem Chap. 3)
        1653             :     !m(t) = m(0) * exp(-t/tau) + Source * tau * (1 - exp(-t/tau))
        1654             :
        1655             :     !convert months -->  seconds (assume 30 days months)
        1656           0 :     TAU_SEC = TAU_MONTHS * DAYSPERMONTH * SECPERDAY
        1657             :
        1658           0 :     C1 = EXP( -TS_SEC / TAU_SEC )
        1659           0 :     C2 = 1.0_hp - C1
        1660             :
        1661             :     ! kg NO/m2
        1662             :     ! NOTE: DEP_RESERVOIR is stored in kg NO/m3, but we just assume
        1663             :     ! that this is kg NO/m2.
        1664           0 :     Inst%DEP_RESERVOIR(I,J) = ( Inst%DEP_RESERVOIR (I,J) * C1 ) &
        1665           0 :                             + DEPN * TAU_SEC * C2
        1666             :
        1667             :     ! 40% runoff.
        1668           0 :     DEP_FERT = Inst%DEP_RESERVOIR(I,J) * 0.6_hp
        1669             :
        1670           0 :   END SUBROUTINE Get_Dep_N
        1671             : !EOC
        1672             : !------------------------------------------------------------------------------
        1673             : !                   Harmonized Emissions Component (HEMCO)                    !
        1674             : !------------------------------------------------------------------------------
        1675             : !BOP
        1676             : !
        1677             : ! !IROUTINE: Source_DryN
        1678             : !
        1679             : ! !DESCRIPTION: Subroutine SOURCE\_DRYN gets dry deposited Nitrogen since
        1680             : !               last emission time step, converts to kg NO/m2/s.
        1681             : !\\
        1682             : !\\
        1683             : ! !INTERFACE:
        1684             : !
        1685           0 :   FUNCTION Source_Dryn( I, J, ExtState, HcoState, Inst ) RESULT( DRYN )
        1686             : !
        1687             : ! !INPUT PARAMETERS:
        1688             : !
        1689             :     INTEGER,         INTENT(IN) :: I
        1690             :     INTEGER,         INTENT(IN) :: J
        1691             :     TYPE(Ext_State), POINTER    :: ExtState   ! Module options
        1692             :     TYPE(HCO_State), POINTER    :: HcoState   ! Output obj
        1693             :     TYPE(MyInst),    POINTER    :: Inst
        1694             : !
        1695             : ! !RETURN VALUE:
        1696             : !
        1697             :     REAL(hp)                      :: DRYN       ! Dry dep. N since prev timestep
        1698             : !
        1699             : ! !REVISION HISTORY:
        1700             : !  See https://github.com/geoschem/hemco for complete history
        1701             : !EOP
        1702             : !------------------------------------------------------------------------------
        1703             : !BOC
        1704             : !
        1705             : ! !LOCAL VARIABLES:
        1706             : !
        1707             :     REAL(hp),  PARAMETER   :: CM2_PER_M2  = 1.0e+4_hp
        1708             :     REAL(hp)               :: NTS
        1709             :
        1710             :     ! Divide through by number of chemistry timesteps
        1711             :     ! because DRY_TOTN is summed over chemistry timesteps
        1712             :     ! need to get average
        1713             :
        1714             :     !Molecules/cm2/s --> kg NO/m2/s
        1715           0 :     NTS  = HcoState%TS_EMIS / HcoState%TS_CHEM
        1716           0 :     DRYN = ExtState%DRY_TOTN%Arr%Val(I,J) * CM2_PER_M2 / NTS / &
        1717           0 :            HcoState%Phys%Avgdr * HcoState%Spc(Inst%IDTNO)%MW_g / 1000.0_hp
        1718             :
        1719           0 :   END FUNCTION Source_DryN
        1720             : !EOC
        1721             : !------------------------------------------------------------------------------
        1722             : !                   Harmonized Emissions Component (HEMCO)                    !
        1723             : !------------------------------------------------------------------------------
        1724             : !BOP
        1725             : !
        1726             : ! !IROUTINE: Source_WetN
        1727             : !
        1728             : ! !DESCRIPTION: Subroutine Source\_WetN gets wet deposited Nitrogen since
        1729             : !  last emission time step, converts to kg NO/m2/s.
        1730             : !\\
        1731             : !\\
        1732             : ! !INTERFACE:
        1733             : !
        1734           0 :   FUNCTION Source_WetN( I, J, ExtState, HcoState ) RESULT(WETN )
        1735             : !
        1736             : ! !INPUT PARAMETERS:
        1737             : !
        1738             :     INTEGER,         INTENT(IN) :: I
        1739             :     INTEGER,         INTENT(IN) :: J
        1740             :     TYPE(Ext_State), POINTER    :: ExtState   ! Module options
        1741             :     TYPE(HCO_State), POINTER    :: HcoState   ! Output obj
        1742             : !
        1743             : ! !RETURN VALUE:
        1744             : !
        1745             :     REAL(hp)                      :: WETN       ! Dry dep. N since prev timestep
        1746             : !
        1747             : ! !REVISION HISTORY:
        1748             : !  See https://github.com/geoschem/hemco for complete history
        1749             : !EOP
        1750             : !------------------------------------------------------------------------------
        1751             : !BOC
        1752             : !
        1753             : ! !LOCAL VARIABLES:
        1754             : !
        1755             :     REAL(hp)               :: NTS,    AREA_M2
        1756             :
        1757             :     ! Divide through by number of transport timesteps
        1758             :     ! because WET_TOTN is summed over transport timesteps
        1759             :     ! need to get average
        1760             :
        1761           0 :     NTS     = HcoState%TS_EMIS / HcoState%TS_DYN
        1762           0 :     AREA_M2 = HcoState%Grid%AREA_M2%Val(I,J)
        1763             :
        1764             :     ! Total N wet dep
        1765           0 :     WETN = ExtState%WET_TOTN%Arr%Val(I,J) / AREA_M2 / NTS
        1766             :
        1767           0 :   END FUNCTION Source_WetN
        1768             : !EOC
        1769             : !------------------------------------------------------------------------------
        1770             : !                   Harmonized Emissions Component (HEMCO)                    !
        1771             : !------------------------------------------------------------------------------
        1772             : !BOP
        1773             : !
        1774             : ! !IROUTINE: SoilTemp
        1775             : !
        1776             : ! !DESCRIPTION: Function SoilTemp computes the temperature-dependent term
        1777             : !  of the soil NOx emissions in ng N/m2/s and converts to molec/cm2/s
        1778             : !\\
        1779             : !\\
        1780             : ! !INTERFACE:
        1781             : !
        1782           0 :   FUNCTION SoilTemp( NN, TC, GWET ) RESULT( SOIL_TEMP )
        1783             : !
        1784             : ! !INPUT PARAMETERS:
        1785             : !
        1786             :     INTEGER,  INTENT(IN) :: NN           ! Soil biome type
        1787             :     REAL(hp), INTENT(IN) :: TC           ! Surface air temperature [C]
        1788             :     REAL(hp), INTENT(IN) :: GWET         ! Top soil moisture
        1789             : !
        1790             : ! !RETURN VALUE:
        1791             : !
        1792             :     REAL(hp)             :: SOIL_TEMP    ! Temperature-dependent term of
        1793             :                                          ! soil NOx emissions [unitless]
        1794             : !
        1795             : ! !REMARKS:
        1796             : !    Based on Ormeci et al., [1999] and Otter et al., [1999]
        1797             : !    there exists and entirely exponential relationship between
        1798             : !    temperature and soil NOx emissions at constant soil moisture
        1799             : !    Therefore we use the following relationship based
        1800             : !    on Yienger and Levy et al., [1995] for temperatures 0-30C:
        1801             : !
        1802             : !
        1803             : !         f(T) =  exp( 0.103+/-0.04 * T )
        1804             : !           in ng N/m2/s
        1805             : !
        1806             : !
        1807             : !     where T is the temperature in degrees Celsius....Below
        1808             : !     0 C, we assume emissions are zero because they are insignificant
        1809             : !     for the purposes of this global source. ...
        1810             : !
        1811             : !  References:
        1812             : !  ============================================================================
        1813             : !  (1 ) Ormeci, B., S. L. Sanin, and J. J. Pierce, Laboratory study of
        1814             : !        NO flux from agricultural soil: Effects of soil moisture, pH,
        1815             : !        and temperature, J. Geophys. Res., 104 ,16211629, 1999.
        1816             : !  (2 ) Otter, L. B., W. X. Yang, M. C. Scholes, and F. X. Meixner,
        1817             : !        Nitric oxide emissions from a southern African savanna, J.
        1818             : !        Geophys. Res., 105 , 20,69720,706, 1999.
        1819             : !  (3 ) Yienger, J.J, and H. Levy, Empirical model of global soil-biogenic
        1820             : !        NOx emissions, J. Geophys. Res., 100, D6, 11,447-11464, June 20, 1995.
        1821             : !
        1822             : ! !REVISION HISTORY:
        1823             : !  17 Aug 2009 - R. Yantosca - Initial Version
        1824             : !  See https://github.com/geoschem/hemco for complete history
        1825             : !EOP
        1826             : !------------------------------------------------------------------------------
        1827             : !BOC
        1828             : !
        1829             : ! !LOCAL VARIABLES
        1830             : !
        1831             :     REAL(hp)  :: TMMP
        1832             :
        1833             :     !==============================================================
        1834             :     ! 1) Convert from Surface Temp  --> Soil Temp
        1835             :     !==============================================================
        1836             :
        1837             :     ! Save surface air temp in shadow variable TMMP
        1838           0 :     TMMP   = TC
        1839             :
        1840             :     ! DRY
        1841           0 :     IF ( GWET < 0.3_hp ) THEN
        1842             :
        1843             :        ! Convert surface air temperature to model temperature
        1844             :        ! by adding 5 degrees C to model temperature
        1845           0 :        TMMP = TMMP + 5.0_hp
        1846             :
        1847             :     ! WET
        1848             :     ELSE
        1849             :
        1850           0 :        TMMP = SOILTA(NN) * TMMP + SOILTB(NN)
        1851             :
        1852             :     ENDIF
        1853             :
        1854             :     !==============================================================
        1855             :     ! 2) Compute Temperature Dependence
        1856             :     !==============================================================
        1857             :
        1858             :     ! Compute the soil temperature dependence term according
        1859             :     ! to equations 9b, 9a of Yienger & Levy [1995].
        1860             :     ! We now assume that soil response is exponential 0-30C
        1861             :     ! based on latest observations, caps at 30C
        1862             :
        1863           0 :     IF ( TMMP <= 0.0_hp ) THEN
        1864             :
        1865             :        ! No soil emissions if temp below freezing
        1866             :        SOIL_TEMP = 0.0_hp
        1867             :
        1868             :     ELSE
        1869             :
        1870             :        ! Caps temperature response at 30C
        1871           0 :        IF ( TMMP >= 30.0_hp ) TMMP = 30.0_hp
        1872             :
        1873           0 :        SOIL_TEMP =  EXP( 0.103_hp * TMMP )
        1874             :
        1875             :     ENDIF
        1876             :
        1877           0 :   END FUNCTION SoilTemp
        1878             : !EOC
        1879             : !----------------------------------------------------------------------------
        1880             : !                   Harmonized Emissions Component (HEMCO)                    !
        1881             : !------------------------------------------------------------------------------
        1882             : !BOP
        1883             : !
        1884             : ! !IROUTINE: SoilWet
        1885             : !
        1886             : ! !DESCRIPTION: Function SoilWet returns the soil moisture scaling
        1887             : !  of soil NOx emissions (values from 0-1).
        1888             : !\\
        1889             : !\\
        1890             : ! !INTERFACE:
        1891             : !
        1892           0 :   FUNCTION SoilWet( GWET, ARID, NONARID ) RESULT( WETSCALE )
        1893             : !
        1894             : ! !INPUT PARAMETERS:
        1895             : !
        1896             :     ! Top soil wetness [unitless]
        1897             :     REAL(hp), INTENT(IN) :: GWET
        1898             :
        1899             :     ! Fraction of arid & non-arid soil in the gridbox
        1900             :     REAL(hp), INTENT(IN) :: ARID
        1901             :     REAL(hp), INTENT(IN) :: NONARID
        1902             : !
        1903             : ! !RETURN_VALUE:
        1904             : !
        1905             :     ! A scaling term between 0-1 based on soil moisture
        1906             :     REAL(hp)             :: WETSCALE
        1907             : !
        1908             : ! !REMARKS:
        1909             : !  Soil moisture and temperature and now decoupled, the temperature
        1910             : !  term is scaled with a value from 0-1 based on water filled pore space
        1911             : !  WFPS in top-soil.
        1912             : !
        1913             : !  From N.E. Moore thesis:
        1914             : !  The response of SNOx is not monotonic to WFPS. SNOx are low for the
        1915             : !  extreme values of WFPS (0 and 1). For low values, emissions are
        1916             : !  substrate-limited. For high values, emissions are trapped and cannot
        1917             : !  diffuse to the surface [Yan et al., 2005]. SNOx dependence on soil
        1918             : !  moisture is best described as a Poisson function [Parsons et al., 1996;
        1919             : !  Otter et al., 1999; Pierce and Aneja, 2000; Kirkman et al., 2001;
        1920             : !  van Dijk and Meixner, 2001; van Dijk et al., 2002]:
        1921             : !
        1922             : !     scaling = a*x*exp(-b*x^2)
        1923             : !
        1924             : !  where the values of a and b are chosen such that the maximum value
        1925             : !  (unity) occurs for WFPS=0.3, which laboratory and field measurements have
        1926             : !  found to be the optimal value for emissions in most soils. The typical
        1927             : !  range of values are 0.2 (arid) up to 0.45 (floodplain)
        1928             : !  [Yang and Meixner, 1997; Ormeci et al., 1999].
        1929             : !
        1930             : !  Rice paddies no longer have to be scaled as in the Yienger & Levy model.
        1931             : !
        1932             : !  References:
        1933             : !  ============================================================================
        1934             : !  (1 ) Galbally, I. E., and R. Roy, Loss of fixed nitrogen from soils
        1935             : !        by nitric oxide exhalation, Nature, 275 , 734735, 1978.
        1936             : !  (2 ) Kirkman, G. A., W. X. Yang, and F. X. Meixner, Biogenic nitric
        1937             : !        oxide emissions upscaling: An approach for Zimbabwe, Global
        1938             : !        Biogeochemical Cycles, 15 ,1005 1020, 2001.
        1939             : !  (3 ) Ormeci, B., S. L. Sanin, and J. J. Pierce, Laboratory study of NO
        1940             : !        flux from agricultural soil: Effects of soil moisture, pH, and
        1941             : !        temperature, J. Geophys. Res., 104 , 16211629, 1999.
        1942             : !  (4 ) Otter, L. B., W. X. Yang, M. C. Scholes, and F. X. Meixner,
        1943             : !        Nitric oxide emissions from a southern African savanna, J.
        1944             : !        Geophys. Res., 105 , 20,69720,706, 1999.
        1945             : !  (5 ) Parsons, D. A., M. C. Scholes, R. J. Scholes, and J. S. Levine,
        1946             : !        Biogenic NO emissions from savanna soils as a function of fire
        1947             : !        regime, soil type, soil nitrogen, and water status, J. Geophys.
        1948             : !        Res., 101 , 23,68323,688, 1996.
        1949             : !  (6 ) Pierce, J. J., and V. P. Aneja, Nitric oxide emissions from
        1950             : !        engineered soil systems, Journal of Environmental Engineering,
        1951             : !        pp. 225232, 2000.
        1952             : !  (7 ) van Dijk, S. M., and J. H. Duyzer, Nitric oxide emissions from
        1953             : !        forest soils, J. Geophys. Res., 104 , 15,95515,961, 1999.
        1954             : !  (8 ) van Dijk, S. M., and F. X. Meixner, Production and consumption of
        1955             : !        NO in forest and pasture soils from the Amazon basin, Water, Air,
        1956             : !        and Soil Pollution: Focus 1 , pp. 119130, 2001.
        1957             : !  (9 ) Yang, W. X., and F. X. Meixner, Gaseous Nitrogen Emissions from
        1958             : !        Grasslands, CAB Int., Wallingford, UK, 1997, 67-71.
        1959             : !
        1960             : ! !REVISION HISTORY:
        1961             : !  31 Jan 2011 - R. Hudman   - Rewrote scaling scheme
        1962             : !  See https://github.com/geoschem/hemco for complete history
        1963             : !EOP
        1964             : !------------------------------------------------------------------------------
        1965             : !BOC
        1966             : !
        1967             :     !Scale by soil moisture
        1968           0 :     IF ( ARID >= NONARID .AND. ARID > 0.0_hp ) THEN
        1969             :        !Arid, max Poisson = 0.2
        1970           0 :        WETSCALE = 8.24_hp * GWET * EXP( -12.5_hp * GWET * GWET )
        1971             :     ELSE
        1972             :        !Non-arid, max Poisson = 0.3
        1973           0 :        WETSCALE = 5.5_hp * GWET * EXP( -5.55_hp * GWET * GWET )
        1974             :     ENDIF
        1975             :
        1976           0 :   END FUNCTION SoilWet
        1977             : !EOC
        1978             : !------------------------------------------------------------------------------
        1979             : !                   Harmonized Emissions Component (HEMCO)                    !
        1980             : !------------------------------------------------------------------------------
        1981             : !BOP
        1982             : !
        1983             : ! !IROUTINE: SoilCrf
        1984             : !
        1985             : ! !DESCRIPTION: Computes the canopy reduction factor for the soil NOx
        1986             : !  emissions according to Jacob \% Bakwin [1991] (and as used in Wang
        1987             : !  et al [1998]).
        1988             : !\\
        1989             : !\\
        1990             : ! !INTERFACE:
        1991             : !
        1992           0 :   FUNCTION SoilCrf( K, LAI, CPYNOX, WINDSQR, SUNCOS ) RESULT( SOIL_CRF )
        1993             : !
        1994             : ! !INPUT PARAMETERS:
        1995             : !
        1996             :     INTEGER,   INTENT(IN) :: K          ! Soil biome type
        1997             :     REAL(hp),  INTENT(IN) :: LAI        ! Leaf area index [cm2/cm2]
        1998             :     REAL(hp),  INTENT(IN) :: CPYNOX     ! Bulk sfc resistance to NOx [1/s]
        1999             :     REAL(hp),  INTENT(IN) :: WINDSQR    ! Square of sfc wind speed [m2/s2]
        2000             :     REAL(hp),  INTENT(IN) :: SUNCOS     ! Cosine of solar zenith angle
        2001             : !
        2002             : ! !RETURN_VALUE:
        2003             : !
        2004             :     REAL(hp)              :: SOIL_CRF   ! Canopy reduction factor (see below)
        2005             : !
        2006             : ! !REMARKS:
        2007             : !  Also note, CANOPYNOX (the bulk surface resistance to NOx) is computed
        2008             : !  in routine GET_CANOPY_NOx (in "canopy_nox_mod.f") and is passed here
        2009             : !  as an argument.
        2010             : !
        2011             : ! !REVISION HISTORY:
        2012             : !  17 Aug 2009 - R. Yantosca - Initial Version
        2013             : !  See https://github.com/geoschem/hemco for complete history
        2014             : !EOP
        2015             : !------------------------------------------------------------------------------
        2016             : !BOC
        2017             : !
        2018             : ! !DEFINED PARAMETERS:
        2019             : !
        2020             :     ! Ventilation velocity for NOx, day & night values [m/s]
        2021             :     REAL(hp),  PARAMETER :: VFDAY   = 1.0e-2_hp
        2022             :     REAL(hp),  PARAMETER :: VFNIGHT = 0.2e-2_hp
        2023             : !
        2024             : ! !LOCAL VARIABLES:
        2025             : !
        2026             :     REAL(hp) :: VFNEW
        2027             :
        2028             :     ! Pick proper ventilation velocity for day or night
        2029           0 :     IF ( SUNCOS > 0.0_hp ) THEN
        2030             :        VFNEW = VFDAY
        2031             :     ELSE
        2032           0 :        VFNEW = VFNIGHT
        2033             :     ENDIF
        2034             :
        2035             :     ! If the leaf area index and the bulk surface resistance
        2036             :     ! of the canopy to NOx deposition are both nonzero ...
        2037           0 :     IF ( LAI > 0.0_hp .and. CPYNOX > 0.0_hp ) THEN
        2038             :
        2039             :        ! Adjust the ventilation velocity.
        2040             :        ! NOTE: SOILEXC(21) is the canopy wind extinction
        2041             :        ! coefficient for the tropical rainforest biome.
        2042             :        VFNEW    = ( VFNEW        * SQRT( WINDSQR/9.0_hp * 7.0_hp/LAI     )  &
        2043           0 :                 * ( SOILEXC(21)  / SOILEXC(K)                            ) )
        2044             :
        2045             :        ! Soil canopy reduction factor
        2046           0 :        SOIL_CRF = CPYNOX / ( CPYNOX + VFNEW )
        2047             :
        2048             :     ELSE
        2049             :
        2050             :        ! Otherwise set the soil canopy reduction factor to zero
        2051             :        SOIL_CRF = 0.0_hp
        2052             :
        2053             :     ENDIF
        2054             :
        2055           0 :   END FUNCTION SoilCrf
        2056             : !EOC
        2057             : !-----------------------------------------------------------------------------
        2058             : !                   Harmonized Emissions Component (HEMCO)                    !
        2059             : !------------------------------------------------------------------------------
        2060             : !BOP
        2061             : !
        2062             : ! !IROUTINE: FertAdd
        2063             : !
        2064             : ! !DESCRIPTION: Function FertAdd computes fertilizer emissions
        2065             : !\\
        2066             : !\\
        2067             : ! !INTERFACE:
        2068             : !
        2069           0 :   FUNCTION FertAdd( SOILFRT, DEPN ) RESULT( FERT_ADD )
        2070             : !
        2071             : ! !INPUT PARAMETERS:
        2072             : !
        2073             :     REAL(hp), INTENT(IN) :: DEPN       ! N emissions from deposition
        2074             :     REAL(hp), INTENT(IN) :: SOILFRT    ! N emissions from fertilizers
        2075             :                                        !  read in from disk and passed
        2076             :                                        !  here as an argument [ng N/m2/s]
        2077             : !
        2078             : ! !RETURN_VALUE:
        2079             : !
        2080             :     REAL(hp)             :: FERT_ADD   ! Total Fert emissions
        2081             : !
        2082             : ! !REMARKS:
        2083             : !  We use a new spatially explicit data set of chemical and manure fert
        2084             : !  (native resolution 0.5\B0x0.5\B0) from Potter et al., [2010]
        2085             : !  distributed using MODIS EVI seasonality as described in
        2086             : !  N.E. Moore thesis, and Hudman et al., in prep.
        2087             : !
        2088             : !  In previous model, fertilizer emissions were emitted instantaneously as
        2089             : !  2.5% of applied fertilizer, independent of soil moisture/soil temperature,
        2090             : !  so that they were constant over the growing season.
        2091             : !
        2092             : !  Similar to the YL  parameterization, we now treat fertilizer emissions
        2093             : !  as part of the Aw. If we treat the wet biome coefficient as a measure of
        2094             : !  available N multiplied by a mean emission rate, we can treat fertilizer
        2095             : !  N in the same manner.
        2096             : !
        2097             : !  AW = SOILAW(BinewsoilAWS_08112011_emissonlyome) + N available in soil
        2098             : !       x mean emission rate
        2099             : !
        2100             : !  Instead of choosing an emission rate for each box equivalent to 2.5%
        2101             : !  of applied N yearly as done in the YL scheme, we chose the mean emission
        2102             : !  rate so that the total global above canopy SNOx due to fertilizer matches
        2103             : !  observed estimates of fertilizer emissions of 1.8 Tg N yr-1 from Stehfest
        2104             : !  and Bouman [2006].  This treatment allows for interannual and daily
        2105             : !  variability in the strength  of response to temperature and precipitation.
        2106             : !  Note: this scaling must be set for each resolution.
        2107             : !
        2108             : !  References:
        2109             : !  ============================================================================
        2110             : !  (1 ) Potter, P., Ramankutty, N., Bennett, E.,  and Donner, S.:
        2111             : !        Characterizing the Spatial Patterns of Global Fertilizer
        2112             : !        Application and Manure Production, Earth Interactions,
        2113             : !        in press, 2010.
        2114             : !  (2 ) Stehfest, E. and L. Bouwman, N2O and NO emission from
        2115             : !        agricultural fields and soils under natural vegetation:
        2116             : !        summarizing available measurement data and modeling
        2117             : !        of global annual emissions, Nutrient Cycling in Agroecosystems
        2118             : !        (2006), 74:207-228 DOI 10.1007/s10705-006-9000-7.
        2119             : !
        2120             : ! !REVISION HISTORY:
        2121             : !  31 Jan 2011 - R. Hudman   - Rewrote pulsing scheme
        2122             : !  See https://github.com/geoschem/hemco for complete history
        2123             : !EOP
        2124             : !------------------------------------------------------------------------------
        2125             : !BOC
        2126             : !
        2127             :     ! Seconds per year
        2128             :     REAL(hp), PARAMETER      :: SEC_PER_YEAR = 3.1536e+7_hp
        2129             :
        2130             :
        2131             :     ! Soil fert and dep [ kg/m2 ], a measure of N avail. in soil
        2132           0 :     FERT_ADD = ( SOILFRT + DEPN ) / SEC_PER_YEAR
        2133             :
        2134             :     ! Convert [ng N/m2] --> [kg/m2/s]
        2135             :     ! (scale needed to force fert emiss of 1.8 Tg N/yr w/o canopy uptake)
        2136             : !    FERT_ADD = FERT_ADD * FERT_SCALE
        2137             :
        2138           0 :   END FUNCTION FERTADD
        2139             : !EOC
        2140             : !------------------------------------------------------------------------------
        2141             : !                   Harmonized Emissions Component (HEMCO)                    !
        2142             : !------------------------------------------------------------------------------
        2143             : !BOP
        2144             : !
        2145             : ! !IROUTINE: Pulsing
        2146             : !
        2147             : ! !DESCRIPTION: Function Pulsing calculates the increase (or "pulse") of
        2148             : !  soil NOx emission that happens after preciptiation falls on dry soil.
        2149             : !                                                                             .
        2150             : !  According to  Yan et al., [2005] , this pulsing process is thought to
        2151             : !  be due to a release of inorganic nitrogen trapped on top of the dry soil
        2152             : !  and a subsequent reactivation of water-stressed bacteria, which then
        2153             : !  metabolize the excess nitrogen. This can happen in seasonally dry
        2154             : !  grasslands and savannahs or over freshly fertilized fields.
        2155             : !\\
        2156             : !\\
        2157             : ! !INTERFACE:
        2158             : !
        2159           0 :   FUNCTION Pulsing( GWET,      TS_EMIS,             &
        2160             :                     GWET_PREV, PFACTOR, DRYPERIOD ) &
        2161             :                     RESULT( THE_PULSING )
        2162             : !
        2163             : ! !INPUT PARAMETERS:
        2164             : !
        2165             :     REAL(hp), INTENT(IN)    :: GWET        ! Soil Moisture
        2166             :     REAL*4,   INTENT(IN)    :: TS_EMIS     ! Emissions timestep [s]
        2167             :
        2168             : ! !INPUT/OUTPUT PARAMETERS:
        2169             : !
        2170             :     REAL(sp), INTENT(INOUT) :: GWET_PREV   ! Soil Moisture Prev timestep
        2171             :     REAL(sp), INTENT(INOUT) :: PFACTOR     ! Pulsing Factor
        2172             :     REAL(sp), INTENT(INOUT) :: DRYPERIOD   ! Dry period length in hours
        2173             : !
        2174             : ! !RETURN VALUE:
        2175             : !
        2176             :     REAL(hp)                :: THE_PULSING ! Factor to multiply baseline
        2177             :                                            ! emissions by to account for
        2178             :                                            ! soil pulsing of all types
        2179             : !
        2180             : ! !REMARKS:
        2181             : !  Soil NOx emissions consist of baseline emissions plus discrete "pulsing"
        2182             : !  episodes.  We follow thw Yan et al., [2005] algorithm, where the pulse
        2183             : !  (relative to the flux prewetting) is determined by the antecedent dry
        2184             : !  period, with a simple logarithmic relationship,
        2185             : !
        2186             : !     PFACTOR = 13.01 ln ( DRYPERIOD ) -  53.6
        2187             : !
        2188             : !  where PFACTOR is the magnitude of peak flux relative to prewetting flux,
        2189             : !  and DRYPERIOD  is the length of the antecedent dry period in hours.
        2190             : !
        2191             : !  The pulse decays with
        2192             : !
        2193             : !     PFACTOR = PFACTOR * EXP( -0.068e+0_hp * DTSRCE )
        2194             : !
        2195             : !  References:
        2196             : !  ============================================================================
        2197             : !  (1 ) Yan, X., T. Ohara, and H. Akimoto (2005), Statistical modeling of
        2198             : !        global soil NOx emissions, Global Biogeochem. Cycles, 19, GB3019,
        2199             : !        doi:10.1029/2004GB002276.Section 2.3.3
        2200             : !
        2201             : ! !REVISION HISTORY:
        2202             : !  31 Jan 2011 - R. Hudman   - Rewrote pulsing scheme
        2203             : !  See https://github.com/geoschem/hemco for complete history
        2204             : !EOP
        2205             : !------------------------------------------------------------------------------
        2206             : !BOC
        2207             : !
        2208             : ! !LOCAL VARIABLES:
        2209             : !
        2210             :     REAL(hp)  :: DTSRCE, GDIFF
        2211             :
        2212             :     !=================================================================
        2213             :     ! PULSING begins here!
        2214             :     !=================================================================
        2215             :
        2216             :     ! Emission timestep [s --> hours]
        2217           0 :     DTSRCE = TS_EMIS / 3600.0_hp
        2218             :
        2219             :     ! If soil moisture less than 0.3 and no pulse is taking place
        2220           0 :     IF ( GWET < 0.3_hp .and. PFACTOR == 1.0_hp) THEN
        2221             :
        2222             :        ! Get change in soil moisture since previous timestep
        2223           0 :        GDIFF = ( GWET - GWET_PREV )
        2224             :
        2225             :        ! If change in soil moisture is > 0.01 (rains)
        2226           0 :        IF ( GDIFF > 0.01_hp ) THEN
        2227             :
        2228             :           ! Initialize new pulse factor (dry period hours)
        2229           0 :           IF ( DRYPERIOD > 0.0_hp ) THEN
        2230           0 :              PFACTOR = 13.01_hp * LOG( DRYPERIOD ) - 53.6_hp
        2231             :           ELSE
        2232           0 :              PFACTOR = -53.6_hp
        2233             :           ENDIF
        2234             :
        2235             :           ! If dry period < ~3 days then no pulse
        2236           0 :           IF ( PFACTOR < 1.0_hp ) PFACTOR = 1.0_hp
        2237             :
        2238             :           ! Reinitialize dry period
        2239           0 :           DRYPERIOD = 0.0_hp
        2240             :
        2241             :        ! If no rain (i.e.,  change in soil moisture is < 0.01)
        2242             :        ELSE
        2243             :
        2244             :           ! Add one timestep to dry period
        2245           0 :           DRYPERIOD = DRYPERIOD + DTSRCE
        2246             :
        2247             :        ENDIF
        2248             :
        2249             :     ! If box is already pulsing , then decay pulse one timestep
        2250           0 :     ELSEIF ( PFACTOR /= 1.0_hp ) THEN
        2251             :
        2252             :        ! Decay pulse
        2253           0 :        PFACTOR   = PFACTOR * EXP( -0.068_hp * DTSRCE )
        2254             :
        2255             :        ! Update dry period
        2256           0 :        IF ( GWET < 0.3_hp ) DRYPERIOD = DRYPERIOD + DTSRCE
        2257             :
        2258             :        ! If end of pulse
        2259           0 :        IF ( PFACTOR < 1.0_hp ) PFACTOR = 1.0_hp
        2260             :
        2261             :     ENDIF
        2262             :
        2263             :     ! Update soil moisture holder for previous timestep
        2264           0 :     GWET_PREV = GWET
        2265             :
        2266             :     ! Return the pulsing factor
        2267           0 :     THE_PULSING = PFACTOR
        2268             :
        2269           0 :   END FUNCTION Pulsing
        2270             : !EOC
        2271             : !------------------------------------------------------------------------------
        2272             : !                   Harmonized Emissions Component (HEMCO)                    !
        2273             : !------------------------------------------------------------------------------
        2274             : !BOP
        2275             : !
        2276             : ! !IROUTINE: InstGet
        2277             : !
        2278             : ! !DESCRIPTION: Subroutine InstGet returns a pointer to the desired instance.
        2279             : !\\
        2280             : !\\
        2281             : ! !INTERFACE:
        2282             : !
        2283           0 :   SUBROUTINE InstGet ( Instance, Inst, RC, PrevInst )
        2284             : !
        2285             : ! !INPUT PARAMETERS:
        2286             : !
        2287             :     INTEGER                             :: Instance
        2288             :     TYPE(MyInst),     POINTER           :: Inst
        2289             :     INTEGER                             :: RC
        2290             :     TYPE(MyInst),     POINTER, OPTIONAL :: PrevInst
        2291             : !
        2292             : ! !REVISION HISTORY:
        2293             : !  18 Feb 2016 - C. Keller   - Initial version
        2294             : !  See https://github.com/geoschem/hemco for complete history
        2295             : !EOP
        2296             : !------------------------------------------------------------------------------
        2297             : !BOC
        2298             :     TYPE(MyInst),     POINTER    :: PrvInst
        2299             :
        2300             :     !=================================================================
        2301             :     ! InstGet begins here!
        2302             :     !=================================================================
        2303             :
        2304             :     ! Get instance. Also archive previous instance.
        2305           0 :     PrvInst => NULL()
        2306           0 :     Inst    => AllInst
        2307           0 :     DO WHILE ( ASSOCIATED(Inst) )
        2308           0 :        IF ( Inst%Instance == Instance ) EXIT
        2309           0 :        PrvInst => Inst
        2310           0 :        Inst    => Inst%NextInst
        2311             :     END DO
        2312           0 :     IF ( .NOT. ASSOCIATED( Inst ) ) THEN
        2313           0 :        RC = HCO_FAIL
        2314           0 :        RETURN
        2315             :     ENDIF
        2316             :
        2317             :     ! Pass output arguments
        2318           0 :     IF ( PRESENT(PrevInst) ) PrevInst => PrvInst
        2319             :
        2320             :     ! Cleanup & Return
        2321           0 :     PrvInst => NULL()
        2322           0 :     RC = HCO_SUCCESS
        2323             :
        2324             :   END SUBROUTINE InstGet
        2325             : !EOC
        2326             : !------------------------------------------------------------------------------
        2327             : !                   Harmonized Emissions Component (HEMCO)                    !
        2328             : !------------------------------------------------------------------------------
        2329             : !BOP
        2330             : !
        2331             : ! !IROUTINE: InstCreate
        2332             : !
        2333             : ! !DESCRIPTION: Subroutine InstCreate adds a new instance to the list of
        2334             : !  instances, assigns a unique instance number to this new instance, and
        2335             : !  archives this instance number to output argument Instance.
        2336             : !\\
        2337             : !\\
        2338             : ! !INTERFACE:
        2339             : !
        2340           0 :   SUBROUTINE InstCreate ( ExtNr, Instance, Inst, RC )
        2341             : !
        2342             : ! !INPUT PARAMETERS:
        2343             : !
        2344             :     INTEGER,       INTENT(IN)       :: ExtNr
        2345             : !
        2346             : ! !OUTPUT PARAMETERS:
        2347             : !
        2348             :     INTEGER,       INTENT(  OUT)    :: Instance
        2349             :     TYPE(MyInst),  POINTER          :: Inst
        2350             : !
        2351             : ! !INPUT/OUTPUT PARAMETERS:
        2352             : !
        2353             :     INTEGER,       INTENT(INOUT)    :: RC
        2354             : !
        2355             : ! !REVISION HISTORY:
        2356             : !  18 Feb 2016 - C. Keller   - Initial version
        2357             : !  See https://github.com/geoschem/hemco for complete history
        2358             : !EOP
        2359             : !------------------------------------------------------------------------------
        2360             : !BOC
        2361             :     TYPE(MyInst), POINTER          :: TmpInst
        2362             :     INTEGER                        :: nnInst
        2363             :
        2364             :     !=================================================================
        2365             :     ! InstCreate begins here!
        2366             :     !=================================================================
        2367             :
        2368             :     ! ----------------------------------------------------------------
        2369             :     ! Generic instance initialization
        2370             :     ! ----------------------------------------------------------------
        2371             :
        2372             :     ! Initialize
        2373           0 :     Inst => NULL()
        2374             :
        2375             :     ! Get number of already existing instances
        2376           0 :     TmpInst => AllInst
        2377           0 :     nnInst = 0
        2378           0 :     DO WHILE ( ASSOCIATED(TmpInst) )
        2379           0 :        nnInst  =  nnInst + 1
        2380           0 :        TmpInst => TmpInst%NextInst
        2381             :     END DO
        2382             :
        2383             :     ! Create new instance
        2384           0 :     ALLOCATE(Inst)
        2385           0 :     Inst%Instance = nnInst + 1
        2386           0 :     Inst%ExtNr    = ExtNr
        2387             :
        2388             :     ! Attach to instance list
        2389           0 :     Inst%NextInst => AllInst
        2390           0 :     AllInst       => Inst
        2391             :
        2392             :     ! Update output instance
        2393           0 :     Instance = Inst%Instance
        2394             :
        2395             :     ! ----------------------------------------------------------------
        2396             :     ! Type specific initialization statements follow below
        2397             :     ! ----------------------------------------------------------------
        2398             :
        2399             :     ! Make sure pointers are not dangling
        2400           0 :     Inst%DRYCOEFF  => NULL()
        2401           0 :     Inst%CLIMARID  => NULL()
        2402           0 :     Inst%CLIMNARID => NULL()
        2403           0 :     Inst%SOILFERT  => NULL()
        2404           0 :     Inst%LANDTYPE  => NULL()
        2405             :
        2406             :     ! Return w/ success
        2407           0 :     RC = HCO_SUCCESS
        2408             :
        2409           0 :   END SUBROUTINE InstCreate
        2410             : !EOC
        2411             : !------------------------------------------------------------------------------
        2412             : !                   Harmonized Emissions Component (HEMCO)                    !
        2413             : !------------------------------------------------------------------------------
        2414             : !BOP
        2415             : !
        2416             : ! !IROUTINE: InstRemove
        2417             : !
        2418             : ! !DESCRIPTION: Subroutine InstRemove removes an instance from the list of
        2419             : ! instances.
        2420             : !\\
        2421             : !\\
        2422             : ! !INTERFACE:
        2423             : !
        2424           0 :   SUBROUTINE InstRemove ( ExtState )
        2425             : !
        2426             : ! !INPUT PARAMETERS:
        2427             : !
        2428             :     TYPE(Ext_State),  POINTER       :: ExtState      ! Extension state
        2429             : !
        2430             : ! !REVISION HISTORY:
        2431             : !  18 Feb 2016 - C. Keller   - Initial version
        2432             : !  See https://github.com/geoschem/hemco for complete history
        2433             : !EOP
        2434             : !------------------------------------------------------------------------------
        2435             : !BOC
        2436             :     INTEGER                     :: RC
        2437             :     INTEGER                     :: I
        2438             :     TYPE(MyInst), POINTER       :: PrevInst
        2439             :     TYPE(MyInst), POINTER       :: Inst
        2440             :
        2441             :     !=================================================================
        2442             :     ! InstRemove begins here!
        2443             :     !=================================================================
        2444             :
        2445             :     ! Get instance. Also archive previous instance.
        2446           0 :     PrevInst => NULL()
        2447           0 :     Inst     => NULL()
        2448           0 :     CALL InstGet ( ExtState%SoilNOx, Inst, RC, PrevInst=PrevInst )
        2449             :
        2450             :     ! Instance-specific deallocation
        2451           0 :     IF ( ASSOCIATED(Inst) ) THEN
        2452             :
        2453             :        !---------------------------------------------------------------------
        2454             :        ! Deallocate fields of Inst before popping off from the list
        2455             :        ! in order to avoid memory leaks (Bob Yantosca (17 Aug 2022)
        2456             :        !---------------------------------------------------------------------
        2457           0 :        IF ( ALLOCATED( Inst%PFACTOR ) ) THEN
        2458           0 :           DEALLOCATE( Inst%PFACTOR  )
        2459             :        ENDIF
        2460             :
        2461           0 :        IF ( ALLOCATED( Inst%GWET_PREV ) ) THEN
        2462           0 :           DEALLOCATE( Inst%GWET_PREV )
        2463             :        ENDIF
        2464             :
        2465           0 :        IF ( ALLOCATED( Inst%CANOPYNOX ) ) THEN
        2466           0 :           DEALLOCATE( Inst%CANOPYNOX )
        2467             :        ENDIF
        2468             :
        2469           0 :        IF ( ALLOCATED( Inst%DEP_RESERVOIR ) ) THEN
        2470           0 :           DEALLOCATE( Inst%DEP_RESERVOIR )
        2471             :        ENDIF
        2472             :
        2473           0 :        IF ( ALLOCATED( Inst%SpcScalVal ) ) THEN
        2474           0 :           DEALLOCATE( Inst%SpcScalVal )
        2475             :        ENDIF
        2476             :
        2477           0 :        IF ( ALLOCATED( Inst%SpcScalFldNme ) ) THEN
        2478           0 :           DEALLOCATE( Inst%SpcScalFldNme )
        2479             :        ENDIF
        2480             :
        2481           0 :        IF ( ASSOCIATED( Inst%FertNO_Diag ) ) THEN
        2482           0 :           DEALLOCATE( Inst%FertNO_Diag )
        2483             :        ENDIF
        2484           0 :        Inst%FertNO_Diag => NULL()
        2485             :
        2486           0 :        IF ( ASSOCIATED( Inst%CLIMARID ) ) THEN
        2487           0 :           DEALLOCATE( Inst%CLIMARID )
        2488             :        ENDIF
        2489           0 :        Inst%CLIMARID => NULL()
        2490             :
        2491           0 :        IF ( ASSOCIATED( Inst%CLIMNARID ) ) THEN
        2492           0 :           DEALLOCATE( Inst%CLIMNARID )
        2493             :        ENDIF
        2494           0 :        Inst%CLIMNARID => NULL()
        2495             :
        2496           0 :        IF ( ASSOCIATED( Inst%SOILFERT ) ) THEN
        2497           0 :           DEALLOCATE( Inst%SOILFERT )
        2498             :        ENDIF
        2499           0 :        Inst%SOILFERT => NULL()
        2500             :
        2501             :        ! Deallocate LANDTYPE vector
        2502           0 :        IF ( ASSOCIATED(Inst%LANDTYPE) ) THEN
        2503           0 :           DO I = 1,NBIOM
        2504           0 :              IF ( ASSOCIATED( Inst%LANDTYPE(I)%VAL ) ) THEN
        2505           0 :                 DEALLOCATE( Inst%LANDTYPE(I)%Val )
        2506             :              ENDIF
        2507           0 :              Inst%LANDTYPE(I)%Val => NULL()
        2508             :           ENDDO
        2509           0 :           DEALLOCATE ( Inst%LANDTYPE )
        2510             :
        2511             :        ENDIF
        2512           0 :        Inst%LANDTYPE => NULL()
        2513             :
        2514             :        ! Eventually deallocate DRYCOEFF.
        2515             :        ! Make sure ExtState DRYCOEFF pointer is not dangling!
        2516           0 :        IF ( ASSOCIATED ( Inst%DRYCOEFF ) ) THEN
        2517           0 :           DEALLOCATE ( Inst%DRYCOEFF )
        2518           0 :           ExtState%DRYCOEFF => NULL()
        2519             :        ENDIF
        2520           0 :        Inst%DRYCOEFF => NULL()
        2521             :
        2522             :        !---------------------------------------------------------------------
        2523             :        ! Pop off instance from list
        2524             :        !---------------------------------------------------------------------
        2525           0 :        IF ( ASSOCIATED(PrevInst) ) THEN
        2526           0 :           PrevInst%NextInst => Inst%NextInst
        2527             :        ELSE
        2528           0 :           AllInst => Inst%NextInst
        2529             :        ENDIF
        2530           0 :        DEALLOCATE(Inst)
        2531             :     ENDIF
        2532             :
        2533             :     ! Free pointers before exiting
        2534           0 :     PrevInst => NULL()
        2535           0 :     Inst     => NULL()
        2536             :
        2537           0 :    END SUBROUTINE InstRemove
        2538             : !EOC
        2539           0 : END MODULE HCOX_SoilNOx_Mod
        2540             : !EOM

Generated by: LCOV version 1.14