diff --git a/bld/namelist_files/namelist_definition_ctsm.xml b/bld/namelist_files/namelist_definition_ctsm.xml
index 69a243bd27..736597a6ac 100644
--- a/bld/namelist_files/namelist_definition_ctsm.xml
+++ b/bld/namelist_files/namelist_definition_ctsm.xml
@@ -996,6 +996,11 @@ Method for distributing soil thickness across hillslope columns
Toggle to turn on the plant hydraulic stress model
+
+Toggle to turn on soil NOx in the N cycle
+
+
How LUNA and Photosynthesis (if needed) will get Leaf nitrogen content
diff --git a/src/biogeochem/CNBalanceCheckMod.F90 b/src/biogeochem/CNBalanceCheckMod.F90
index 7b88422843..9c96d314aa 100644
--- a/src/biogeochem/CNBalanceCheckMod.F90
+++ b/src/biogeochem/CNBalanceCheckMod.F90
@@ -10,7 +10,7 @@ module CNBalanceCheckMod
use shr_log_mod , only : errMsg => shr_log_errMsg
use decompMod , only : bounds_type, subgrid_level_gridcell, subgrid_level_column
use abortutils , only : endrun
- use clm_varctl , only : iulog, use_nitrif_denitrif, use_fates_bgc
+ use clm_varctl , only : iulog, use_nitrif_denitrif, use_fates_bgc, use_soil_nox
use clm_time_manager , only : get_step_size_real
use CNVegNitrogenFluxType , only : cnveg_nitrogenflux_type
use CNVegNitrogenStateType , only : cnveg_nitrogenstate_type
@@ -547,6 +547,14 @@ subroutine NBalanceCheck(this, bounds, num_soilc, filter_soilc, &
smin_no3_leached => soilbiogeochem_nitrogenflux_inst%smin_no3_leached_col , & ! Input: [real(r8) (:) ] (gN/m2/s) soil mineral NO3 pool loss to leaching
smin_no3_runoff => soilbiogeochem_nitrogenflux_inst%smin_no3_runoff_col , & ! Input: [real(r8) (:) ] (gN/m2/s) soil mineral NO3 pool loss to runoff
f_n2o_nit => soilbiogeochem_nitrogenflux_inst%f_n2o_nit_col , & ! Input: [real(r8) (:) ] (gN/m2/s) flux of N2o from nitrification
+
+ f_n2o_denit => soilbiogeochem_nitrogenflux_inst%f_n2o_denit_col , & ! Input: [real(r8) (:) ] (gN/m2/s) flux of N2o from denitrification
+ f_n2_denit => soilbiogeochem_nitrogenflux_inst%f_n2_denit_col , & ! Input: [real(r8) (:) ] (gN/m2/s) flux of N2 from denitrification
+ f_nox_denit => soilbiogeochem_nitrogenflux_inst%f_nox_denit_col , & ! Input:[real(r8) (:) ] (gN/m2/s) flux of NOx from denitrification
+ f_nox_nit => soilbiogeochem_nitrogenflux_inst%f_nox_nit_col , & ! Input:[real(r8) (:) ] (gN/m2/s) flux of NOx from nitrification
+ f_nox_denit_atmos => soilbiogeochem_nitrogenflux_inst%f_nox_denit_atmos_col , & ! Input:[real(r8) (:) ] (gN/m2/s) flux of NOx from denitrification
+ f_nox_nit_atmos => soilbiogeochem_nitrogenflux_inst%f_nox_nit_atmos_col , & ! Input: [real(r8) (:) ] (gN/m2/s) flux of NOx from nitrification
+
som_n_leached => soilbiogeochem_nitrogenflux_inst%som_n_leached_col , & ! Input: [real(r8) (:) ] (gN/m2/s) total SOM N loss from vertical transport
col_fire_nloss => cnveg_nitrogenflux_inst%fire_nloss_col , & ! Input: [real(r8) (:) ] (gN/m2/s) total column-level fire N loss
@@ -599,9 +607,15 @@ subroutine NBalanceCheck(this, bounds, num_soilc, filter_soilc, &
col_ninputs_partial(c) = col_ninputs(c)
! calculate total column-level outputs
-
- col_noutputs(c) = denit(c)
-
+ if (use_soil_nox) then
+ !denit=f_n2o_denit + f_n2_denit + f_nox_denit
+ col_noutputs(c) = f_n2_denit(c)+f_n2o_denit(c)
+ ! only f_nox_denit_atmos leaves the system, NOx denit trapped in the canopy stays in the system.
+ col_noutputs(c) = col_noutputs(c)+f_nox_denit_atmos(c)
+ else
+ col_noutputs(c) = denit(c)
+ endif
+
if( .not.col%is_fates(c) ) then
col_noutputs(c) = col_noutputs(c) + col_fire_nloss(c) + gru_conv_nflux(c)
@@ -626,8 +640,10 @@ subroutine NBalanceCheck(this, bounds, num_soilc, filter_soilc, &
col_noutputs(c) = col_noutputs(c) + sminn_leached(c)
else
col_noutputs(c) = col_noutputs(c) + f_n2o_nit(c)
-
col_noutputs(c) = col_noutputs(c) + smin_no3_leached(c) + smin_no3_runoff(c)
+ if (use_soil_nox) then
+ col_noutputs(c) = col_noutputs(c) + f_nox_nit_atmos(c)
+ endif
end if
col_noutputs(c) = col_noutputs(c) - som_n_leached(c)
@@ -652,7 +668,11 @@ subroutine NBalanceCheck(this, bounds, num_soilc, filter_soilc, &
if (abs(col_errnb(c)) > this%nwarning) then
write(iulog,*) 'nbalance warning at c =', c, col_errnb(c), col_endnb(c)
write(iulog,*)'inputs,ffix,nfix,ndep = ',ffix_to_sminn(c)*dt,nfix_to_sminn(c)*dt,ndep_to_sminn(c)*dt
- write(iulog,*)'outputs,lch,roff,dnit = ',smin_no3_leached(c)*dt, smin_no3_runoff(c)*dt,f_n2o_nit(c)*dt
+ write(iulog,*)'outputs,lch,roff,dnit = ',smin_no3_leached(c)*dt, smin_no3_runoff(c)*dt,denit(c)*dt
+ if (use_soil_nox) then
+ write(iulog,*)'outputs,fnoxnit,fnoxdenit,fn2 = ',f_nox_nit(c)*dt, f_nox_denit(c)*dt,f_n2_denit(c)*dt
+ write(iulog,*)'outputs,fn2onit,fn2odenit= ',f_n2o_nit(c)*dt, f_n2o_denit(c)*dt
+ end if
end if
end do ! end of columns loop
@@ -673,9 +693,9 @@ subroutine NBalanceCheck(this, bounds, num_soilc, filter_soilc, &
write(iulog,*)'inputs,ffix,nfix,ndep = ',ffix_to_sminn(c)*dt,nfix_to_sminn(c)*dt,ndep_to_sminn(c)*dt
end if
if(col%is_fates(c))then
- write(iulog,*)'outputs,lch,roff,dnit,plnt = ',smin_no3_leached(c)*dt, smin_no3_runoff(c)*dt,f_n2o_nit(c)*dt,sminn_to_plant(c)*dt
+ write(iulog,*)'outputs,lch,roff,dnit,plnt = ',smin_no3_leached(c)*dt, smin_no3_runoff(c)*dt,denit(c)*dt,sminn_to_plant(c)*dt
else
- write(iulog,*)'outputs,lch,roff,dnit = ',smin_no3_leached(c)*dt, smin_no3_runoff(c)*dt,f_n2o_nit(c)*dt
+ write(iulog,*)'outputs,lch,roff,dnit = ',smin_no3_leached(c)*dt, smin_no3_runoff(c)*dt,denit(c)*dt
end if
call endrun(subgrid_index=c, subgrid_level=subgrid_level_column, msg=errMsg(sourcefile, __LINE__))
end if
diff --git a/src/biogeochem/CNDriverMod.F90 b/src/biogeochem/CNDriverMod.F90
index d14e619b03..6f5e5593ca 100644
--- a/src/biogeochem/CNDriverMod.F90
+++ b/src/biogeochem/CNDriverMod.F90
@@ -45,6 +45,7 @@ module CNDriverMod
use CLMFatesInterfaceMod , only : hlm_fates_interface_type
use CropReprPoolsMod , only : nrepr
use SoilHydrologyType , only: soilhydrology_type
+ use DryDepVelocity , only : drydepvel_type
!
! !PUBLIC TYPES:
implicit none
@@ -103,7 +104,8 @@ subroutine CNDriverNoLeaching(bounds,
wateratm2lndbulk_inst, canopystate_inst, soilstate_inst, temperature_inst, &
soil_water_retention_curve, crop_inst, ch4_inst, &
dgvs_inst, photosyns_inst, saturated_excess_runoff_inst, energyflux_inst, &
- nutrient_competition_method, cnfire_method, dribble_crophrv_xsmrpool_2atm)
+ nutrient_competition_method, cnfire_method, dribble_crophrv_xsmrpool_2atm, &
+ drydepvel_inst)
!
! !DESCRIPTION:
! The core CN code is executed here. Calculates fluxes for maintenance
@@ -207,6 +209,7 @@ subroutine CNDriverNoLeaching(bounds,
class(fire_method_type) , intent(inout) :: cnfire_method
logical , intent(in) :: dribble_crophrv_xsmrpool_2atm
type(hlm_fates_interface_type) , intent(inout) :: clm_fates
+ type(drydepvel_type) , intent(inout) :: drydepvel_inst
!
! !LOCAL VARIABLES:
real(r8):: cn_decomp_pools(bounds%begc:bounds%endc,1:nlevdecomp,1:ndecomp_pools)
@@ -479,7 +482,7 @@ subroutine CNDriverNoLeaching(bounds,
cnveg_carbonflux_inst,cnveg_nitrogenstate_inst,cnveg_nitrogenflux_inst, &
soilbiogeochem_carbonflux_inst,&
soilbiogeochem_state_inst,soilbiogeochem_nitrogenstate_inst, &
- soilbiogeochem_nitrogenflux_inst,canopystate_inst)
+ soilbiogeochem_nitrogenflux_inst,canopystate_inst, drydepvel_inst)
call t_stopf('soilbiogeochemcompetition')
! distribute the available N between the competing patches on the basis of
diff --git a/src/biogeochem/CNVegetationFacade.F90 b/src/biogeochem/CNVegetationFacade.F90
index 9ff26dd56d..aa6f576a93 100644
--- a/src/biogeochem/CNVegetationFacade.F90
+++ b/src/biogeochem/CNVegetationFacade.F90
@@ -97,6 +97,7 @@ module CNVegetationFacade
use SoilWaterRetentionCurveMod , only : soil_water_retention_curve_type
use CLMFatesInterfaceMod , only : hlm_fates_interface_type
use SoilHydrologyType , only : soilhydrology_type
+ use DryDepVelocity , only : drydepvel_type
!
implicit none
private
@@ -128,7 +129,7 @@ module CNVegetationFacade
type(cn_balance_type) :: cn_balance_inst
class(fire_method_type), allocatable :: cnfire_method
type(dgvs_type) :: dgvs_inst
-
+ type(drydepvel_type) :: drydepvel_inst
! Control variables
logical, private :: reseed_dead_plants ! Flag to indicate if should reseed dead plants when starting up the model
logical, private :: dribble_crophrv_xsmrpool_2atm = .False. ! Flag to indicate if should harvest xsmrpool to the atmosphere
@@ -963,7 +964,7 @@ subroutine EcosystemDynamicsPreDrainage(this, bounds, &
wateratm2lndbulk_inst, canopystate_inst, soilstate_inst, temperature_inst, &
soil_water_retention_curve, crop_inst, ch4_inst, &
photosyns_inst, saturated_excess_runoff_inst, energyflux_inst, &
- nutrient_competition_method, fireemis_inst)
+ nutrient_competition_method, fireemis_inst, drydepvel_inst)
!
! !DESCRIPTION:
! Do the main science for biogeochemistry that needs to be done before hydrology-drainage
@@ -1020,6 +1021,7 @@ subroutine EcosystemDynamicsPreDrainage(this, bounds, &
class(nutrient_competition_method_type) , intent(inout) :: nutrient_competition_method
type(fireemis_type) , intent(inout) :: fireemis_inst
type(hlm_fates_interface_type) , intent(inout) :: clm_fates
+ type(drydepvel_type) , intent(inout) :: drydepvel_inst
!
! !LOCAL VARIABLES:
@@ -1054,7 +1056,8 @@ subroutine EcosystemDynamicsPreDrainage(this, bounds, &
wateratm2lndbulk_inst, canopystate_inst, soilstate_inst, temperature_inst, &
soil_water_retention_curve, crop_inst, ch4_inst, &
this%dgvs_inst, photosyns_inst, saturated_excess_runoff_inst, energyflux_inst, &
- nutrient_competition_method, this%cnfire_method, this%dribble_crophrv_xsmrpool_2atm)
+ nutrient_competition_method, this%cnfire_method, this%dribble_crophrv_xsmrpool_2atm, &
+ drydepvel_inst)
! fire carbon emissions
call CNFireEmisUpdate(bounds, num_bgc_vegp, filter_bgc_vegp, &
diff --git a/src/biogeochem/DryDepVelocity.F90 b/src/biogeochem/DryDepVelocity.F90
index b13d28c765..aa4ccec191 100644
--- a/src/biogeochem/DryDepVelocity.F90
+++ b/src/biogeochem/DryDepVelocity.F90
@@ -84,6 +84,7 @@ Module DryDepVelocity
real(r8), pointer, public :: velocity_patch (:,:) ! Dry Deposition Velocity
real(r8), pointer, private :: rs_drydep_patch (:) ! Stomatal resistance associated with dry deposition velocity for Ozone
+ real(r8), pointer, public :: crf_drydep_patch (:) ! Canopy reduction factor for NO associated with dry deposition velocity
contains
@@ -136,6 +137,7 @@ subroutine InitAllocate(this, bounds)
if ( n_drydep > 0 )then
allocate(this%velocity_patch(begp:endp, n_drydep)); this%velocity_patch(:,:) = nan
allocate(this%rs_drydep_patch(begp:endp)) ; this%rs_drydep_patch(:) = nan
+ allocate(this%crf_drydep_patch(begp:endp)) ; this%crf_drydep_patch(:) = nan
end if
end subroutine InitAllocate
@@ -180,6 +182,11 @@ subroutine InitHistory(this, bounds)
avgflag='A', long_name='Stomatal Resistance Associated with Ozone Dry Deposition Velocity', &
ptr_patch=this%rs_drydep_patch, default='inactive' )
+ this%crf_drydep_patch(begp:endp)= spval
+ call hist_addfld1d ( fname='CRF_DRYDEP_NO', units='unitless', &
+ avgflag='A', long_name='Canopy Reduction Factor Associated with NO Dry Deposition Velocity', &
+ ptr_patch=this%crf_drydep_patch, default='inactive')
+
end subroutine InitHistory
!-----------------------------------------------------------------------
@@ -253,7 +260,10 @@ subroutine depvel_compute( bounds, &
!mvm 11/30/2013
real(r8) :: rlu_lai ! constant to calculate rlu over bulk canopy
-
+ real(r8) :: ratio_stai_lai !ratio of stomata area index to lai [unitless]
+ real(r8) :: ks !factor to adjust SAI for soil canopy reduction [unitless] Yan et al (2005) https://doi.org/10.1029/2004GB002276
+ real(r8) :: kc !factor to adjust LAI for soil canopy reduction [unitless] Yan et al (2005) https://doi.org/10.1029/2004GB002276
+
logical :: has_dew
logical :: has_rain
real(r8), parameter :: rain_threshold = 1.e-7_r8 ! of the order of 1cm/day expressed in m/s
@@ -268,7 +278,7 @@ subroutine depvel_compute( bounds, &
real(r8) :: rc !combined surface resistance
real(r8) :: cts !correction to flu rcl and rgs for frost
real(r8) :: rlux_o3 !to calculate O3 leaf resistance in dew/rain conditions
-
+
! constants
real(r8), parameter :: slope = 0._r8 ! Used to calculate rdc in (lower canopy resistance)
integer, parameter :: wveg_unset = -1 ! Unset Wesley vegetation type
@@ -307,7 +317,9 @@ subroutine depvel_compute( bounds, &
annlai => canopystate_inst%annlai_patch , & ! Input: [real(r8) (:,:) ] 12 months of monthly lai from input data set
velocity => drydepvel_inst%velocity_patch , & ! Output: [real(r8) (:,:) ] cm/sec
- rs_drydep => drydepvel_inst%rs_drydep_patch & ! Output: [real(r8) (:) ] stomatal resistance associated with Ozone dry deposition velocity (s/m)
+ rs_drydep => drydepvel_inst%rs_drydep_patch , & ! Output: [real(r8) (:) ] stomatal resistance associated with Ozone dry deposition velocity (s/m)
+ crf_drydep => drydepvel_inst%crf_drydep_patch & ! Output: [real(r8) (:) ] canopy reduction factor for NO associated with dry deposition velocity [unitless]
+
)
!_________________________________________________________________
@@ -654,6 +666,52 @@ subroutine depvel_compute( bounds, &
(1._r8/(rdc+rclx(ispec)))+(1._r8/(rac(index_season,wesveg)+rgsx(ispec))))
rc = max( 10._r8, rc)
!
+
+ ! CANOPY REDUCTION FACTOR FOR soil NO fluxes
+ ! (Probably there may be a cleaner way to code this)
+ if ( drydep_list(ispec) == 'NO') then
+ !factors for canopy reduction
+ ks=11.6_r8
+ kc=0.32_r8
+
+ !no canopy reduction factor (ie zero update) over bare land or elai < 0.
+ if( clmveg == noveg) then
+ crf_drydep(pi)=1._r8
+ else
+ ratio_stai_lai = 0._r8 ! Initialize
+
+ !needleaf forest (coniferous forest)
+ if (clmveg == ndllf_evr_tmp_tree .or. clmveg == ndllf_evr_brl_tree .or. clmveg == ndllf_dcd_brl_tree ) ratio_stai_lai = 0.03_r8
+
+ !broadleaf evergreen (rainforest)
+ if (clmveg == nbrdlf_evr_trp_tree .or. clmveg == nbrdlf_evr_tmp_tree ) ratio_stai_lai = 0.015_r8
+
+ !broadleaf decidouos
+ if (clmveg == nbrdlf_dcd_trp_tree .or. clmveg == nbrdlf_dcd_tmp_tree .or. clmveg == nbrdlf_dcd_brl_tree) then
+ if (index_season == 1 .or. index_season == 5) then
+ ! spring/summer
+ ratio_stai_lai = 0.015_r8
+ else
+ !winter/fall
+ ratio_stai_lai = 0.005_r8
+ endif
+ endif
+
+ !shrubs
+ if (clmveg == nbrdlf_evr_shrub .or. clmveg == nbrdlf_dcd_tmp_shrub .or. clmveg == nbrdlf_dcd_brl_shrub ) ratio_stai_lai = 0.005_r8
+
+ !grass
+ if (clmveg == nc3_arctic_grass .or. clmveg == nc3_nonarctic_grass .or. clmveg == nc4_grass ) ratio_stai_lai = 0.005_r8
+
+ !crops
+ if (clmveg == nc3crop .or. clmveg == nc3irrig ) ratio_stai_lai = 0.008_r8
+ if (is_prognostic_crop(clmveg)) ratio_stai_lai = 0.008_r8
+
+ crf_drydep(pi)=(exp(-ks*elai(pi)*ratio_stai_lai)+exp(-kc*elai(pi)))/2
+ endif !no veg
+ endif ! drydep_list(ispec)
+
+
! assume no surface resistance for SO2 over water
!
if ( drydep_list(ispec) == 'SO2' .and. wesveg == 7 ) then
diff --git a/src/biogeophys/WaterStateType.F90 b/src/biogeophys/WaterStateType.F90
index 35441d65d9..157743c8dc 100644
--- a/src/biogeophys/WaterStateType.F90
+++ b/src/biogeophys/WaterStateType.F90
@@ -40,6 +40,9 @@ module WaterStateType
real(r8), pointer :: h2osoi_vol_col (:,:) ! col volumetric soil water (0<=h2osoi_vol<=watsat) [m3/m3] (nlevgrnd)
real(r8), pointer :: h2osoi_vol_prs_grc (:,:) ! grc volumetric soil water prescribed (0<=h2osoi_vol<=watsat) [m3/m3] (nlevgrnd)
real(r8), pointer :: h2osfc_col (:) ! col surface water (mm H2O)
+ real(r8), pointer :: h2osoi_vol_old_col (:,:) ! col volumetric soil water for previous time step (0<=h2osoi_vol<=watsat) [m3/m3] (nlevgrnd)
+ real(r8), pointer :: ldry_col (:,:) ! col dry period rain pulses for NOx from nitrification as state variable [hours]
+
real(r8), pointer :: snocan_patch (:) ! patch canopy snow water (mm H2O)
real(r8), pointer :: liqcan_patch (:) ! patch canopy liquid water (mm H2O)
@@ -142,6 +145,15 @@ subroutine InitAllocate(this, bounds, tracer_vars)
container = tracer_vars, &
bounds = bounds, subgrid_level = subgrid_level_column, &
dim2beg = -nlevsno+1, dim2end = nlevmaxurbgrnd)
+ call AllocateVar2d(var = this%h2osoi_vol_old_col, name = 'h2osoi_vol_old_col', &
+ container = tracer_vars, &
+ bounds = bounds, subgrid_level = subgrid_level_column, &
+ dim2beg = 1, dim2end = nlevmaxurbgrnd)
+ call AllocateVar2d(var = this%ldry_col, name = 'ldry_col', &
+ container = tracer_vars, &
+ bounds = bounds, subgrid_level = subgrid_level_column, &
+ dim2beg = 1, dim2end = nlevmaxurbgrnd)
+
call AllocateVar1d(var = this%snocan_patch, name = 'snocan_patch', &
container = tracer_vars, &
bounds = bounds, subgrid_level = subgrid_level_patch)
diff --git a/src/main/clm_driver.F90 b/src/main/clm_driver.F90
index a436b75a97..6a652a0565 100644
--- a/src/main/clm_driver.F90
+++ b/src/main/clm_driver.F90
@@ -1075,7 +1075,7 @@ subroutine clm_drv(doalb, nextsw_cday, declinp1, declin, rstwr, nlend, rdate, ro
water_inst%wateratm2lndbulk_inst, canopystate_inst, soilstate_inst, temperature_inst, &
soil_water_retention_curve, crop_inst, ch4_inst, &
photosyns_inst, saturated_excess_runoff_inst, energyflux_inst, &
- nutrient_competition_method, fireemis_inst)
+ nutrient_competition_method, fireemis_inst,drydepvel_inst)
call t_stopf('ecosysdyn')
end if
@@ -1590,7 +1590,10 @@ subroutine clm_drv_init(bounds, &
eflx_bot => energyflux_inst%eflx_bot_col , & ! Output: [real(r8) (:) ] heat flux from beneath soil/ice column (W/m**2)
cisun_z => photosyns_inst%cisun_z_patch , & ! Output: [real(r8) (:) ] intracellular sunlit leaf CO2 (Pa)
- cisha_z => photosyns_inst%cisha_z_patch & ! Output: [real(r8) (:) ] intracellular shaded leaf CO2 (Pa)
+ cisha_z => photosyns_inst%cisha_z_patch , & ! Output: [real(r8) (:) ] intracellular shaded leaf CO2 (Pa)
+
+ h2osoi_vol => waterstatebulk_inst%h2osoi_vol_col , & ! Input: [real(r8) (:,:) ] volumetric soil water (m3/m3)
+ h2osoi_vol_old => waterstatebulk_inst%h2osoi_vol_old_col & ! Output: [real(r8) (:,:) ] volumetric soil water for previous time step m3/m3
)
! Initialize intracellular CO2 (Pa) parameters each timestep for use in VOCEmission
@@ -1602,6 +1605,8 @@ subroutine clm_drv_init(bounds, &
do c = bounds%begc,bounds%endc
! Reset flux from beneath soil/ice column
eflx_bot(c) = 0._r8
+ !save volumetric water from previous time step-soil nox
+ h2osoi_vol_old(c, :)=h2osoi_vol(c, :)
end do
! Initialize fraction of vegetation not covered by snow
diff --git a/src/main/clm_varcon.F90 b/src/main/clm_varcon.F90
index 985660142a..9d864e2a75 100644
--- a/src/main/clm_varcon.F90
+++ b/src/main/clm_varcon.F90
@@ -199,6 +199,7 @@ module clm_varcon
!real(r8), parameter :: nitrif_n2o_loss_frac = 0.02_r8 ! fraction of N lost as N2O in nitrification (Parton et al., 2001)
real(r8), public, parameter :: nitrif_n2o_loss_frac = 6.e-4_r8 ! fraction of N lost as N2O in nitrification (Li et al., 2000)
+
real(r8), public, parameter :: frac_minrlztn_to_no3 = 0.2_r8 ! fraction of N mineralized that is dieverted to the nitrification stream (Parton et al., 2001)
!------------------------------------------------------------------
diff --git a/src/main/clm_varctl.F90 b/src/main/clm_varctl.F90
index 9a0d2901b9..5a892e69a5 100644
--- a/src/main/clm_varctl.F90
+++ b/src/main/clm_varctl.F90
@@ -456,6 +456,14 @@ module clm_varctl
logical, public :: use_hydrstress = .false. ! true => use plant hydraulic stress calculation
+ !----------------------------------------------------------
+ ! SOIL EMISSIONS (Flow of Soil NOx) switch
+ !----------------------------------------------------------
+ !
+ logical, public :: use_soil_nox = .false. ! true => add soil NOx emissions in N cycle
+
+
+
!----------------------------------------------------------
! glacier_mec control variables: default values (may be overwritten by namelist)
!----------------------------------------------------------
diff --git a/src/main/controlMod.F90 b/src/main/controlMod.F90
index 089503dc8b..09234d6713 100644
--- a/src/main/controlMod.F90
+++ b/src/main/controlMod.F90
@@ -319,7 +319,7 @@ subroutine control_init(dtime)
use_lch4, use_nitrif_denitrif, use_extralakelayers, &
use_vichydro, use_cn, use_cndv, use_crop, use_fertilizer, &
use_grainproduct, use_snicar_frc, use_vancouver, use_mexicocity, use_noio, &
- use_nguardrail, crop_residue_removal_frac, flush_gdd20, use_nvmovement
+ use_nguardrail, crop_residue_removal_frac, flush_gdd20, use_nvmovement, use_soil_nox
! SNICAR
namelist /clm_inparm/ &
@@ -747,7 +747,7 @@ subroutine control_spmd()
call mpi_bcast (use_mexicocity, 1, MPI_LOGICAL, 0, mpicom, ier)
call mpi_bcast (use_noio, 1, MPI_LOGICAL, 0, mpicom, ier)
call mpi_bcast (use_SSRE, 1, MPI_LOGICAL, 0, mpicom, ier)
-
+ call mpi_bcast (use_soil_nox, 1, MPI_LOGICAL, 0, mpicom, ier)
! initial file variables
call mpi_bcast (nrevsn, len(nrevsn), MPI_CHARACTER, 0, mpicom, ier)
call mpi_bcast (finidat, len(finidat), MPI_CHARACTER, 0, mpicom, ier)
@@ -1047,6 +1047,8 @@ subroutine control_print ()
write(iulog,*) ' use_mexicocity = ', use_mexicocity
write(iulog,*) ' use_noio = ', use_noio
write(iulog,*) ' use_SSRE = ', use_SSRE
+ write(iulog,*) ' use_soil_nox = ', use_soil_nox
+
write(iulog,*) 'input data files:'
write(iulog,*) ' PFT physiology and parameters file = ',trim(paramfile)
if (fsurdat == ' ') then
diff --git a/src/soilbiogeochem/SoilBiogeochemCompetitionMod.F90 b/src/soilbiogeochem/SoilBiogeochemCompetitionMod.F90
index 57bc82984e..3e44f014da 100644
--- a/src/soilbiogeochem/SoilBiogeochemCompetitionMod.F90
+++ b/src/soilbiogeochem/SoilBiogeochemCompetitionMod.F90
@@ -8,7 +8,7 @@ module SoilBiogeochemCompetitionMod
use shr_kind_mod , only : r8 => shr_kind_r8
use shr_log_mod , only : errMsg => shr_log_errMsg
use clm_varcon , only : dzsoi_decomp
- use clm_varctl , only : use_nitrif_denitrif
+ use clm_varctl , only : use_nitrif_denitrif, use_soil_nox
use abortutils , only : endrun
use decompMod , only : bounds_type
use SoilBiogeochemStateType , only : soilbiogeochem_state_type
@@ -30,6 +30,7 @@ module SoilBiogeochemCompetitionMod
use TemperatureType , only : temperature_type
use SoilStateType , only : soilstate_type
use CanopyStateType , only : CanopyState_type
+ use DryDepVelocity , only : drydepvel_type
!
implicit none
private
@@ -59,7 +60,7 @@ module SoilBiogeochemCompetitionMod
! !PRIVATE DATA MEMBERS:
real(r8) :: dt ! decomp timestep (seconds)
real(r8) :: bdnr ! bulk denitrification rate (1/s)
-
+ real(r8) :: dthr ! decomp timestep (hours) for soil nox rain pulsing
character(len=*), parameter, private :: sourcefile = &
__FILE__
!-----------------------------------------------------------------------
@@ -143,7 +144,7 @@ subroutine SoilBiogeochemCompetitionInit ( bounds)
! set time steps
dt = get_step_size_real()
-
+ dthr=dt/3600._r8 !time step in hours for soil nox nitrif rain pulsing
! set space-and-time parameters from parameter file
bdnr = params_inst%bdnr * (dt/secspday)
@@ -172,7 +173,8 @@ subroutine SoilBiogeochemCompetition (bounds, num_bgc_soilc, filter_bgc_soilc,nu
cnveg_carbonflux_inst,cnveg_nitrogenstate_inst,cnveg_nitrogenflux_inst, &
soilbiogeochem_carbonflux_inst, &
soilbiogeochem_state_inst, soilbiogeochem_nitrogenstate_inst, &
- soilbiogeochem_nitrogenflux_inst,canopystate_inst)
+ soilbiogeochem_nitrogenflux_inst,canopystate_inst, &
+ drydepvel_inst)
!
! !USES:
use clm_varctl , only: cnallocate_carbon_only, iulog
@@ -207,7 +209,8 @@ subroutine SoilBiogeochemCompetition (bounds, num_bgc_soilc, filter_bgc_soilc,nu
type(soilbiogeochem_nitrogenstate_type) , intent(inout) :: soilbiogeochem_nitrogenstate_inst
type(soilbiogeochem_nitrogenflux_type) , intent(inout) :: soilbiogeochem_nitrogenflux_inst
type(canopystate_type) , intent(inout) :: canopystate_inst
-!
+ type(drydepvel_type) , intent(in) :: drydepvel_inst !for NO soil canopy reduction
+
!
! !LOCAL VARIABLES:
integer :: c,p,l,pi,j,k ! indices
@@ -241,6 +244,10 @@ subroutine SoilBiogeochemCompetition (bounds, num_bgc_soilc, filter_bgc_soilc,nu
real(r8) :: residual_smin_no3(bounds%begc:bounds%endc)
real(r8) :: residual_plant_ndemand(bounds%begc:bounds%endc)
real(r8) :: sminn_to_plant_new(bounds%begc:bounds%endc)
+ real(r8) :: h2osoi_diff(bounds%begc:bounds%endc,1:nlevdecomp) !soil NO rain pulse
+ real(r8), pointer :: crf_drydep_col(:) !soil NO canopy_reduction
+ real(r8), target :: crf_drydep_col_target(bounds%begc:bounds%endc)
+
!-----------------------------------------------------------------------
associate( &
@@ -272,6 +279,22 @@ subroutine SoilBiogeochemCompetition (bounds, num_bgc_soilc, filter_bgc_soilc,nu
n2_n2o_ratio_denit_vr => soilbiogeochem_nitrogenflux_inst%n2_n2o_ratio_denit_vr_col , & ! Output: [real(r8) (:,:) ] ratio of N2 to N2O production by denitrification [gN/gN]
f_n2o_denit_vr => soilbiogeochem_nitrogenflux_inst%f_n2o_denit_vr_col , & ! Output: [real(r8) (:,:) ] flux of N2O from denitrification [gN/m3/s]
f_n2o_nit_vr => soilbiogeochem_nitrogenflux_inst%f_n2o_nit_vr_col , & ! Output: [real(r8) (:,:) ] flux of N2O from nitrification [gN/m3/s]
+!!!!SOIL NOx
+ f_n2_denit_vr => soilbiogeochem_nitrogenflux_inst%f_n2_denit_vr_col , & ! Output: [real(r8) (:,:) ] flux of N2 from denitrification [gN/m3/s]
+ nox_n2o_ratio_vr => soilbiogeochem_nitrogenflux_inst%nox_n2o_ratio_vr_col , & ! Output: [real(r8) (:,:) ] ratio of NOx to N2O production by nitrification and denitrification [gN/gN]
+ pot_f_nox_denit_vr => soilbiogeochem_nitrogenflux_inst%pot_f_nox_denit_vr_col , & ! Output: [real(r8) (:,:) ] potential flux of NOx from denitrification [gN/m3/s]
+ pot_f_nox_nit_vr => soilbiogeochem_nitrogenflux_inst%pot_f_nox_nit_vr_col , & ! Output: [real(r8) (:,:) ] potential flux of NOx from nitrification [gN/m3/s]
+ f_nox_denit_vr => soilbiogeochem_nitrogenflux_inst%f_nox_denit_vr_col , & ! Output: [real(r8) (:,:) ] flux of NOx from denitrification [gN/m3/s]
+ f_nox_nit_vr => soilbiogeochem_nitrogenflux_inst%f_nox_nit_vr_col , & ! Output: [real(r8) (:,:) ] flux of NOx from nitrification [gN/m3/s]
+ crf_drydep =>drydepvel_inst%crf_drydep_patch ,& ! Input: [real(r8) (:) ] canopy reduction factor for NOx associated with dry deposition velocity [unitless])
+ f_nox_denit_atmos_vr => soilbiogeochem_nitrogenflux_inst%f_nox_denit_atmos_vr_col , & ! Output: [real(r8) (:,:) ] flux of NOx from denitrification [gN/m3/s]
+ f_nox_nit_atmos_vr => soilbiogeochem_nitrogenflux_inst%f_nox_nit_atmos_vr_col , & ! Output: [real(r8) (:,:) ] flux of NOx from nitrification [gN/m3/s]
+ pfactor_vr => soilbiogeochem_nitrogenstate_inst%pfactor_vr_col , & ! Output: [real(r8) (:,:) ] rain pulse factor for NOx nitrification flux [unitless]
+ ldry_vr => soilbiogeochem_nitrogenstate_inst%ldry_vr_col , & ! In/Output: [real(r8) (:,:) ] dry period counter for rain pulse in NOx nitrification flux [hours]
+ h2osoi_diff_vr => soilbiogeochem_nitrogenflux_inst%h2osoi_diff_vr_col , & ! Output: [real(r8) (:,:) ] delta volumetric water content from day and previous day [m3/m3]
+ h2osoi_vol_old => waterstatebulk_inst%h2osoi_vol_old_col , & ! Input: [real(r8) (:,:) ] volumetric soil water previous time step (0<=h2osoi_vol<=watsat) [m3/m3] (nlevgrnd)
+ h2osoi_vol => waterstatebulk_inst%h2osoi_vol_col , & ! Input: [real(r8) (:,:) ] volumetric soil water (0<=h2osoi_vol<=watsat) [m3/m3]
+!!!!!
supplement_to_sminn_vr => soilbiogeochem_nitrogenflux_inst%supplement_to_sminn_vr_col , & ! Output: [real(r8) (:,:) ]
sminn_to_plant_vr => soilbiogeochem_nitrogenflux_inst%sminn_to_plant_vr_col , & ! Output: [real(r8) (:,:) ]
potential_immob_vr => soilbiogeochem_nitrogenflux_inst%potential_immob_vr_col , & ! Input: [real(r8) (:,:) ]
@@ -281,6 +304,13 @@ subroutine SoilBiogeochemCompetition (bounds, num_bgc_soilc, filter_bgc_soilc,nu
sminn_to_plant_fun_nh4_vr => soilbiogeochem_nitrogenflux_inst%sminn_to_plant_fun_nh4_vr_col & ! Iutput: [real(r8) (:) ] Total layer nh4 uptake of FUN (gN/m2/s)
)
+ !pft to column average for soil NO canopy reduction
+ crf_drydep_col => crf_drydep_col_target
+ call p2c(bounds, num_bgc_soilc, filter_bgc_soilc, &
+ crf_drydep(bounds%begp:bounds%endp), &
+ crf_drydep_col(bounds%begc:bounds%endc))
+
+
! calcualte nitrogen uptake profile
! nuptake_prof(:,:) = nan
! call SoilBiogelchemNitrogenUptakeProfile(bounds, &
@@ -721,6 +751,77 @@ subroutine SoilBiogeochemCompetition (bounds, num_bgc_soilc, filter_bgc_soilc,nu
f_n2o_nit_vr(c,j) = f_nit_vr(c,j) * nitrif_n2o_loss_frac
f_n2o_denit_vr(c,j) = f_denit_vr(c,j) / (1._r8 + n2_n2o_ratio_denit_vr(c,j))
+ if (use_soil_nox) then
+
+ !need to have released N2 to balance the cycle later
+ f_n2_denit_vr(c,j) = f_denit_vr(c,j)/ (1._r8 + 1._r8/n2_n2o_ratio_denit_vr(c,j))
+
+
+ !soil NOx emissions nitrification
+ pot_f_nox_nit_vr(c,j) = f_n2o_nit_vr(c,j) * nox_n2o_ratio_vr(c,j)
+
+ ! RAIN PULSES FOR NOX NITRIFICATION Hudman et al (2012) https://doi.org/10.5194/acp-12-7779-2012 Equation 4
+
+ !An increase of 0.5% (v/v) in 7 cm of surface soil is equivalent to about 3.5 mm of rainfall, which
+ !is the rainfall amount often reported to cause a pulse
+ ! Calculate the change in volumetric soil moisture since previous timestep
+ h2osoi_diff(c,j)=h2osoi_vol(c,j)-h2osoi_vol_old(c,j)
+ if (h2osoi_vol(c,j) < 0.30_r8 .and. pfactor_vr(c,j) .eq. 1._r8) then
+
+ ! If change in soil moisture is > 0.005 (0.5% v/v) and there is precipitation, then start pulsing
+ if ( h2osoi_diff(c,j) > 0.005_r8 ) then
+ ! Initialize new pulse factor (dry period hours)
+ if ( ldry_vr(c,j) > 72._r8 ) then
+ pfactor_vr(c,j) = 13._r8 *LOG(ldry_vr(c,j)) - 53.6_r8
+ else
+ pfactor_vr(c,j)=1._r8
+ endif
+ ! If dry period < ~3 days then no pulse. "13._r8 *LOG(ldry(c,j))" becomes > 53.6 at > 72 hours
+ ! Reinitialize dry period
+ ldry_vr(c,j) = 0._r8
+ else ! If no rain (i.e., change in soil moisture is < 0.5%)
+ ! Add one timestep (in hours) to dry period
+ ldry_vr(c,j) = ldry_vr(c,j) + dthr
+ endif
+ else
+ pfactor_vr(c,j) = pfactor_vr(c,j) * exp( -0.068_r8 * dthr )
+
+ ! Update dry period
+ If ( h2osoi_vol(c,j) < 0.30_r8 ) then
+ ldry_vr(c,j)=ldry_vr(c,j) + dthr
+ endif
+ endif
+
+ ! If end of pulse
+ if ( pfactor_vr(c,j) < 1._r8 ) then
+ pfactor_vr(c,j) = 1._r8
+ endif
+ ! Apply the pulse factor to NOx nitrification emission
+ pot_f_nox_nit_vr(c,j) = pot_f_nox_nit_vr(c,j) * pfactor_vr(c,j)
+
+ !soil NOx denitrification (DayCent code)
+ pot_f_nox_denit_vr(c,j) = f_n2o_denit_vr(c,j) * nox_n2o_ratio_vr(c,j)
+
+
+ ! final NOx cannot be larger than available nitrogen from nitrification, considering the N2O nitrified
+ f_nox_nit_vr(c,j) = min(pot_f_nox_nit_vr(c,j), max(0._r8,f_nit_vr(c,j)-f_n2o_nit_vr(c,j)))
+ ! final NOx from denitrification cannot be larger than available nitrogen from denitrification
+ f_nox_denit_vr(c,j) = min(pot_f_nox_denit_vr(c,j), max(0._r8, f_n2o_denit_vr(c,j)))
+
+ !some of the NOx denit needs to come from f_n2 as f_NOx+f_N2+f_N2O cannot be larger than Fdenit)
+ f_n2_denit_vr(c,j)=f_n2_denit_vr(c,j)-f_nox_denit_vr(c,j)
+
+ !remove NOx trapped in canopy, ie, NO that makes it into the atmosphere. N in canopy is assumed to returned to the soil for now.
+ !NOx flux above canopy
+ !sanity check, make sure CRF ranges from 0 to 1
+ crf_drydep_col(c)= min(1._r8, max(0._r8, crf_drydep_col(c)))
+
+ f_nox_nit_atmos_vr(c,j) =f_nox_nit_vr(c,j)* crf_drydep_col(c)
+ f_nox_denit_atmos_vr(c,j)=f_nox_denit_vr(c,j)*crf_drydep_col(c)
+ endif
+
+
+
! this code block controls the addition of N to sminn pool
! to eliminate any N limitation, when Carbon_Only is set. This lets the
@@ -744,6 +845,7 @@ subroutine SoilBiogeochemCompetition (bounds, num_bgc_soilc, filter_bgc_soilc,nu
end if
sminn_to_plant_vr(c,j) = smin_no3_to_plant_vr(c,j) + smin_nh4_to_plant_vr(c,j)
end if
+
! sum up no3 and nh4 fluxes
fpi_vr(c,j) = fpi_no3_vr(c,j) + fpi_nh4_vr(c,j)
@@ -752,6 +854,11 @@ subroutine SoilBiogeochemCompetition (bounds, num_bgc_soilc, filter_bgc_soilc,nu
end do
end do
+
+
+
+
+
if ( local_use_fun ) then
call t_startf( 'CNFUN' )
call CNFUN(bounds,num_bgc_soilc,filter_bgc_soilc,num_bgc_vegp,filter_bgc_vegp,waterstatebulk_inst,&
diff --git a/src/soilbiogeochem/SoilBiogeochemNStateUpdate1Mod.F90 b/src/soilbiogeochem/SoilBiogeochemNStateUpdate1Mod.F90
index 8655b6c72d..77a2935359 100644
--- a/src/soilbiogeochem/SoilBiogeochemNStateUpdate1Mod.F90
+++ b/src/soilbiogeochem/SoilBiogeochemNStateUpdate1Mod.F90
@@ -8,7 +8,7 @@ module SoilBiogeochemNStateUpdate1Mod
use shr_kind_mod , only: r8 => shr_kind_r8
use clm_time_manager , only : get_step_size_real
use clm_varpar , only : nlevdecomp, ndecomp_cascade_transitions
- use clm_varctl , only : iulog, use_nitrif_denitrif, use_crop
+ use clm_varctl , only : iulog, use_nitrif_denitrif, use_crop, use_soil_nox
use SoilBiogeochemDecompCascadeConType , only : use_soil_matrixcn
use clm_varcon , only : nitrif_n2o_loss_frac
use SoilBiogeochemStateType , only : soilbiogeochem_state_type
@@ -243,17 +243,37 @@ subroutine SoilBiogeochemNStateUpdate1(num_bgc_soilc, filter_bgc_soilc, &
ns%smin_no3_vr_col(c,j) = ns%smin_no3_vr_col(c,j) - nf%sminn_to_plant_fun_no3_vr_col(c,j)*dt
end if
-
+ if (use_soil_nox) then
! Account for nitrification fluxes
- ns%smin_nh4_vr_col(c,j) = ns%smin_nh4_vr_col(c,j) - nf%f_nit_vr_col(c,j) * dt
+ ns%smin_nh4_vr_col(c,j) = ns%smin_nh4_vr_col(c,j) - nf%f_nit_vr_col(c,j) * dt
- ns%smin_no3_vr_col(c,j) = ns%smin_no3_vr_col(c,j) + nf%f_nit_vr_col(c,j) * dt &
- * (1._r8 - nitrif_n2o_loss_frac)
+ ! Adding nitrified NH4 to NO3
+ ns%smin_no3_vr_col(c,j) = ns%smin_no3_vr_col(c,j) + nf%f_nit_vr_col(c,j) * dt
+
+ !Deducting N2O leakage from nitrification
+ ns%smin_no3_vr_col(c,j) = ns%smin_no3_vr_col(c,j) - nf%f_n2o_nit_vr_col(c,j) * dt
+ ! Deducting NOx leakage from nitrification. Only NOx_nit not trapped in the canopy leaves the system
+ ns%smin_no3_vr_col(c,j) = ns%smin_no3_vr_col(c,j) - nf%f_nox_nit_atmos_vr_col(c,j) * dt
! Account for denitrification fluxes
- ns%smin_no3_vr_col(c,j) = ns%smin_no3_vr_col(c,j) - nf%f_denit_vr_col(c,j) * dt
+ ! Deducting denitrification fluxes (denit accounts for N2O denit, N2 and NOx denit)
+ ns%smin_no3_vr_col(c,j) = ns%smin_no3_vr_col(c,j) - nf%f_n2o_denit_vr_col(c,j) * dt
+ ns%smin_no3_vr_col(c,j) = ns%smin_no3_vr_col(c,j) - nf%f_n2_denit_vr_col(c,j) * dt
+ ! Not all NOx denit leaves the system, only f_nox_denit_atmos
+ ns%smin_no3_vr_col(c,j) = ns%smin_no3_vr_col(c,j) - nf%f_nox_denit_atmos_vr_col(c,j)* dt
+ else
+ !ORIGINAL CODE
+ ! Account for nitrification fluxes
+ ns%smin_nh4_vr_col(c,j) = ns%smin_nh4_vr_col(c,j) - nf%f_nit_vr_col(c,j) * dt
+
+ ns%smin_no3_vr_col(c,j) = ns%smin_no3_vr_col(c,j) + nf%f_nit_vr_col(c,j) * dt &
+ * (1._r8 - nitrif_n2o_loss_frac)
+ ! Account for denitrification fluxes
+ ns%smin_no3_vr_col(c,j) = ns%smin_no3_vr_col(c,j) - nf%f_denit_vr_col(c,j) * dt
+ endif
+
! flux that prevents N limitation (when Carbon_only is set; put all into NH4)
ns%smin_nh4_vr_col(c,j) = ns%smin_nh4_vr_col(c,j) + nf%supplement_to_sminn_vr_col(c,j)*dt
diff --git a/src/soilbiogeochem/SoilBiogeochemNitrifDenitrifMod.F90 b/src/soilbiogeochem/SoilBiogeochemNitrifDenitrifMod.F90
index 44d013cece..f995846c8d 100644
--- a/src/soilbiogeochem/SoilBiogeochemNitrifDenitrifMod.F90
+++ b/src/soilbiogeochem/SoilBiogeochemNitrifDenitrifMod.F90
@@ -12,11 +12,11 @@ module SoilBiogeochemNitrifDenitrifMod
use clm_varpar , only : nlevdecomp
use clm_varcon , only : rpi, grav
use clm_varcon , only : d_con_g, d_con_w, secspday
- use clm_varctl , only : use_lch4
+ use clm_varctl , only : use_lch4, use_soil_nox
use abortutils , only : endrun
use decompMod , only : bounds_type
use SoilStatetype , only : soilstate_type
- use WaterStateBulkType , only : waterstatebulk_type
+ use WaterStateBulkType , only : waterstatebulk_type
use TemperatureType , only : temperature_type
use SoilBiogeochemCarbonFluxType , only : soilbiogeochem_carbonflux_type
use SoilBiogeochemNitrogenStateType , only : soilbiogeochem_nitrogenstate_type
@@ -47,7 +47,7 @@ module SoilBiogeochemNitrifDenitrifMod
type(params_type), private :: params_inst
- logical, public :: no_frozen_nitrif_denitrif = .false. ! stop nitrification and denitrification in frozen soils
+ logical, public :: no_frozen_nitrif_denitrif = .true. ! stop nitrification and denitrification in frozen soils
character(len=*), parameter, private :: sourcefile = &
__FILE__
@@ -189,6 +189,10 @@ subroutine SoilBiogeochemNitrifDenitrif(bounds, num_bgc_soilc, filter_bgc_soilc,
real(r8) :: organic_max ! organic matter content (kg/m3) where
! soil is assumed to act like peat
character(len=32) :: subname='nitrif_denitrif' ! subroutine name
+
+ real(r8) :: Dr ! relative gas diffusivity in soil compare to air
+ real(r8) :: fno3_co2 ! for denitrification rate
+
!-----------------------------------------------------------------------
associate( &
@@ -243,7 +247,10 @@ subroutine SoilBiogeochemNitrifDenitrif(bounds, num_bgc_soilc, filter_bgc_soilc,
pot_f_nit_vr => soilbiogeochem_nitrogenflux_inst%pot_f_nit_vr_col , & ! Output: [real(r8) (:,:) ] (gN/m3/s) potential soil nitrification flux
pot_f_denit_vr => soilbiogeochem_nitrogenflux_inst%pot_f_denit_vr_col , & ! Output: [real(r8) (:,:) ] (gN/m3/s) potential soil denitrification flux
- n2_n2o_ratio_denit_vr => soilbiogeochem_nitrogenflux_inst%n2_n2o_ratio_denit_vr_col & ! Output: [real(r8) (:,:) ] ratio of N2 to N2O production by denitrification [gN/gN]
+ n2_n2o_ratio_denit_vr => soilbiogeochem_nitrogenflux_inst%n2_n2o_ratio_denit_vr_col , & ! Output: [real(r8) (:,:) ] ratio of N2 to N2O production by denitrification [gN/gN]
+ afps_vr => soilbiogeochem_nitrogenflux_inst%afps_vr_col , & !Output: [real(r8) (:,:) ] portion of soil pore space that is filled with air
+ gross_nmin_vr => soilbiogeochem_nitrogenflux_inst%gross_nmin_vr_col , & ! Input: [real(r8) (:,:) ] ! to add missing term in nitrification rate
+ nox_n2o_ratio_vr => soilbiogeochem_nitrogenflux_inst%nox_n2o_ratio_vr_col & !Output: [real(r8) (:,:) ] ratio of NOx to N2O production, for soil NOx emissions [gN/gN]
)
surface_tension_water = params_inst%surface_tension_water
@@ -342,9 +349,10 @@ subroutine SoilBiogeochemNitrifDenitrif(bounds, num_bgc_soilc, filter_bgc_soilc,
! note that k_nitr_max_perday is converted from 1/day to 1/s
k_nitr_vr(c,j) = k_nitr_max_perday/secspday * k_nitr_t_vr(c,j) * k_nitr_h2o_vr(c,j) * k_nitr_ph_vr(c,j)
- ! first-order decay of ammonium pool with scalar defined above
+ ! first-order decay of ammonium pool with scalar defined above. It is missing the N from SOM turn over term but that affects the N plant uptake and decreases carbon cycle substantially
+ ! pot_f_nit_vr(c,j) = max(0.1_r8 * gross_nmin_vr(c,j) + smin_nh4_vr(c,j) * k_nitr_vr(c,j), 0._r8)
pot_f_nit_vr(c,j) = max(smin_nh4_vr(c,j) * k_nitr_vr(c,j), 0._r8)
-
+
! limit to oxic fraction of soils
pot_f_nit_vr(c,j) = pot_f_nit_vr(c,j) * (1._r8 - anaerobic_frac(c,j))
@@ -413,9 +421,25 @@ subroutine SoilBiogeochemNitrifDenitrif(bounds, num_bgc_soilc, filter_bgc_soilc,
wfps_vr(c,j) = max(min(h2osoi_vol(c,j)/watsat(c, j), 1._r8), 0._r8) * 100._r8
fr_WFPS(c,j) = max(0.1_r8, 0.015_r8 * wfps_vr(c,j) - 0.32_r8)
+ if (use_soil_nox) then
+ ! final ratio expression
+ !DAYCENT4.5 EQUATION modified n2_n2o_ratio_denit to be consistent with daycent code
+ !Parton et al 2001 Equation 5 https://doi.org/10.1029/2001JD900101
+ fno3_co2 = max(0.16_r8*ratio_k1(c,j), ratio_k1(c,j)*exp(-0.8_r8 * ratio_no3_co2(c,j)))
+ n2_n2o_ratio_denit_vr(c,j) = max(0.1_r8,fno3_co2*fr_WFPS(c,j))
+ ! Dr from Zhao et al. (2017) Equation A2 https://doi.org/10.5194/acp-17-9781-2017
+ afps_vr(c,j) = 1._r8-max(min(h2osoi_vol(c,j)/watsat(c,j), 1._r8), 0._r8)
+ Dr = 0.209_r8*afps_vr(c,j)**(4._r8/3._r8)
+
+ nox_n2o_ratio_vr(c,j) = 15.2_r8 + (35.5_r8 * atan(0.68_r8 * rpi * (10.0_r8 * Dr - 1.86_r8))) / rpi
+
+ nox_n2o_ratio_vr(c,j) = max(0._r8, nox_n2o_ratio_vr(c,j))
+
+ else
! final ratio expression
- n2_n2o_ratio_denit_vr(c,j) = max(0.16_r8*ratio_k1(c,j), ratio_k1(c,j)*exp(-0.8_r8 * ratio_no3_co2(c,j))) * fr_WFPS(c,j)
-
+ n2_n2o_ratio_denit_vr(c,j) = max(0.16*ratio_k1(c,j), ratio_k1(c,j)*exp(-0.8 * ratio_no3_co2(c,j))) * fr_WFPS(c,j)
+ endif
+
end do
end do
diff --git a/src/soilbiogeochem/SoilBiogeochemNitrogenFluxType.F90 b/src/soilbiogeochem/SoilBiogeochemNitrogenFluxType.F90
index 94d5b4324f..a0f43d77c7 100644
--- a/src/soilbiogeochem/SoilBiogeochemNitrogenFluxType.F90
+++ b/src/soilbiogeochem/SoilBiogeochemNitrogenFluxType.F90
@@ -7,7 +7,7 @@ module SoilBiogeochemNitrogenFluxType
use clm_varpar , only : nlevdecomp_full, nlevdecomp, ndecomp_pools_vr
use clm_varcon , only : spval, ispval, dzsoi_decomp
use decompMod , only : bounds_type
- use clm_varctl , only : use_nitrif_denitrif, use_crop, use_fates
+ use clm_varctl , only : use_nitrif_denitrif, use_crop, use_fates, use_soil_nox
use CNSharedParamsMod , only : use_fun
use SoilBiogeochemDecompCascadeConType , only : decomp_cascade_con, use_soil_matrixcn
use abortutils , only : endrun
@@ -66,7 +66,29 @@ module SoilBiogeochemNitrogenFluxType
real(r8), pointer :: f_n2o_denit_col (:) ! col flux of N2o from denitrification [gN/m^2/s]
real(r8), pointer :: f_n2o_nit_vr_col (:,:) ! col flux of N2o from nitrification [gN/m^3/s]
real(r8), pointer :: f_n2o_nit_col (:) ! col flux of N2o from nitrification [gN/m^2/s]
-
+!for NOx fluxes
+ real(r8), pointer :: f_n2_denit_vr_col (:,:) ! col flux of N2 from denitrification [gN/m^3/s]
+ real(r8), pointer :: f_n2_denit_col (:) ! col flux of N2 from denitrification [gN/m^2/s]
+ real(r8), pointer :: nox_n2o_ratio_vr_col (:,:) ! col ratio of NOx to N2O production by nitrification and denitrification [gN/gN]
+ real(r8), pointer :: f_nox_denit_vr_col (:,:) ! col flux of NOx from denitrifiation [gN/m3/s]
+ real(r8), pointer :: f_nox_denit_col (:) ! col (integrated) flux of NOx from denitrifiation [gN/m2/s]
+ real(r8), pointer :: f_nox_nit_vr_col (:,:) ! col flux of NOx from nitrification [gN/m3/s]
+ real(r8), pointer :: f_nox_nit_col (:) ! col (integrated) flux of NOx from nitrification [gN/m2/s]
+ real(r8), pointer :: h2osoi_diff_vr_col (:,:) ! col volumtric water content different [m3/m3]. FOR DIAGNOSTICS
+ real(r8), pointer :: pot_f_nox_denit_vr_col (:,:) ! col potential flux of NOx from denitrifiation [gN/m3/s]
+ real(r8), pointer :: pot_f_nox_denit_col (:) ! col potential (integrated) flux of NOx from denitrifiation [gN/m2/s]
+ real(r8), pointer :: pot_f_nox_nit_vr_col (:,:) ! col potential flux of NOx from nitrification [gN/m3/s]
+ real(r8), pointer :: pot_f_nox_nit_col (:) ! col potential (integrated) flux of NOx from nitrification [gN/m2/s]
+ real(r8), pointer :: f_nox_denit_atmos_vr_col (:,:) ! col flux of NOx from denitrifiation with canopy reduction [gN/m3/s]
+ real(r8), pointer :: f_nox_denit_atmos_col (:) ! col (integrated) flux of NOx from denitrifiation with canopy reduction [gN/m2/s]
+ real(r8), pointer :: f_nox_nit_atmos_vr_col (:,:) ! col flux of NOx from nitrification with canopy reduction [gN/m3/s]
+ real(r8), pointer :: f_nox_nit_atmos_col (:) ! col (integrated) flux of NOx from nitrification with canopy reduction [gN/m2/s]
+ real(r8), pointer :: soil_nox_atmos_col (:) ! col (integrated) flux of NOx from nitrification and denitrification with canopy reduction [gN/m2/s]. This is the output to the atmosphere
+ real(r8), pointer :: soil_nox_atmos_crop_col (:) ! col (integrated) flux of NOx from nitrification and denitrification from croplands with canopy reduction [gN/m2/s]. This is the output to the atmosphere
+ real(r8), pointer :: soil_nox_atmos_nat_col (:) ! col (integrated) flux of NOx from nitrification and denitrification from natural ecosystems with canopy reduction [gN/m2/s]. This is the output to the atmosphere
+ real(r8), pointer :: afps_vr_col (:,:) ! air soil porosity
+
+
! immobilization / uptake fluxes
real(r8), pointer :: actual_immob_no3_vr_col (:,:) ! col vertically-resolved actual immobilization of NO3 (gN/m3/s)
real(r8), pointer :: actual_immob_nh4_vr_col (:,:) ! col vertically-resolved actual immobilization of NH4 (gN/m3/s)
@@ -235,7 +257,28 @@ subroutine InitAllocate(this, bounds)
allocate(this%f_n2o_denit_vr_col (begc:endc,1:nlevdecomp_full)) ; this%f_n2o_denit_vr_col (:,:) = nan
allocate(this%f_n2o_nit_col (begc:endc)) ; this%f_n2o_nit_col (:) = nan
allocate(this%f_n2o_nit_vr_col (begc:endc,1:nlevdecomp_full)) ; this%f_n2o_nit_vr_col (:,:) = nan
-
+!for soil NO fluxes
+ allocate(this%nox_n2o_ratio_vr_col (begc:endc,1:nlevdecomp_full)) ; this%nox_n2o_ratio_vr_col (:,:) = nan
+ allocate(this%f_nox_denit_col (begc:endc)) ; this%f_nox_denit_col (:) = nan
+ allocate(this%f_nox_denit_vr_col (begc:endc,1:nlevdecomp_full)) ; this%f_nox_denit_vr_col (:,:) = nan
+ allocate(this%f_nox_nit_col (begc:endc)) ; this%f_nox_nit_col (:) = nan
+ allocate(this%f_nox_nit_vr_col (begc:endc,1:nlevdecomp_full)) ; this%f_nox_nit_vr_col (:,:) = nan
+ allocate(this%h2osoi_diff_vr_col (begc:endc,1:nlevdecomp_full)) ; this%h2osoi_diff_vr_col (:,:) = nan !for nox nit rain pulse
+ allocate(this%pot_f_nox_denit_col (begc:endc)) ; this%pot_f_nox_denit_col (:) = nan
+ allocate(this%pot_f_nox_denit_vr_col (begc:endc,1:nlevdecomp_full)) ; this%pot_f_nox_denit_vr_col (:,:) = nan
+ allocate(this%pot_f_nox_nit_col (begc:endc)) ; this%pot_f_nox_nit_col (:) = nan
+ allocate(this%pot_f_nox_nit_vr_col (begc:endc,1:nlevdecomp_full)) ; this%pot_f_nox_nit_vr_col (:,:) = nan
+ allocate(this%f_n2_denit_col (begc:endc)) ; this%f_n2_denit_col (:) = nan
+ allocate(this%f_n2_denit_vr_col (begc:endc,1:nlevdecomp_full)) ; this%f_n2_denit_vr_col (:,:) = nan
+ allocate(this%f_nox_denit_atmos_col (begc:endc)) ; this%f_nox_denit_atmos_col (:) = nan
+ allocate(this%f_nox_denit_atmos_vr_col (begc:endc,1:nlevdecomp_full)) ; this%f_nox_denit_atmos_vr_col (:,:) = nan
+ allocate(this%f_nox_nit_atmos_col (begc:endc)) ; this%f_nox_nit_atmos_col (:) = nan
+ allocate(this%f_nox_nit_atmos_vr_col (begc:endc,1:nlevdecomp_full)) ; this%f_nox_nit_atmos_vr_col (:,:) = nan
+ allocate(this%soil_nox_atmos_col (begc:endc)) ; this%soil_nox_atmos_col (:) = nan
+ allocate(this%soil_nox_atmos_crop_col (begc:endc)) ; this%soil_nox_atmos_crop_col (:) = nan
+ allocate(this%soil_nox_atmos_nat_col (begc:endc)) ; this%soil_nox_atmos_nat_col (:) = nan
+ allocate(this%afps_vr_col (begc:endc,1:nlevdecomp_full)) ; this%afps_vr_col (:,:) = nan
+
allocate(this%smin_no3_massdens_vr_col (begc:endc,1:nlevdecomp_full)) ; this%smin_no3_massdens_vr_col (:,:) = nan
allocate(this%soil_bulkdensity_col (begc:endc,1:nlevdecomp_full)) ; this%soil_bulkdensity_col (:,:) = nan
allocate(this%k_nitr_t_vr_col (begc:endc,1:nlevdecomp_full)) ; this%k_nitr_t_vr_col (:,:) = nan
@@ -857,6 +900,56 @@ subroutine InitHistory(this, bounds)
call hist_addfld1d (fname='F_N2O_DENIT', units='gN/m^2/s', &
avgflag='A', long_name='denitrification N2O flux', &
ptr_col=this%f_n2o_denit_col)
+
+ if (use_soil_nox) then
+
+ this%nox_n2o_ratio_vr_col(begc:endc,:) = spval
+ call hist_addfld_decomp (fname='NOx_N2O_RATIO', units='gN/gN', type2d='levdcmp', &
+ avgflag='A', long_name='NOx to N2O ratio', &
+ ptr_col=this%nox_n2o_ratio_vr_col)
+
+ this%f_nox_nit_col(begc:endc) = spval
+ call hist_addfld1d (fname='F_NOx_NIT', units='gN/m^2/s', &
+ avgflag='A', long_name='nitrification NOx flux', &
+ ptr_col=this%f_nox_nit_col)
+
+ this%f_nox_denit_col(begc:endc) = spval
+ call hist_addfld1d (fname='F_NOx_DENIT', units='gN/m^2/s', &
+ avgflag='A', long_name='denitrification NOx flux', &
+ ptr_col=this%f_nox_denit_col)
+
+ this%f_n2_denit_col(begc:endc) = spval
+ call hist_addfld1d (fname='F_N2_DENIT', units='gN/m^2/s', &
+ avgflag='A', long_name='denitrification N2 flux', &
+ ptr_col=this%f_n2_denit_col)
+
+ this%f_nox_nit_atmos_col(begc:endc) = spval
+ call hist_addfld1d (fname='F_NOx_NIT_ATMOS', units='gN/m^2/s', &
+ avgflag='A', long_name='nitrification NOx flux with canopy reduction', &
+ ptr_col=this%f_nox_nit_atmos_col)
+
+ this%f_nox_denit_atmos_col(begc:endc) = spval
+ call hist_addfld1d (fname='F_NOx_DENIT_ATMOS', units='gN/m^2/s', &
+ avgflag='A', long_name='denitrification NOx flux with canopy reduction', &
+ ptr_col=this%f_nox_denit_atmos_col)
+
+ this%soil_nox_atmos_col(begc:endc) = spval
+ call hist_addfld1d (fname='SOIL_NOx', units='gN/m^2/s', &
+ avgflag='A', long_name='total soil NOx (nit+denit) that is passed to the atmosphere)', &
+ ptr_col=this%soil_nox_atmos_col)
+
+ this%soil_nox_atmos_crop_col(begc:endc) = spval
+ call hist_addfld1d (fname='SOIL_NOx_CROP', units='gN/m^2/s', &
+ avgflag='A', long_name='soil NOx (nit+denit) from croplands that is passed to the atmosphere)', &
+ ptr_col=this%soil_nox_atmos_crop_col)
+
+ this%soil_nox_atmos_nat_col(begc:endc) = spval
+ call hist_addfld1d (fname='SOIL_NOx_NAT', units='gN/m^2/s', &
+ avgflag='A', long_name='soil NOx (nit+denit) from natural ecosystems that is passed to the atmosphere)', &
+ ptr_col=this%soil_nox_atmos_nat_col)
+ endif
+
+
end if
if (use_crop) then
@@ -954,7 +1047,18 @@ subroutine SetValues ( this, &
this%smin_nh4_to_plant_vr_col(i,j) = value_column
this%f_n2o_denit_vr_col(i,j) = value_column
this%f_n2o_nit_vr_col(i,j) = value_column
-
+ !soil NOx fluxes
+ this%f_n2_denit_vr_col(i,j) = value_column
+ this%f_nox_denit_vr_col(i,j) = value_column
+ this%f_nox_nit_vr_col(i,j) = value_column
+ this%nox_n2o_ratio_vr_col(i,j) = value_column
+ this%pot_f_nox_denit_vr_col(i,j) = value_column
+ this%pot_f_nox_nit_vr_col(i,j) = value_column
+ this%f_nox_denit_atmos_vr_col(i,j) = value_column
+ this%f_nox_nit_atmos_vr_col(i,j) = value_column
+ this%h2osoi_diff_vr_col(i,j) = value_column
+ this%afps_vr_col(i,j) = value_column
+
this%smin_no3_massdens_vr_col(i,j) = value_column
this%k_nitr_t_vr_col(i,j) = value_column
this%k_nitr_ph_vr_col(i,j) = value_column
@@ -1011,6 +1115,17 @@ subroutine SetValues ( this, &
this%f_n2o_nit_col(i) = value_column
this%smin_no3_leached_col(i) = value_column
this%smin_no3_runoff_col(i) = value_column
+ !soil NOx flux
+ this%f_n2_denit_col(i) = value_column
+ this%f_nox_denit_col(i) = value_column
+ this%f_nox_nit_col(i) = value_column
+ this%f_nox_denit_atmos_col(i) = value_column
+ this%f_nox_nit_atmos_col(i) = value_column
+ this%pot_f_nox_denit_col(i) = value_column
+ this%pot_f_nox_nit_col(i) = value_column
+ this%soil_nox_atmos_col(i) = value_column
+ this%soil_nox_atmos_crop_col(i) = value_column
+ this%soil_nox_atmos_nat_col(i) = value_column
else
this%sminn_to_denit_excess_col(i) = value_column
this%sminn_leached_col(i) = value_column
@@ -1081,6 +1196,7 @@ subroutine Summary(this, bounds, num_bgc_soilc, filter_bgc_soilc)
! !USES:
use clm_varpar , only: nlevdecomp, ndecomp_cascade_transitions,ndecomp_pools
use clm_varctl , only: use_nitrif_denitrif
+ use landunit_varcon , only : istcrop
!
! !ARGUMENTS:
class (soilbiogeochem_nitrogenflux_type) :: this
@@ -1193,7 +1309,7 @@ subroutine Summary(this, bounds, num_bgc_soilc, filter_bgc_soilc)
this%f_n2o_denit_col(c) = &
this%f_n2o_denit_col(c) + &
this%f_n2o_denit_vr_col(c,j) * dzsoi_decomp(j)
-
+
! leaching/runoff flux
this%smin_no3_leached_col(c) = &
this%smin_no3_leached_col(c) + &
@@ -1203,12 +1319,49 @@ subroutine Summary(this, bounds, num_bgc_soilc, filter_bgc_soilc)
this%smin_no3_runoff_col(c) + &
this%smin_no3_runoff_vr_col(c,j) * dzsoi_decomp(j)
+ !soil NOx fluxes
+ this%f_n2_denit_col(c) = &
+ this%f_n2_denit_col(c) + &
+ this%f_n2_denit_vr_col(c,j) * dzsoi_decomp(j)
+
+ this%f_nox_nit_col(c) = &
+ this%f_nox_nit_col(c) + &
+ this%f_nox_nit_vr_col(c,j) * dzsoi_decomp(j)
+
+ this%f_nox_denit_col(c) = &
+ this%f_nox_denit_col(c) + &
+ this%f_nox_denit_vr_col(c,j) * dzsoi_decomp(j)
+
+ this%f_nox_nit_atmos_col(c) = &
+ this%f_nox_nit_atmos_col(c) + &
+ this%f_nox_nit_atmos_vr_col(c,j) * dzsoi_decomp(j)
+
+ this%f_nox_denit_atmos_col(c) = &
+ this%f_nox_denit_atmos_col(c) + &
+ this%f_nox_denit_atmos_vr_col(c,j) * dzsoi_decomp(j)
+
+ this%pot_f_nox_nit_col(c) = &
+ this%pot_f_nox_nit_col(c) + &
+ this%pot_f_nox_nit_vr_col(c,j) * dzsoi_decomp(j)
+
+ this%pot_f_nox_denit_col(c) = &
+ this%pot_f_nox_denit_col(c) + &
+ this%pot_f_nox_denit_vr_col(c,j) * dzsoi_decomp(j)
+
end do
end do
do fc = 1,num_bgc_soilc
c = filter_bgc_soilc(fc)
this%denit_col(c) = this%f_denit_col(c)
+
+ this%soil_nox_atmos_col(c)=this%f_nox_nit_atmos_col(c)+this%f_nox_denit_atmos_col(c)
+
+ if (lun%itype(col%landunit(c)) == istcrop) then
+ this%soil_nox_atmos_crop_col(c)= this%f_nox_nit_atmos_col(c)+this%f_nox_denit_atmos_col(c)
+ else
+ this%soil_nox_atmos_nat_col(c)= this%f_nox_nit_atmos_col(c)+this%f_nox_denit_atmos_col(c)
+ end if
end do
end if
diff --git a/src/soilbiogeochem/SoilBiogeochemNitrogenStateType.F90 b/src/soilbiogeochem/SoilBiogeochemNitrogenStateType.F90
index d3042b02e8..ec9c514db8 100644
--- a/src/soilbiogeochem/SoilBiogeochemNitrogenStateType.F90
+++ b/src/soilbiogeochem/SoilBiogeochemNitrogenStateType.F90
@@ -66,6 +66,8 @@ module SoilBiogeochemNitrogenStateType
real(r8), pointer :: totn_col (:) ! (gN/m2) total column nitrogen, incl veg
real(r8), pointer :: totecosysn_col (:) ! (gN/m2) total ecosystem nitrogen, incl veg
real(r8), pointer :: totn_grc (:) ! (gN/m2) total gridcell nitrogen
+
+
! Matrix-cn
real(r8), pointer :: matrix_cap_decomp_npools_col (:,:) ! col (gN/m2) N capacity in decomposing (litter, cwd, soil) N pools in dimension (col,npools)
@@ -80,6 +82,10 @@ module SoilBiogeochemNitrogenStateType
type(sparse_matrix_type) :: AKXnacc ! col (gN/m3/yr) accumulated N transfers from j to i (col,i,j) per year in dimension(col,nlev*npools,nlev*npools) in sparse matrix type
type(vector_type) :: matrix_Ninter ! col (gN/m3) vertically-resolved decomposing (litter, cwd, soil) N pools in dimension(col,nlev*npools) in vector type
+ ! rain pulse for soil NOx
+ real(r8), pointer :: ldry_vr_col (:,:) ! dry period for rain pulses for NOx from nitrification [hours]
+ real(r8), pointer :: pfactor_vr_col (:,:) ! pulsing factor [unitless]
+
contains
procedure , public :: Init
@@ -190,7 +196,10 @@ subroutine InitAllocate(this, bounds)
allocate(this%totn_col (begc:endc)) ; this%totn_col (:) = nan
allocate(this%totecosysn_col (begc:endc)) ; this%totecosysn_col (:) = nan
allocate(this%totn_grc (bounds%begg:bounds%endg)) ; this%totn_grc (:) = nan
-
+ !soil nox-rain pulse
+ allocate(this%ldry_vr_col (begc:endc,1:nlevdecomp_full)) ; this%ldry_vr_col (:,:) = nan
+ allocate(this%pfactor_vr_col (begc:endc,1:nlevdecomp_full)) ; this%pfactor_vr_col (:,:) = nan
+
end subroutine InitAllocate
!------------------------------------------------------------------------
@@ -505,6 +514,10 @@ subroutine InitCold(this, bounds, &
do j = 1, nlevdecomp_full
this%smin_nh4_vr_col(c,j) = 0._r8
this%smin_no3_vr_col(c,j) = 0._r8
+ !soil nox-rail pulse
+ this%ldry_vr_col(c,j) = 0._r8
+ this%pfactor_vr_col(c,j) = 0._r8
+
end do
this%smin_nh4_col(c) = 0._r8
this%smin_no3_col(c) = 0._r8
@@ -956,6 +969,9 @@ subroutine SetValues ( this, num_column, filter_column, value_column )
if (use_nitrif_denitrif) then
this%smin_no3_vr_col(i,j) = value_column
this%smin_nh4_vr_col(i,j) = value_column
+ !soilnox-rain pulse
+ this%ldry_vr_col(i,j) = value_column
+ this%pfactor_vr_col(i,j) = value_column
end if
end do
end do