module SoilBiogeochemNitrifDenitrifMod !----------------------------------------------------------------------- ! !DESCRIPTION: ! Calculate nitrification and denitrification rates ! ! ! !USES: use shr_kind_mod , only : r8 => shr_kind_r8 use shr_const_mod , only : SHR_CONST_TKFRZ use shr_log_mod , only : errMsg => shr_log_errMsg 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 abortutils , only : endrun use decompMod , only : bounds_type use SoilStatetype , only : soilstate_type use WaterStateType , only : waterstate_type use TemperatureType , only : temperature_type use SoilBiogeochemCarbonFluxType , only : soilbiogeochem_carbonflux_type use SoilBiogeochemNitrogenStateType , only : soilbiogeochem_nitrogenstate_type use SoilBiogeochemNitrogenFluxType , only : soilbiogeochem_nitrogenflux_type use ch4Mod , only : ch4_type use ColumnType , only : col ! implicit none private ! public :: readParams ! Read in parameters from params file public :: nitrifReadNML ! Read in namelist public :: SoilBiogeochemNitrifDenitrif ! Calculate nitrification and ! type, private :: params_type real(r8) :: k_nitr_max ! maximum nitrification rate constant (1/s) real(r8) :: surface_tension_water ! surface tension of water(J/m^2), Arah an and Vinten 1995 real(r8) :: rij_kro_a ! Arah and Vinten 1995) real(r8) :: rij_kro_alpha ! parameter to calculate anoxic fraction of soil (Arah and Vinten 1995) real(r8) :: rij_kro_beta ! (Arah and Vinten 1995) real(r8) :: rij_kro_gamma ! (Arah and Vinten 1995) real(r8) :: rij_kro_delta ! (Arah and Vinten 1995) real(r8) :: denitrif_respiration_coefficient ! Multiplier for heterotrophic respiration for max denitrif rates real(r8) :: denitrif_respiration_exponent ! Exponents for heterotrophic respiration for max denitrif rates real(r8) :: denitrif_nitrateconc_coefficient ! Multiplier for nitrate concentration for max denitrif rates real(r8) :: denitrif_nitrateconc_exponent ! Exponent for nitrate concentration for max denitrif rates end type params_type type(params_type), private :: params_inst logical, public :: no_frozen_nitrif_denitrif = .false. ! stop nitrification and denitrification in frozen soils character(len=*), parameter, private :: sourcefile = & __FILE__ !----------------------------------------------------------------------- contains !----------------------------------------------------------------------- subroutine readParams ( ncid ) ! use ncdio_pio, only: file_desc_t,ncd_io ! ! !ARGUMENTS: type(file_desc_t),intent(inout) :: ncid ! pio netCDF file id ! ! !LOCAL VARIABLES: character(len=32) :: subname = 'CNNitrifDenitrifParamsType' character(len=100) :: errCode = '-Error reading in parameters file:' logical :: readv ! has variable been read in or not real(r8) :: tempr ! temporary to read in constant character(len=100) :: tString ! temp. var for reading !----------------------------------------------------------------------- ! ! read in constants ! tString='surface_tension_water' call ncd_io(trim(tString),tempr, 'read', ncid, readvar=readv) if ( .not. readv ) call endrun(msg=trim(errCode)//trim(tString)//errMsg(sourcefile, __LINE__)) params_inst%surface_tension_water=tempr tString='rij_kro_a' call ncd_io(trim(tString),tempr, 'read', ncid, readvar=readv) if ( .not. readv ) call endrun(msg=trim(errCode)//trim(tString)//errMsg(sourcefile, __LINE__)) params_inst%rij_kro_a=tempr tString='rij_kro_alpha' call ncd_io(trim(tString),tempr, 'read', ncid, readvar=readv) if ( .not. readv ) call endrun(msg=trim(errCode)//trim(tString)//errMsg(sourcefile, __LINE__)) params_inst%rij_kro_alpha=tempr tString='rij_kro_beta' call ncd_io(trim(tString),tempr, 'read', ncid, readvar=readv) if ( .not. readv ) call endrun(msg=trim(errCode)//trim(tString)//errMsg(sourcefile, __LINE__)) params_inst%rij_kro_beta=tempr tString='rij_kro_gamma' call ncd_io(trim(tString),tempr, 'read', ncid, readvar=readv) if ( .not. readv ) call endrun(msg=trim(errCode)//trim(tString)//errMsg(sourcefile, __LINE__)) params_inst%rij_kro_gamma=tempr tString='rij_kro_delta' call ncd_io(trim(tString),tempr, 'read', ncid, readvar=readv) if ( .not. readv ) call endrun(msg=trim(errCode)//trim(tString)//errMsg(sourcefile, __LINE__)) params_inst%rij_kro_delta=tempr end subroutine readParams !----------------------------------------------------------------------- subroutine nitrifReadNML( NLFilename ) ! ! !DESCRIPTION: ! Read the namelist for nitrification/denitrification ! ! !USES: use fileutils , only : getavu, relavu, opnfil use shr_nl_mod , only : shr_nl_find_group_name use spmdMod , only : masterproc, mpicom use shr_mpi_mod , only : shr_mpi_bcast use clm_varctl , only : iulog ! ! !ARGUMENTS: character(len=*), intent(in) :: NLFilename ! Namelist filename ! ! !LOCAL VARIABLES: integer :: ierr ! error code integer :: unitn ! unit for namelist file character(len=*), parameter :: subname = 'ReadNML' character(len=*), parameter :: nmlname = 'nitrif_inparm' !----------------------------------------------------------------------- real(r8) :: k_nitr_max_perday, denitrif_respiration_coefficient, & denitrif_respiration_exponent, denitrif_nitrateconc_coefficient, & denitrif_nitrateconc_exponent namelist /nitrif_inparm/ k_nitr_max_perday, denitrif_respiration_coefficient, & denitrif_respiration_exponent, denitrif_nitrateconc_coefficient, & denitrif_nitrateconc_exponent ! Initialize options to default values, in case they are not specified in ! the namelist denitrif_respiration_coefficient = 0.1_r8 denitrif_respiration_exponent = 1.3_r8 denitrif_nitrateconc_coefficient = 1.15_r8 denitrif_nitrateconc_exponent = 0.57_r8 k_nitr_max_perday = 0.1_r8 if (masterproc) then unitn = getavu() write(iulog,*) 'Read in '//nmlname//' namelist' call opnfil (NLFilename, unitn, 'F') call shr_nl_find_group_name(unitn, nmlname, status=ierr) if (ierr == 0) then read(unitn, nml=nitrif_inparm, iostat=ierr) if (ierr /= 0) then call endrun(msg="ERROR reading "//nmlname//"namelist"//errmsg(sourcefile, __LINE__)) end if else call endrun(msg="ERROR could NOT find "//nmlname//"namelist"//errmsg(sourcefile, __LINE__)) end if call relavu( unitn ) end if call shr_mpi_bcast (k_nitr_max_perday , mpicom) call shr_mpi_bcast (denitrif_respiration_coefficient , mpicom) call shr_mpi_bcast (denitrif_respiration_exponent , mpicom) call shr_mpi_bcast (denitrif_nitrateconc_coefficient , mpicom) call shr_mpi_bcast (denitrif_nitrateconc_exponent , mpicom) params_inst%k_nitr_max = k_nitr_max_perday / secspday ! Change units to per second params_inst%denitrif_respiration_coefficient = denitrif_respiration_coefficient params_inst%denitrif_respiration_exponent = denitrif_respiration_exponent params_inst%denitrif_nitrateconc_coefficient = denitrif_nitrateconc_coefficient params_inst%denitrif_nitrateconc_exponent = denitrif_nitrateconc_exponent if (masterproc) then write(iulog,*) ' ' write(iulog,*) nmlname//' settings:' write(iulog,nml=nitrif_inparm) write(iulog,*) ' ' end if end subroutine nitrifReadNML !----------------------------------------------------------------------- subroutine SoilBiogeochemNitrifDenitrif(bounds, num_soilc, filter_soilc, & soilstate_inst, waterstate_inst, temperature_inst, ch4_inst, & soilbiogeochem_carbonflux_inst, soilbiogeochem_nitrogenstate_inst, soilbiogeochem_nitrogenflux_inst) ! ! !DESCRIPTION: ! calculate nitrification and denitrification rates ! ! !USES: use clm_time_manager , only : get_curr_date, get_step_size use CNSharedParamsMod , only : anoxia_wtsat, CNParamsShareInst ! ! !ARGUMENTS: type(bounds_type) , intent(in) :: bounds integer , intent(in) :: num_soilc ! number of soil columns in filter integer , intent(in) :: filter_soilc(:) ! filter for soil columns type(soilstate_type) , intent(in) :: soilstate_inst type(waterstate_type) , intent(in) :: waterstate_inst type(temperature_type) , intent(in) :: temperature_inst type(ch4_type) , intent(in) :: ch4_inst type(soilbiogeochem_carbonflux_type) , intent(in) :: soilbiogeochem_carbonflux_inst type(soilbiogeochem_nitrogenstate_type) , intent(in) :: soilbiogeochem_nitrogenstate_inst type(soilbiogeochem_nitrogenflux_type) , intent(inout) :: soilbiogeochem_nitrogenflux_inst ! ! !LOCAL VARIABLES: integer :: c, fc, reflev, j real(r8) :: soil_hr_vr(bounds%begc:bounds%endc,1:nlevdecomp) ! total soil respiration rate (g C / m3 / s) real(r8) :: g_per_m3__to__ug_per_gsoil real(r8) :: g_per_m3_sec__to__ug_per_gsoil_day real(r8) :: mu, sigma real(r8) :: t real(r8) :: pH(bounds%begc:bounds%endc) !debug-- put these type structure for outing to hist files real(r8) :: co2diff_con(2) ! diffusion constants for CO2 real(r8) :: eps real(r8) :: f_a real(r8) :: surface_tension_water ! (J/m^2), Arah and Vinten 1995 real(r8) :: rij_kro_a ! Arah and Vinten 1995 real(r8) :: rij_kro_alpha ! Arah and Vinten 1995 real(r8) :: rij_kro_beta ! Arah and Vinten 1995 real(r8) :: rij_kro_gamma ! Arah and Vinten 1995 real(r8) :: rij_kro_delta ! Arah and Vinten 1995 real(r8) :: rho_w = 1.e3_r8 ! (kg/m3) real(r8) :: r_max real(r8) :: r_min(bounds%begc:bounds%endc,1:nlevdecomp) real(r8) :: ratio_diffusivity_water_gas(bounds%begc:bounds%endc,1:nlevdecomp) real(r8) :: om_frac real(r8) :: anaerobic_frac_sat, r_psi_sat, r_min_sat ! scalar values in sat portion for averaging 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 !----------------------------------------------------------------------- associate( & watsat => soilstate_inst%watsat_col , & ! Input: [real(r8) (:,:) ] volumetric soil water at saturation (porosity) (nlevgrnd) watfc => soilstate_inst%watfc_col , & ! Input: [real(r8) (:,:) ] volumetric soil water at field capacity (nlevsoi) bd => soilstate_inst%bd_col , & ! Input: [real(r8) (:,:) ] bulk density of dry soil material [kg/m3] bsw => soilstate_inst%bsw_col , & ! Input: [real(r8) (:,:) ] Clapp and Hornberger "b" (nlevgrnd) cellorg => soilstate_inst%cellorg_col , & ! Input: [real(r8) (:,:) ] column 3D org (kg/m3 organic matter) (nlevgrnd) sucsat => soilstate_inst%sucsat_col , & ! Input: [real(r8) (:,:) ] minimum soil suction (mm) soilpsi => soilstate_inst%soilpsi_col , & ! Input: [real(r8) (:,:) ] soil water potential in each soil layer (MPa) h2osoi_vol => waterstate_inst%h2osoi_vol_col , & ! Input: [real(r8) (:,:) ] volumetric soil water (0<=h2osoi_vol<=watsat) [m3/m3] (nlevgrnd) h2osoi_liq => waterstate_inst%h2osoi_liq_col , & ! Input: [real(r8) (:,:) ] liquid water (kg/m2) (new) (-nlevsno+1:nlevgrnd) t_soisno => temperature_inst%t_soisno_col , & ! Input: [real(r8) (:,:) ] soil temperature (Kelvin) (-nlevsno+1:nlevgrnd) o2_decomp_depth_unsat => ch4_inst%o2_decomp_depth_unsat_col , & ! Input: [real(r8) (:,:) ] O2 consumption during decomposition in each soil layer (nlevsoi) (mol/m3/s) conc_o2_unsat => ch4_inst%conc_o2_unsat_col , & ! Input: [real(r8) (:,:) ] O2 conc in each soil layer (mol/m3) (nlevsoi) o2_decomp_depth_sat => ch4_inst%o2_decomp_depth_sat_col , & ! Input: [real(r8) (:,:) ] O2 consumption during decomposition in each soil layer (nlevsoi) (mol/m3/s) conc_o2_sat => ch4_inst%conc_o2_sat_col , & ! Input: [real(r8) (:,:) ] O2 conc in each soil layer (mol/m3) (nlevsoi) finundated => ch4_inst%finundated_col , & ! Input: [real(r8) (:) ] fractional inundated area in soil column (excluding dedicated wetland columns) smin_nh4_vr => soilbiogeochem_nitrogenstate_inst%smin_nh4_vr_col , & ! Input: [real(r8) (:,:) ] (gN/m3) soil mineral NH4 pool smin_no3_vr => soilbiogeochem_nitrogenstate_inst%smin_no3_vr_col , & ! Input: [real(r8) (:,:) ] (gN/m3) soil mineral NO3 pool phr_vr => soilbiogeochem_carbonflux_inst%phr_vr_col , & ! Input: [real(r8) (:,:) ] potential hr (not N-limited) w_scalar => soilbiogeochem_carbonflux_inst%w_scalar_col , & ! Input: [real(r8) (:,:) ] soil water scalar for decomp t_scalar => soilbiogeochem_carbonflux_inst%t_scalar_col , & ! Input: [real(r8) (:,:) ] temperature scalar for decomp denit_resp_coef => params_inst%denitrif_respiration_coefficient , & ! Input: [real(r8) ] coefficient for max denitrification rate based on respiration denit_resp_exp => params_inst%denitrif_respiration_exponent , & ! Input: [real(r8) ] exponent for max denitrification rate based on respiration denit_nitrate_coef => params_inst%denitrif_nitrateconc_coefficient , & ! Input: [real(r8) ] coefficient for max denitrification rate based on nitrate concentration denit_nitrate_exp => params_inst%denitrif_nitrateconc_exponent , & ! Input: [real(r8) ] exponent for max denitrification rate based on nitrate concentration k_nitr_max => params_inst%k_nitr_max , & ! Input: r_psi => soilbiogeochem_nitrogenflux_inst%r_psi_col , & ! Output: [real(r8) (:,:) ] anaerobic_frac => soilbiogeochem_nitrogenflux_inst%anaerobic_frac_col , & ! Output: [real(r8) (:,:) ] ! ! subsets of the n flux calcs (for diagnostic/debugging purposes) smin_no3_massdens_vr => soilbiogeochem_nitrogenflux_inst%smin_no3_massdens_vr_col , & ! Output: [real(r8) (:,:) ] (ugN / g soil) soil nitrate concentration k_nitr_t_vr => soilbiogeochem_nitrogenflux_inst%k_nitr_t_vr_col , & ! Output: [real(r8) (:,:) ] k_nitr_ph_vr => soilbiogeochem_nitrogenflux_inst%k_nitr_ph_vr_col , & ! Output: [real(r8) (:,:) ] k_nitr_h2o_vr => soilbiogeochem_nitrogenflux_inst%k_nitr_h2o_vr_col , & ! Output: [real(r8) (:,:) ] k_nitr_vr => soilbiogeochem_nitrogenflux_inst%k_nitr_vr_col , & ! Output: [real(r8) (:,:) ] wfps_vr => soilbiogeochem_nitrogenflux_inst%wfps_vr_col , & ! Output: [real(r8) (:,:) ] fmax_denit_carbonsubstrate_vr => soilbiogeochem_nitrogenflux_inst%fmax_denit_carbonsubstrate_vr_col , & ! Output: [real(r8) (:,:) ] fmax_denit_nitrate_vr => soilbiogeochem_nitrogenflux_inst%fmax_denit_nitrate_vr_col , & ! Output: [real(r8) (:,:) ] f_denit_base_vr => soilbiogeochem_nitrogenflux_inst%f_denit_base_vr_col , & ! Output: [real(r8) (:,:) ] diffus => soilbiogeochem_nitrogenflux_inst%diffus_col , & ! Output: [real(r8) (:,:) ] diffusivity (unitless fraction of total diffusivity) ratio_k1 => soilbiogeochem_nitrogenflux_inst%ratio_k1_col , & ! Output: [real(r8) (:,:) ] ratio_no3_co2 => soilbiogeochem_nitrogenflux_inst%ratio_no3_co2_col , & ! Output: [real(r8) (:,:) ] soil_co2_prod => soilbiogeochem_nitrogenflux_inst%soil_co2_prod_col , & ! Output: [real(r8) (:,:) ] (ug C / g soil / day) fr_WFPS => soilbiogeochem_nitrogenflux_inst%fr_WFPS_col , & ! Output: [real(r8) (:,:) ] soil_bulkdensity => soilbiogeochem_nitrogenflux_inst%soil_bulkdensity_col , & ! Output: [real(r8) (:,:) ] (kg soil / m3) bulk density of soil (including water) 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] ) surface_tension_water = params_inst%surface_tension_water ! Set parameters from simple-structure model to calculate anoxic fratction (Arah and Vinten 1995) rij_kro_a = params_inst%rij_kro_a rij_kro_alpha = params_inst%rij_kro_alpha rij_kro_beta = params_inst%rij_kro_beta rij_kro_gamma = params_inst%rij_kro_gamma rij_kro_delta = params_inst%rij_kro_delta organic_max = CNParamsShareInst%organic_max pH(bounds%begc:bounds%endc) = 6.5 !!! set all soils with the same pH as placeholder here co2diff_con(1) = 0.1325_r8 co2diff_con(2) = 0.0009_r8 do j = 1, nlevdecomp do fc = 1,num_soilc c = filter_soilc(fc) !---------------- calculate soil anoxia state ! calculate gas diffusivity of soil at field capacity here ! use expression from methane code, but neglect OM for now f_a = 1._r8 - watfc(c,j) / watsat(c,j) eps = watsat(c,j)-watfc(c,j) ! Air-filled fraction of total soil volume ! use diffusivity calculation including peat if (use_lch4) then if (organic_max > 0._r8) then om_frac = min(cellorg(c,j)/organic_max, 1._r8) ! Use first power, not square as in iniTimeConst else om_frac = 1._r8 end if diffus (c,j) = (d_con_g(2,1) + d_con_g(2,2)*t_soisno(c,j)) * 1.e-4_r8 * & (om_frac * f_a**(10._r8/3._r8) / watsat(c,j)**2 + & (1._r8-om_frac) * eps**2 * f_a**(3._r8 / bsw(c,j)) ) ! calculate anoxic fraction of soils ! use rijtema and kroess model after Riley et al., 2000 ! caclulated r_psi as a function of psi r_min(c,j) = 2 * surface_tension_water / (rho_w * grav * abs(soilpsi(c,j))) r_max = 2 * surface_tension_water / (rho_w * grav * 0.1_r8) r_psi(c,j) = sqrt(r_min(c,j) * r_max) ratio_diffusivity_water_gas(c,j) = (d_con_g(2,1) + d_con_g(2,2)*t_soisno(c,j) ) * 1.e-4_r8 / & ((d_con_w(2,1) + d_con_w(2,2)*t_soisno(c,j) + d_con_w(2,3)*t_soisno(c,j)**2) * 1.e-9_r8) if (o2_decomp_depth_unsat(c,j) > 0._r8) then anaerobic_frac(c,j) = exp(-rij_kro_a * r_psi(c,j)**(-rij_kro_alpha) * & o2_decomp_depth_unsat(c,j)**(-rij_kro_beta) * & conc_o2_unsat(c,j)**rij_kro_gamma * (h2osoi_vol(c,j) + ratio_diffusivity_water_gas(c,j) * & watsat(c,j))**rij_kro_delta) else anaerobic_frac(c,j) = 0._r8 endif if (anoxia_wtsat) then ! Average saturated fraction values into anaerobic_frac(c,j). r_min_sat = 2._r8 * surface_tension_water / (rho_w * grav * abs(grav * 1.e-6_r8 * sucsat(c,j))) r_psi_sat = sqrt(r_min_sat * r_max) if (o2_decomp_depth_sat(c,j) > 0._r8) then anaerobic_frac_sat = exp(-rij_kro_a * r_psi_sat**(-rij_kro_alpha) * & o2_decomp_depth_sat(c,j)**(-rij_kro_beta) * & conc_o2_sat(c,j)**rij_kro_gamma * (watsat(c,j) + ratio_diffusivity_water_gas(c,j) * & watsat(c,j))**rij_kro_delta) else anaerobic_frac_sat = 0._r8 endif anaerobic_frac(c,j) = (1._r8 - finundated(c))*anaerobic_frac(c,j) + finundated(c)*anaerobic_frac_sat end if else ! NITRIF_DENITRIF requires Methane model to be active, ! otherwise diffusivity will be zeroed out here. EBK CDK 10/18/2011 anaerobic_frac(c,j) = 0._r8 diffus (c,j) = 0._r8 !call endrun(msg=' ERROR: NITRIF_DENITRIF requires Methane model to be active'//errMsg(sourcefile, __LINE__) ) end if !---------------- nitrification ! follows CENTURY nitrification scheme (Parton et al., (2001, 1996)) ! assume nitrification temp function equal to the HR scalar k_nitr_t_vr(c,j) = min(t_scalar(c,j), 1._r8) ! ph function from Parton et al., (2001, 1996) k_nitr_ph_vr(c,j) = 0.56 + atan(rpi * 0.45 * (-5.+ pH(c)))/rpi ! moisture function-- assume the same moisture function as limits heterotrophic respiration ! Parton et al. base their nitrification- soil moisture rate constants based on heterotrophic rates-- can we do the same? k_nitr_h2o_vr(c,j) = w_scalar(c,j) ! nitrification constant is a set scalar * temp, moisture, and ph scalars k_nitr_vr(c,j) = k_nitr_max * 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 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)) ! limit to non-frozen soil layers if ( t_soisno(c,j) <= SHR_CONST_TKFRZ .and. no_frozen_nitrif_denitrif) then pot_f_nit_vr(c,j) = 0._r8 endif !---------------- denitrification ! first some input variables an unit conversions soil_hr_vr(c,j) = phr_vr(c,j) ! CENTURY papers give denitrification in units of per gram soil; need to convert from volumetric to mass-based units here soil_bulkdensity(c,j) = bd(c,j) + h2osoi_liq(c,j)/col%dz(c,j) g_per_m3__to__ug_per_gsoil = 1.e3_r8 / soil_bulkdensity(c,j) g_per_m3_sec__to__ug_per_gsoil_day = g_per_m3__to__ug_per_gsoil * secspday smin_no3_massdens_vr(c,j) = max(smin_no3_vr(c,j), 0._r8) * g_per_m3__to__ug_per_gsoil soil_co2_prod(c,j) = (soil_hr_vr(c,j) * (g_per_m3_sec__to__ug_per_gsoil_day)) !! maximum potential denitrification rates based on heterotrophic respiration rates or nitrate concentrations, !! from (del Grosso et al., 2000) fmax_denit_carbonsubstrate_vr(c,j) = (denit_resp_coef * (soil_co2_prod(c,j)**denit_resp_exp)) & / g_per_m3_sec__to__ug_per_gsoil_day ! fmax_denit_nitrate_vr(c,j) = (denit_nitrate_coef * smin_no3_massdens_vr(c,j)**denit_nitrate_exp) & / g_per_m3_sec__to__ug_per_gsoil_day ! find limiting denitrification rate f_denit_base_vr(c,j) = max(min(fmax_denit_carbonsubstrate_vr(c,j), fmax_denit_nitrate_vr(c,j)),0._r8) ! limit to non-frozen soil layers if ( t_soisno(c,j) <= SHR_CONST_TKFRZ .and. no_frozen_nitrif_denitrif ) then f_denit_base_vr(c,j) = 0._r8 endif ! limit to anoxic fraction of soils pot_f_denit_vr(c,j) = f_denit_base_vr(c,j) * anaerobic_frac(c,j) ! now calculate the ratio of N2O to N2 from denitrifictaion, following Del Grosso et al., 2000 ! diffusivity constant (figure 6b) ratio_k1(c,j) = max(1.7_r8, 38.4_r8 - 350._r8 * diffus(c,j)) ! ratio function (figure 7c) if ( soil_co2_prod(c,j) > 1.0e-9_r8 ) then ratio_no3_co2(c,j) = smin_no3_massdens_vr(c,j) / soil_co2_prod(c,j) else ! fucntion saturates at large no3/co2 ratios, so set as some nominally large number ratio_no3_co2(c,j) = 100._r8 endif ! total water limitation function (Del Grosso et al., 2000, figure 7a) 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_lch4) then if (anoxia_wtsat) then fr_WFPS(c,j) = fr_WFPS(c,j)*(1._r8 - finundated(c)) + finundated(c)*1.18_r8 end if end if ! final ratio expression 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) end do end do end associate end subroutine SoilBiogeochemNitrifDenitrif end module SoilBiogeochemNitrifDenitrifMod