module SoilBiogeochemDecompCascadeCNMod !----------------------------------------------------------------------- ! !DESCRIPTION: ! Sets the coeffiecients used in the decomposition cascade submodel. ! This uses the CN parameters as in CLMCN 4.0 ! ! !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 : nlevsoi, nlevgrnd, nlevdecomp, ndecomp_cascade_transitions, ndecomp_pools use clm_varpar , only : i_met_lit, i_cel_lit, i_lig_lit, i_cwd use clm_varctl , only : iulog, spinup_state, anoxia, use_lch4, use_vertsoilc, use_fates use clm_varcon , only : zsoi use decompMod , only : bounds_type use abortutils , only : endrun use CNSharedParamsMod , only : CNParamsShareInst, anoxia_wtsat, nlev_soildecomp_standard use SoilBiogeochemDecompCascadeConType , only : decomp_cascade_con use SoilBiogeochemStateType , only : soilbiogeochem_state_type use SoilBiogeochemCarbonFluxType , only : soilbiogeochem_carbonflux_type use SoilStateType , only : soilstate_type use CanopyStateType , only : canopystate_type use TemperatureType , only : temperature_type use ch4Mod , only : ch4_type use ColumnType , only : col ! implicit none private ! ! !PUBLIC MEMBER FUNCTIONS: public :: readParams public :: init_decompcascade_cn public :: decomp_rate_constants_cn type, private :: params_type real(r8):: cn_s1_cn !C:N for SOM 1 real(r8):: cn_s2_cn !C:N for SOM 2 real(r8):: cn_s3_cn !C:N for SOM 3 real(r8):: cn_s4_cn !C:N for SOM 4 real(r8):: rf_l1s1_cn !respiration fraction litter 1 -> SOM 1 real(r8):: rf_l2s2_cn !respiration fraction litter 2 -> SOM 2 real(r8):: rf_l3s3_cn !respiration fraction litter 3 -> SOM 3 real(r8):: rf_s1s2_cn !respiration fraction SOM 1 -> SOM 2 real(r8):: rf_s2s3_cn !respiration fraction SOM 2 -> SOM 3 real(r8):: rf_s3s4_cn !respiration fraction SOM 3 -> SOM 4 real(r8) :: cwd_fcel_cn !cellulose fraction for CWD real(r8) :: cwd_flig_cn ! real(r8) :: k_l1_cn !decomposition rate for litter 1 real(r8) :: k_l2_cn !decomposition rate for litter 2 real(r8) :: k_l3_cn !decomposition rate for litter 3 real(r8) :: k_s1_cn !decomposition rate for SOM 1 real(r8) :: k_s2_cn !decomposition rate for SOM 2 real(r8) :: k_s3_cn !decomposition rate for SOM 3 real(r8) :: k_s4_cn !decomposition rate for SOM 4 real(r8) :: k_frag_cn !fragmentation rate for CWD real(r8) :: minpsi_cn !minimum soil water potential for heterotrophic resp real(r8) :: maxpsi_cn !maximum soil water potential for heterotrophic resp integer :: nsompools = 4 real(r8), allocatable :: spinup_vector(:) ! multipliers for soil decomp during accelerated spinup end type params_type ! type(params_type), private :: params_inst character(len=*), parameter, private :: sourcefile = & __FILE__ !----------------------------------------------------------------------- contains !----------------------------------------------------------------------- subroutine readParams ( ncid ) ! ! !USES: use ncdio_pio , only : file_desc_t,ncd_io ! ! !ARGUMENTS: implicit none type(file_desc_t),intent(inout) :: ncid ! pio netCDF file id ! ! !CALLED FROM: readParamsMod.F90::CNParamsReadFile ! ! !REVISION HISTORY: ! Dec 3 2012 : Created by S. Muszala ! ! !LOCAL VARIABLES: character(len=32) :: subname = 'SoilBiogeochemDecompCnParamsType' 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 !EOP !----------------------------------------------------------------------- ! These are not read off of netcdf file allocate(params_inst%spinup_vector(params_inst%nsompools)) params_inst%spinup_vector(:) = (/ 1.0_r8, 1.0_r8, 5.0_r8, 70.0_r8 /) ! Read off of netcdf file tString='cn_s1' 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%cn_s1_cn=tempr tString='cn_s2' 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%cn_s2_cn=tempr tString='cn_s3' 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%cn_s3_cn=tempr tString='cn_s4' 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%cn_s4_cn=tempr tString='rf_l1s1' 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%rf_l1s1_cn=tempr tString='rf_l2s2' 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%rf_l2s2_cn=tempr tString='rf_l3s3' 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%rf_l3s3_cn=tempr tString='rf_s1s2' 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%rf_s1s2_cn=tempr tString='rf_s2s3' 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%rf_s2s3_cn=tempr tString='rf_s3s4' 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%rf_s3s4_cn=tempr tString='cwd_fcel' 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%cwd_fcel_cn=tempr tString='k_l1' 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%k_l1_cn=tempr tString='k_l2' 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%k_l2_cn=tempr tString='k_l3' 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%k_l3_cn=tempr tString='k_s1' 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%k_s1_cn=tempr tString='k_s2' 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%k_s2_cn=tempr tString='k_s3' 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%k_s3_cn=tempr tString='k_s4' 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%k_s4_cn=tempr tString='k_frag' 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%k_frag_cn=tempr tString='minpsi_hr' 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%minpsi_cn=tempr tString='maxpsi_hr' 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%maxpsi_cn=tempr tString='cwd_flig' 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%cwd_flig_cn=tempr end subroutine readParams !----------------------------------------------------------------------- subroutine init_decompcascade_cn(bounds, soilbiogeochem_state_inst) ! ! !DESCRIPTION: ! initialize rate constants and decomposition pathways for the BGC model originally implemented in CLM-CN ! written by C. Koven based on original CLM4 decomposition cascade by P. Thornton ! ! !ARGUMENTS: type(bounds_type) , intent(in) :: bounds type(soilbiogeochem_state_type) , intent(inout) :: soilbiogeochem_state_inst ! !-- properties of each pathway along decomposition cascade !-- properties of each decomposing pool real(r8) :: rf_l1s1 !respiration fraction litter 1 -> SOM 1 real(r8) :: rf_l2s2 !respiration fraction litter 2 -> SOM 2 real(r8) :: rf_l3s3 !respiration fraction litter 3 -> SOM 3 real(r8) :: rf_s1s2 !respiration fraction SOM 1 -> SOM 2 real(r8) :: rf_s2s3 !respiration fraction SOM 2 -> SOM 3 real(r8) :: rf_s3s4 !respiration fraction SOM 3 -> SOM 4 real(r8) :: cwd_fcel real(r8) :: cwd_flig real(r8) :: cn_s1 real(r8) :: cn_s2 real(r8) :: cn_s3 real(r8) :: cn_s4 integer :: i_litr1 integer :: i_litr2 integer :: i_litr3 integer :: i_soil1 integer :: i_soil2 integer :: i_soil3 integer :: i_soil4 integer :: i_atm integer :: i_l1s1 integer :: i_l2s2 integer :: i_l3s3 integer :: i_s1s2 integer :: i_s2s3 integer :: i_s3s4 integer :: i_s4atm integer :: i_cwdl2 integer :: i_cwdl3 !----------------------------------------------------------------------- associate( & rf_decomp_cascade => soilbiogeochem_state_inst%rf_decomp_cascade_col , & ! Output: [real(r8) (:,:,:) ] respired fraction in decomposition step (frac) pathfrac_decomp_cascade => soilbiogeochem_state_inst%pathfrac_decomp_cascade_col , & ! Output: [real(r8) (:,:,:) ] what fraction of C leaving a given pool passes through a given transition (frac) cascade_donor_pool => decomp_cascade_con%cascade_donor_pool , & ! Output: [integer (:) ] which pool is C taken from for a given decomposition step cascade_receiver_pool => decomp_cascade_con%cascade_receiver_pool , & ! Output: [integer (:) ] which pool is C added to for a given decomposition step floating_cn_ratio_decomp_pools => decomp_cascade_con%floating_cn_ratio_decomp_pools , & ! Output: [logical (:) ] TRUE => pool has fixed C:N ratio is_litter => decomp_cascade_con%is_litter , & ! Output: [logical (:) ] TRUE => pool is a litter pool is_soil => decomp_cascade_con%is_soil , & ! Output: [logical (:) ] TRUE => pool is a soil pool is_cwd => decomp_cascade_con%is_cwd , & ! Output: [logical (:) ] TRUE => pool is a cwd pool initial_cn_ratio => decomp_cascade_con%initial_cn_ratio , & ! Output: [real(r8) (:) ] c:n ratio for initialization of pools initial_stock => decomp_cascade_con%initial_stock , & ! Output: [real(r8) (:) ] initial concentration for seeding at spinup is_metabolic => decomp_cascade_con%is_metabolic , & ! Output: [logical (:) ] TRUE => pool is metabolic material is_cellulose => decomp_cascade_con%is_cellulose , & ! Output: [logical (:) ] TRUE => pool is cellulose is_lignin => decomp_cascade_con%is_lignin , & ! Output: [logical (:) ] TRUE => pool is lignin spinup_factor => decomp_cascade_con%spinup_factor & ! Output: [real(r8) (:) ] factor for AD spinup associated with each pool ) !------- time-constant coefficients ---------- ! ! set soil organic matter compartment C:N ratios (from Biome-BGC v4.2.0) cn_s1=params_inst%cn_s1_cn cn_s2=params_inst%cn_s2_cn cn_s3=params_inst%cn_s3_cn cn_s4=params_inst%cn_s4_cn ! set respiration fractions for fluxes between compartments ! (from Biome-BGC v4.2.0) rf_l1s1=params_inst%rf_l1s1_cn rf_l2s2=params_inst%rf_l2s2_cn rf_l3s3=params_inst%rf_l3s3_cn rf_s1s2=params_inst%rf_s1s2_cn rf_s2s3=params_inst%rf_s2s3_cn rf_s3s4=params_inst%rf_s3s4_cn ! set the cellulose and lignin fractions for coarse woody debris cwd_fcel=params_inst%cwd_fcel_cn cwd_flig=params_inst%cwd_flig_cn !------------------- list of pools and their attributes ------------ i_litr1 = i_met_lit floating_cn_ratio_decomp_pools(i_litr1) = .true. decomp_cascade_con%decomp_pool_name_restart(i_litr1) = 'litr1' decomp_cascade_con%decomp_pool_name_history(i_litr1) = 'LITR1' decomp_cascade_con%decomp_pool_name_long(i_litr1) = 'litter 1' decomp_cascade_con%decomp_pool_name_short(i_litr1) = 'L1' is_litter(i_litr1) = .true. is_soil(i_litr1) = .false. is_cwd(i_litr1) = .false. initial_cn_ratio(i_litr1) = 90._r8 initial_stock(i_litr1) = 0._r8 is_metabolic(i_litr1) = .true. is_cellulose(i_litr1) = .false. is_lignin(i_litr1) = .false. i_litr2 = i_cel_lit floating_cn_ratio_decomp_pools(i_litr2) = .true. decomp_cascade_con%decomp_pool_name_restart(i_litr2) = 'litr2' decomp_cascade_con%decomp_pool_name_history(i_litr2) = 'LITR2' decomp_cascade_con%decomp_pool_name_long(i_litr2) = 'litter 2' decomp_cascade_con%decomp_pool_name_short(i_litr2) = 'L2' is_litter(i_litr2) = .true. is_soil(i_litr2) = .false. is_cwd(i_litr2) = .false. initial_cn_ratio(i_litr2) = 90._r8 initial_stock(i_litr2) = 0._r8 is_metabolic(i_litr2) = .false. is_cellulose(i_litr2) = .true. is_lignin(i_litr2) = .false. i_litr3 = i_lig_lit floating_cn_ratio_decomp_pools(i_litr3) = .true. decomp_cascade_con%decomp_pool_name_restart(i_litr3) = 'litr3' decomp_cascade_con%decomp_pool_name_history(i_litr3) = 'LITR3' decomp_cascade_con%decomp_pool_name_long(i_litr3) = 'litter 3' decomp_cascade_con%decomp_pool_name_short(i_litr3) = 'L3' is_litter(i_litr3) = .true. is_soil(i_litr3) = .false. is_cwd(i_litr3) = .false. initial_cn_ratio(i_litr3) = 90._r8 initial_stock(i_litr3) = 0._r8 is_metabolic(i_litr3) = .false. is_cellulose(i_litr3) = .false. is_lignin(i_litr3) = .true. if (.not. use_fates) then floating_cn_ratio_decomp_pools(i_cwd) = .true. decomp_cascade_con%decomp_pool_name_restart(i_cwd) = 'cwd' decomp_cascade_con%decomp_pool_name_history(i_cwd) = 'CWD' decomp_cascade_con%decomp_pool_name_long(i_cwd) = 'coarse woody debris' decomp_cascade_con%decomp_pool_name_short(i_cwd) = 'CWD' is_litter(i_cwd) = .false. is_soil(i_cwd) = .false. is_cwd(i_cwd) = .true. initial_cn_ratio(i_cwd) = 500._r8 initial_stock(i_cwd) = 0._r8 is_metabolic(i_cwd) = .false. is_cellulose(i_cwd) = .false. is_lignin(i_cwd) = .false. end if if ( .not. use_fates ) then i_soil1 = 5 else i_soil1 = 4 endif floating_cn_ratio_decomp_pools(i_soil1) = .false. decomp_cascade_con%decomp_pool_name_restart(i_soil1) = 'soil1' decomp_cascade_con%decomp_pool_name_history(i_soil1) = 'SOIL1' decomp_cascade_con%decomp_pool_name_long(i_soil1) = 'soil 1' decomp_cascade_con%decomp_pool_name_short(i_soil1) = 'S1' is_litter(i_soil1) = .false. is_soil(i_soil1) = .true. is_cwd(i_soil1) = .false. initial_cn_ratio(i_soil1) = cn_s1 initial_stock(i_soil1) = 0._r8 is_metabolic(i_soil1) = .false. is_cellulose(i_soil1) = .false. is_lignin(i_soil1) = .false. if ( .not. use_fates ) then i_soil2 = 6 else i_soil2 = 5 endif floating_cn_ratio_decomp_pools(i_soil2) = .false. decomp_cascade_con%decomp_pool_name_restart(i_soil2) = 'soil2' decomp_cascade_con%decomp_pool_name_history(i_soil2) = 'SOIL2' decomp_cascade_con%decomp_pool_name_long(i_soil2) = 'soil 2' decomp_cascade_con%decomp_pool_name_short(i_soil2) = 'S2' is_litter(i_soil2) = .false. is_soil(i_soil2) = .true. is_cwd(i_soil2) = .false. initial_cn_ratio(i_soil2) = cn_s2 initial_stock(i_soil2) = 0._r8 is_metabolic(i_soil2) = .false. is_cellulose(i_soil2) = .false. is_lignin(i_soil2) = .false. if ( .not. use_fates ) then i_soil3 = 7 else i_soil3 = 6 endif floating_cn_ratio_decomp_pools(i_soil3) = .false. decomp_cascade_con%decomp_pool_name_restart(i_soil3) = 'soil3' decomp_cascade_con%decomp_pool_name_history(i_soil3) = 'SOIL3' decomp_cascade_con%decomp_pool_name_long(i_soil3) = 'soil 3' decomp_cascade_con%decomp_pool_name_short(i_soil3) = 'S3' is_litter(i_soil3) = .false. is_soil(i_soil3) = .true. is_cwd(i_soil3) = .false. initial_cn_ratio(i_soil3) = cn_s3 initial_stock(i_soil3) = 0._r8 is_metabolic(i_soil3) = .false. is_cellulose(i_soil3) = .false. is_lignin(i_soil3) = .false. if ( .not. use_fates ) then i_soil4 = 8 else i_soil4 = 7 endif floating_cn_ratio_decomp_pools(i_soil4) = .false. decomp_cascade_con%decomp_pool_name_restart(i_soil4) = 'soil4' decomp_cascade_con%decomp_pool_name_history(i_soil4) = 'SOIL4' decomp_cascade_con%decomp_pool_name_long(i_soil4) = 'soil 4' decomp_cascade_con%decomp_pool_name_short(i_soil4) = 'S4' is_litter(i_soil4) = .false. is_soil(i_soil4) = .true. is_cwd(i_soil4) = .false. initial_cn_ratio(i_soil4) = cn_s4 initial_stock(i_soil4) = 10._r8 is_metabolic(i_soil4) = .false. is_cellulose(i_soil4) = .false. is_lignin(i_soil4) = .false. i_atm = 0 !! for terminal pools (i.e. 100% respiration) floating_cn_ratio_decomp_pools(i_atm) = .false. decomp_cascade_con%decomp_pool_name_restart(i_atm) = 'atmosphere' decomp_cascade_con%decomp_pool_name_history(i_atm) = 'atmosphere' decomp_cascade_con%decomp_pool_name_long(i_atm) = 'atmosphere' decomp_cascade_con%decomp_pool_name_short(i_atm) = '' is_litter(i_atm) = .true. is_soil(i_atm) = .false. is_cwd(i_atm) = .false. initial_cn_ratio(i_atm) = 0._r8 initial_stock(i_atm) = 0._r8 is_metabolic(i_atm) = .false. is_cellulose(i_atm) = .false. is_lignin(i_atm) = .false. spinup_factor(i_litr1) = 1._r8 spinup_factor(i_litr2) = 1._r8 spinup_factor(i_litr3) = 1._r8 if (.not. use_fates) then spinup_factor(i_cwd) = 1._r8 end if spinup_factor(i_soil1) = params_inst%spinup_vector(1) spinup_factor(i_soil2) = params_inst%spinup_vector(2) spinup_factor(i_soil3) = params_inst%spinup_vector(3) spinup_factor(i_soil4) = params_inst%spinup_vector(4) !---------------- list of transitions and their time-independent coefficients ---------------! i_l1s1 = 1 decomp_cascade_con%cascade_step_name(i_l1s1) = 'L1S1' rf_decomp_cascade(bounds%begc:bounds%endc,1:nlevdecomp,i_l1s1) = rf_l1s1 cascade_donor_pool(i_l1s1) = i_litr1 cascade_receiver_pool(i_l1s1) = i_soil1 pathfrac_decomp_cascade(bounds%begc:bounds%endc,1:nlevdecomp,i_l1s1) = 1.0_r8 i_l2s2 = 2 decomp_cascade_con%cascade_step_name(i_l2s2) = 'L2S2' rf_decomp_cascade(bounds%begc:bounds%endc,1:nlevdecomp,i_l2s2) = rf_l2s2 cascade_donor_pool(i_l2s2) = i_litr2 cascade_receiver_pool(i_l2s2) = i_soil2 pathfrac_decomp_cascade(bounds%begc:bounds%endc,1:nlevdecomp,i_l2s2) = 1.0_r8 i_l3s3 = 3 decomp_cascade_con%cascade_step_name(i_l3s3) = 'L3S3' rf_decomp_cascade(bounds%begc:bounds%endc,1:nlevdecomp,i_l3s3) = rf_l3s3 cascade_donor_pool(i_l3s3) = i_litr3 cascade_receiver_pool(i_l3s3) = i_soil3 pathfrac_decomp_cascade(bounds%begc:bounds%endc,1:nlevdecomp,i_l3s3) = 1.0_r8 i_s1s2 = 4 decomp_cascade_con%cascade_step_name(i_s1s2) = 'S1S2' rf_decomp_cascade(bounds%begc:bounds%endc,1:nlevdecomp,i_s1s2) = rf_s1s2 cascade_donor_pool(i_s1s2) = i_soil1 cascade_receiver_pool(i_s1s2) = i_soil2 pathfrac_decomp_cascade(bounds%begc:bounds%endc,1:nlevdecomp,i_s1s2) = 1.0_r8 i_s2s3 = 5 decomp_cascade_con%cascade_step_name(i_s2s3) = 'S2S3' rf_decomp_cascade(bounds%begc:bounds%endc,1:nlevdecomp,i_s2s3) = rf_s2s3 cascade_donor_pool(i_s2s3) = i_soil2 cascade_receiver_pool(i_s2s3) = i_soil3 pathfrac_decomp_cascade(bounds%begc:bounds%endc,1:nlevdecomp,i_s2s3) = 1.0_r8 i_s3s4 = 6 decomp_cascade_con%cascade_step_name(i_s3s4) = 'S3S4' rf_decomp_cascade(bounds%begc:bounds%endc,1:nlevdecomp,i_s3s4) = rf_s3s4 cascade_donor_pool(i_s3s4) = i_soil3 cascade_receiver_pool(i_s3s4) = i_soil4 pathfrac_decomp_cascade(bounds%begc:bounds%endc,1:nlevdecomp,i_s3s4) = 1.0_r8 i_s4atm = 7 decomp_cascade_con%cascade_step_name(i_s4atm) = 'S4' rf_decomp_cascade(bounds%begc:bounds%endc,1:nlevdecomp,i_s4atm) = 1. cascade_donor_pool(i_s4atm) = i_soil4 cascade_receiver_pool(i_s4atm) = i_atm pathfrac_decomp_cascade(bounds%begc:bounds%endc,1:nlevdecomp,i_s4atm) = 1.0_r8 if (.not. use_fates) then i_cwdl2 = 8 decomp_cascade_con%cascade_step_name(i_cwdl2) = 'CWDL2' rf_decomp_cascade(bounds%begc:bounds%endc,1:nlevdecomp,i_cwdl2) = 0._r8 cascade_donor_pool(i_cwdl2) = i_cwd cascade_receiver_pool(i_cwdl2) = i_litr2 pathfrac_decomp_cascade(bounds%begc:bounds%endc,1:nlevdecomp,i_cwdl2) = cwd_fcel i_cwdl3 = 9 decomp_cascade_con%cascade_step_name(i_cwdl3) = 'CWDL3' rf_decomp_cascade(bounds%begc:bounds%endc,1:nlevdecomp,i_cwdl3) = 0._r8 cascade_donor_pool(i_cwdl3) = i_cwd cascade_receiver_pool(i_cwdl3) = i_litr3 pathfrac_decomp_cascade(bounds%begc:bounds%endc,1:nlevdecomp,i_cwdl3) = cwd_flig end if end associate end subroutine init_decompcascade_cn !----------------------------------------------------------------------- subroutine decomp_rate_constants_cn(bounds, & num_soilc, filter_soilc, & canopystate_inst, soilstate_inst, temperature_inst, ch4_inst, soilbiogeochem_carbonflux_inst) ! ! !DESCRIPTION: ! calculate rate constants and decomposition pathways for the BGC model ! originally implemented in CLM-CN ! written by C. Koven based on original CLM4 decomposition cascade by P. Thornton ! ! !USES: use clm_time_manager, only : get_step_size use clm_varcon , only : secspday use clm_varpar , only : i_cwd ! ! !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(canopystate_type) , intent(in) :: canopystate_inst type(soilstate_type) , intent(in) :: soilstate_inst type(temperature_type) , intent(in) :: temperature_inst type(ch4_type) , intent(in) :: ch4_inst type(soilbiogeochem_carbonflux_type) , intent(inout) :: soilbiogeochem_carbonflux_inst ! ! !LOCAL VARIABLES: real(r8):: dt ! decomp timestep (seconds) real(r8):: dtd ! decomp timestep (days) real(r8):: frw(bounds%begc:bounds%endc) ! rooting fraction weight real(r8), allocatable:: fr(:,:) ! column-level rooting fraction by soil depth real(r8):: minpsi, maxpsi ! limits for soil water scalar for decomp real(r8):: psi ! temporary soilpsi for water scalar real(r8):: rate_scalar ! combined rate scalar for decomp real(r8):: k_l1 ! decomposition rate constant litter 1 real(r8):: k_l2 ! decomposition rate constant litter 2 real(r8):: k_l3 ! decomposition rate constant litter 3 real(r8):: k_s1 ! decomposition rate constant SOM 1 real(r8):: k_s2 ! decomposition rate constant SOM 2 real(r8):: k_s3 ! decomposition rate constant SOM 3 real(r8):: k_s4 ! decomposition rate constant SOM 4 real(r8):: k_frag ! fragmentation rate constant CWD real(r8):: ck_l1 ! corrected decomposition rate constant litter 1 real(r8):: ck_l2 ! corrected decomposition rate constant litter 2 real(r8):: ck_l3 ! corrected decomposition rate constant litter 3 real(r8):: ck_s1 ! corrected decomposition rate constant SOM 1 real(r8):: ck_s2 ! corrected decomposition rate constant SOM 2 real(r8):: ck_s3 ! corrected decomposition rate constant SOM 3 real(r8):: ck_s4 ! corrected decomposition rate constant SOM 4 real(r8):: ck_frag ! corrected fragmentation rate constant CWD real(r8):: cwdc_loss ! fragmentation rate for CWD carbon (gC/m2/s) real(r8):: cwdn_loss ! fragmentation rate for CWD nitrogen (gN/m2/s) integer :: i_litr1 integer :: i_litr2 integer :: i_litr3 integer :: i_soil1 integer :: i_soil2 integer :: i_soil3 integer :: i_soil4 integer :: c, fc, j, k, l real(r8):: Q10 ! temperature dependence real(r8):: froz_q10 ! separate q10 for frozen soil respiration rates. default to same as above zero rates real(r8):: decomp_depth_efolding ! (meters) e-folding depth for reduction in decomposition [ real(r8):: depth_scalar(bounds%begc:bounds%endc,1:nlevdecomp) real(r8) :: mino2lim ! minimum anaerobic decomposition rate as a ! fraction of potential aerobic rate !----------------------------------------------------------------------- associate( & dz => col%dz , & ! Input: [real(r8) (:,:) ] soil layer thickness (m) soilpsi => soilstate_inst%soilpsi_col , & ! Input: [real(r8) (:,:) ] soil water potential in each soil layer (MPa) alt_indx => canopystate_inst%alt_indx_col , & ! Input: [integer (:) ] current depth of thaw t_soisno => temperature_inst%t_soisno_col , & ! Input: [real(r8) (:,:) ] soil temperature (Kelvin) (-nlevsno+1:nlevgrnd) o2stress_sat => ch4_inst%o2stress_sat_col , & ! Input: [real(r8) (:,:) ] Ratio of oxygen available to that demanded by roots, aerobes, & methanotrophs (nlevsoi) o2stress_unsat => ch4_inst%o2stress_unsat_col , & ! Input: [real(r8) (:,:) ] Ratio of oxygen available to that demanded by roots, aerobes, & methanotrophs (nlevsoi) finundated => ch4_inst%finundated_col , & ! Input: [real(r8) (:) ] fractional inundated area (excluding dedicated wetland columns) t_scalar => soilbiogeochem_carbonflux_inst%t_scalar_col , & ! Output: [real(r8) (:,:) ] soil temperature scalar for decomp w_scalar => soilbiogeochem_carbonflux_inst%w_scalar_col , & ! Output: [real(r8) (:,:) ] soil water scalar for decomp o_scalar => soilbiogeochem_carbonflux_inst%o_scalar_col , & ! Output: [real(r8) (:,:) ] fraction by which decomposition is limited by anoxia decomp_k => soilbiogeochem_carbonflux_inst%decomp_k_col & ! Output: [real(r8) (:,:,:) ] rate constant for decomposition (1./sec) ) mino2lim = CNParamsShareInst%mino2lim ! set time steps dt = real( get_step_size(), r8 ) dtd = dt/secspday ! set initial base rates for decomposition mass loss (1/day) ! (from Biome-BGC v4.2.0, using three SOM pools) ! Value inside log function is the discrete-time values for a ! daily time step model, and the result of the log function is ! the corresponding continuous-time decay rate (1/day), following ! Olson, 1963. k_l1=params_inst%k_l1_cn k_l2=params_inst%k_l2_cn k_l3=params_inst%k_l3_cn k_s1=params_inst%k_s1_cn k_s2=params_inst%k_s2_cn k_s3=params_inst%k_s3_cn k_s4=params_inst%k_s4_cn k_frag=params_inst%k_frag_cn ! calculate the new discrete-time decay rate for model timestep k_l1 = 1.0_r8-exp(-k_l1*dtd) k_l2 = 1.0_r8-exp(-k_l2*dtd) k_l3 = 1.0_r8-exp(-k_l3*dtd) k_s1 = 1.0_r8-exp(-k_s1*dtd) k_s2 = 1.0_r8-exp(-k_s2*dtd) k_s3 = 1.0_r8-exp(-k_s3*dtd) k_s4 = 1.0_r8-exp(-k_s4*dtd) k_frag = 1.0_r8-exp(-k_frag*dtd) minpsi = params_inst%minpsi_cn maxpsi = params_inst%maxpsi_cn Q10 = CNParamsShareInst%Q10 ! set "froz_q10" parameter froz_q10 = CNParamsShareInst%froz_q10 if (use_vertsoilc) then ! Set "decomp_depth_efolding" parameter decomp_depth_efolding = CNParamsShareInst%decomp_depth_efolding end if ! The following code implements the acceleration part of the AD spinup ! algorithm, by multiplying all of the SOM decomposition base rates by 10.0. if ( spinup_state .eq. 1 ) then k_s1 = k_s1 * params_inst%spinup_vector(1) k_s2 = k_s2 * params_inst%spinup_vector(2) k_s3 = k_s3 * params_inst%spinup_vector(3) k_s4 = k_s4 * params_inst%spinup_vector(4) endif i_litr1 = 1 i_litr2 = 2 i_litr3 = 3 if (use_fates) then i_soil1 = 4 i_soil2 = 5 i_soil3 = 6 i_soil4 = 7 else i_soil1 = 5 i_soil2 = 6 i_soil3 = 7 i_soil4 = 8 endif !--- time dependent coefficients-----! if ( nlevdecomp .eq. 1 ) then ! calculate function to weight the temperature and water potential scalars ! for decomposition control. ! the following normalizes values in fr so that they ! sum to 1.0 across top nlevdecomp levels on a column frw(bounds%begc:bounds%endc) = 0._r8 nlev_soildecomp_standard=5 allocate(fr(bounds%begc:bounds%endc,nlev_soildecomp_standard)) do j=1,nlev_soildecomp_standard do fc = 1,num_soilc c = filter_soilc(fc) frw(c) = frw(c) + dz(c,j) end do end do do j = 1,nlev_soildecomp_standard do fc = 1,num_soilc c = filter_soilc(fc) if (frw(c) /= 0._r8) then fr(c,j) = dz(c,j) / frw(c) else fr(c,j) = 0._r8 end if end do end do ! calculate rate constant scalar for soil temperature ! assuming that the base rate constants are assigned for non-moisture ! limiting conditions at 25 C. ! Peter Thornton: 3/13/09 ! Replaced the Lloyd and Taylor function with a Q10 formula, with Q10 = 1.5 ! as part of the modifications made to improve the seasonal cycle of ! atmospheric CO2 concentration in global simulations. This does not impact ! the base rates at 25 C, which are calibrated from microcosm studies. do j = 1,nlev_soildecomp_standard do fc = 1,num_soilc c = filter_soilc(fc) if (j==1) t_scalar(c,:) = 0._r8 !! use separate (possibly equal) t funcs above and below freezing point !! t_scalar(c,1)=t_scalar(c,1) + (1.5**((t_soisno(c,j)-(SHR_CONST_TKFRZ+25._r8))/10._r8))*fr(c,j) if (t_soisno(c,j) >= SHR_CONST_TKFRZ) then t_scalar(c,1)=t_scalar(c,1) + & (Q10**((t_soisno(c,j)-(SHR_CONST_TKFRZ+25._r8))/10._r8))*fr(c,j) else t_scalar(c,1)=t_scalar(c,1) + & (Q10**(-25._r8/10._r8))*(froz_q10**((t_soisno(c,j)-SHR_CONST_TKFRZ)/10._r8))*fr(c,j) endif end do end do ! calculate the rate constant scalar for soil water content. ! Uses the log relationship with water potential given in ! Andren, O., and K. Paustian, 1987. Barley straw decomposition in the field: ! a comparison of models. Ecology, 68(5):1190-1200. ! and supported by data in ! Orchard, V.A., and F.J. Cook, 1983. Relationship between soil respiration ! and soil moisture. Soil Biol. Biochem., 15(4):447-453. do j = 1,nlev_soildecomp_standard do fc = 1,num_soilc c = filter_soilc(fc) if (j==1) w_scalar(c,:) = 0._r8 psi = min(soilpsi(c,j),maxpsi) ! decomp only if soilpsi is higher than minpsi if (psi > minpsi) then w_scalar(c,1) = w_scalar(c,1) + (log(minpsi/psi)/log(minpsi/maxpsi))*fr(c,j) end if end do end do if (use_lch4) then if (anoxia_wtsat) then ! Adjust for saturated fraction if unfrozen. do fc = 1,num_soilc c = filter_soilc(fc) if (alt_indx(c) >= nlev_soildecomp_standard .and. t_soisno(c,1) > SHR_CONST_TKFRZ) then w_scalar(c,1) = w_scalar(c,1)*(1._r8 - finundated(c)) + finundated(c) end if end do end if end if if (use_lch4) then ! Calculate ANOXIA if (anoxia) then ! Check for anoxia w/o LCH4 now done in controlMod. do j = 1,nlev_soildecomp_standard do fc = 1,num_soilc c = filter_soilc(fc) if (j==1) o_scalar(c,:) = 0._r8 if (.not. anoxia_wtsat) then o_scalar(c,1) = o_scalar(c,1) + fr(c,j) * max(o2stress_unsat(c,j), mino2lim) else o_scalar(c,1) = o_scalar(c,1) + fr(c,j) * & (max(o2stress_unsat(c,j), mino2lim)*(1._r8 - finundated(c)) + & max(o2stress_sat(c,j), mino2lim)*finundated(c) ) end if end do end do else o_scalar(bounds%begc:bounds%endc,1:nlevdecomp) = 1._r8 end if else o_scalar(bounds%begc:bounds%endc,1:nlevdecomp) = 1._r8 end if deallocate(fr) else ! calculate rate constant scalar for soil temperature ! assuming that the base rate constants are assigned for non-moisture ! limiting conditions at 25 C. ! Peter Thornton: 3/13/09 ! Replaced the Lloyd and Taylor function with a Q10 formula, with Q10 = 1.5 ! as part of the modifications made to improve the seasonal cycle of ! atmospheric CO2 concentration in global simulations. This does not impact ! the base rates at 25 C, which are calibrated from microcosm studies. do j = 1, nlevdecomp do fc = 1,num_soilc c = filter_soilc(fc) !! use separate (possibly equal) t funcs above and below freezing point !! t_scalar(c,j)= (1.5**((t_soisno(c,j)-(SHR_CONST_TKFRZ+25._r8))/10._r8)) if (t_soisno(c,j) >= SHR_CONST_TKFRZ) then t_scalar(c,j)= (Q10**((t_soisno(c,j)-(SHR_CONST_TKFRZ+25._r8))/10._r8)) else t_scalar(c,j)= (Q10**(-25._r8/10._r8))*(froz_q10**((t_soisno(c,j)-SHR_CONST_TKFRZ)/10._r8)) endif end do end do ! calculate the rate constant scalar for soil water content. ! Uses the log relationship with water potential given in ! Andren, O., and K. Paustian, 1987. Barley straw decomposition in the field: ! a comparison of models. Ecology, 68(5):1190-1200. ! and supported by data in ! Orchard, V.A., and F.J. Cook, 1983. Relationship between soil respiration ! and soil moisture. Soil Biol. Biochem., 15(4):447-453. do j = 1,nlevdecomp do fc = 1,num_soilc c = filter_soilc(fc) psi = min(soilpsi(c,j),maxpsi) ! decomp only if soilpsi is higher than minpsi if (psi > minpsi) then w_scalar(c,j) = (log(minpsi/psi)/log(minpsi/maxpsi)) else w_scalar(c,j) = 0._r8 end if if (use_lch4) then if (anoxia_wtsat .and. t_soisno(c,j) > SHR_CONST_TKFRZ) then ! wet area will have w_scalar of 1 if unfrozen w_scalar(c,j) = w_scalar(c,j)*(1._r8 - finundated(c)) + finundated(c) end if end if end do end do end if if (use_lch4) then ! Calculate ANOXIA ! Check for anoxia w/o LCH4 now done in controlMod. if (anoxia) then do j = 1,nlevdecomp do fc = 1,num_soilc c = filter_soilc(fc) if (.not. anoxia_wtsat) then o_scalar(c,j) = max(o2stress_unsat(c,j), mino2lim) else o_scalar(c,j) = max(o2stress_unsat(c,j), mino2lim) * (1._r8 - finundated(c)) + & max(o2stress_sat(c,j), mino2lim) * finundated(c) end if end do end do else o_scalar(bounds%begc:bounds%endc,1:nlevdecomp) = 1._r8 end if else o_scalar(bounds%begc:bounds%endc,1:nlevdecomp) = 1._r8 end if if (use_vertsoilc) then ! add a term to reduce decomposition rate at depth ! for now used a fixed e-folding depth do j = 1, nlevdecomp do fc = 1,num_soilc c = filter_soilc(fc) depth_scalar(c,j) = exp(-zsoi(j)/decomp_depth_efolding) end do end do end if ! calculate rate constants for all litter and som pools if (use_vertsoilc) then do j = 1,nlevdecomp do fc = 1,num_soilc c = filter_soilc(fc) decomp_k(c,j,i_litr1) = k_l1 * t_scalar(c,j) * w_scalar(c,j) * depth_scalar(c,j) * o_scalar(c,j) / dt decomp_k(c,j,i_litr2) = k_l2 * t_scalar(c,j) * w_scalar(c,j) * depth_scalar(c,j) * o_scalar(c,j) / dt decomp_k(c,j,i_litr3) = k_l3 * t_scalar(c,j) * w_scalar(c,j) * depth_scalar(c,j) * o_scalar(c,j) / dt decomp_k(c,j,i_soil1) = k_s1 * t_scalar(c,j) * w_scalar(c,j) * depth_scalar(c,j) * o_scalar(c,j) / dt decomp_k(c,j,i_soil2) = k_s2 * t_scalar(c,j) * w_scalar(c,j) * depth_scalar(c,j) * o_scalar(c,j) / dt decomp_k(c,j,i_soil3) = k_s3 * t_scalar(c,j) * w_scalar(c,j) * depth_scalar(c,j) * o_scalar(c,j) / dt decomp_k(c,j,i_soil4) = k_s4 * t_scalar(c,j) * w_scalar(c,j) * depth_scalar(c,j) * o_scalar(c,j) / dt end do end do else do j = 1,nlevdecomp do fc = 1,num_soilc c = filter_soilc(fc) decomp_k(c,j,i_litr1) = k_l1 * t_scalar(c,j) * w_scalar(c,j) * o_scalar(c,j) / dt decomp_k(c,j,i_litr2) = k_l2 * t_scalar(c,j) * w_scalar(c,j) * o_scalar(c,j) / dt decomp_k(c,j,i_litr3) = k_l3 * t_scalar(c,j) * w_scalar(c,j) * o_scalar(c,j) / dt decomp_k(c,j,i_soil1) = k_s1 * t_scalar(c,j) * w_scalar(c,j) * o_scalar(c,j) / dt decomp_k(c,j,i_soil2) = k_s2 * t_scalar(c,j) * w_scalar(c,j) * o_scalar(c,j) / dt decomp_k(c,j,i_soil3) = k_s3 * t_scalar(c,j) * w_scalar(c,j) * o_scalar(c,j) / dt decomp_k(c,j,i_soil4) = k_s4 * t_scalar(c,j) * w_scalar(c,j) * o_scalar(c,j) / dt end do end do end if ! do the same for cwd, but only if fates is not enabled (because fates handles CWD on its own structure if (.not. use_fates) then if (use_vertsoilc) then do j = 1,nlevdecomp do fc = 1,num_soilc c = filter_soilc(fc) decomp_k(c,j,i_cwd) = k_frag * t_scalar(c,j) * w_scalar(c,j) * depth_scalar(c,j) * o_scalar(c,j) / dt end do end do else do j = 1,nlevdecomp do fc = 1,num_soilc c = filter_soilc(fc) decomp_k(c,j,i_cwd) = k_frag * t_scalar(c,j) * w_scalar(c,j) * o_scalar(c,j) / dt end do end do end if end if end associate end subroutine decomp_rate_constants_cn end module SoilBiogeochemDecompCascadeCNMod