Skip to content

Commit

Permalink
Merge pull request ESCOMP#2869 from wwieder/BNF_v2
Browse files Browse the repository at this point in the history 
Update temperature cost function for symbiotic Nfix in FUN
  • Loading branch information
@slevis-lmwg
slevis-lmwg authored Dec 16, 2024
2 parents 7417668 + 9b0be80 commit 2db5c6a
Show file tree
Hide file tree
Showing 7 changed files with 111 additions and 24 deletions.
10 changes: 9 additions & 1 deletion bld/CLMBuildNamelist.pm
View file
Original file line number Diff line number Diff line change
Expand Up @@ -3381,7 +3381,7 @@ sub setup_logic_mineral_nitrogen_dynamics {
#
my ($opts, $nl_flags, $definition, $defaults, $nl) = @_;

my @vars = ( "freelivfix_slope_wet", "freelivfix_intercept" );
my @vars = ( "freelivfix_slope_wet", "freelivfix_intercept", "nfix_method" );
if ( &value_is_true($nl_flags->{'use_cn'}) && &value_is_true($nl->get_value('use_fun')) ) {
foreach my $var ( @vars ) {
add_default($opts, $nl_flags->{'inputdata_rootdir'}, $definition, $defaults, $nl, $var,
Expand All @@ -3394,6 +3394,14 @@ sub setup_logic_mineral_nitrogen_dynamics {
}
}
}

if ( &value_is_true($nl_flags->{'use_cn'}) && &value_is_true($nl->get_value('use_fun')) ) {
my $var = $nl->get_value('nfix_method');
if ( $var ne "'Houlton'" && $var ne "'Bytnerowicz'" ) {
$log->fatal_error("$var is incorrect entry for the namelist variable nfix_method; expected Houlton or Bytnerowicz");
}
}

}


Expand Down
10 changes: 6 additions & 4 deletions bld/namelist_files/namelist_defaults_ctsm.xml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -544,10 +544,10 @@ attributes from the config_cache.xml file (with keys converted to upper-case).
<!-- The default filenames are given relative to the root directory
for the CLM2 data in the CESM distribution -->
<!-- Plant function types (relative to {csmdata}) -->
<paramfile phys="clm6_0" >lnd/clm2/paramdata/ctsm60_params.c241017.nc</paramfile>
<paramfile phys="clm5_1" >lnd/clm2/paramdata/ctsm51_params.c241017.nc</paramfile>
<paramfile phys="clm5_0" >lnd/clm2/paramdata/clm50_params.c241017.nc</paramfile>
<paramfile phys="clm4_5" >lnd/clm2/paramdata/clm45_params.c241017.nc</paramfile>
<paramfile phys="clm6_0" >lnd/clm2/paramdata/ctsm60_params.c241119.nc</paramfile>
<paramfile phys="clm5_1" >lnd/clm2/paramdata/ctsm60_params.c241119.nc</paramfile>
<paramfile phys="clm5_0" >lnd/clm2/paramdata/clm50_params.c241119.nc</paramfile>
<paramfile phys="clm4_5" >lnd/clm2/paramdata/clm45_params.c241119.nc</paramfile>

<!-- ================================================================== -->
<!-- FATES default parameter file -->
Expand Down Expand Up @@ -2024,6 +2024,8 @@ lnd/clm2/surfdata_esmf/NEON/ctsm5.3.0/surfdata_1x1_NEON_TOOL_hist_2000_78pfts_c2
<freelivfix_intercept use_fun=".true." use_cn=".true.">0.0117d00</freelivfix_intercept>
<freelivfix_slope_wet use_fun=".true." use_cn=".true.">0.0006d00</freelivfix_slope_wet>

<nfix_method>Houlton</nfix_method>

<br_root>0.83d-06</br_root>

<!-- Scalar of leaf respiration to vcmax (used for SP mode and with luna) (could vary with physics version) -->
Expand Down
5 changes: 5 additions & 0 deletions bld/namelist_files/namelist_definition_ctsm.xml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -388,6 +388,11 @@ Slope of free living Nitrogen fixation with annual ET
Intercept of free living Nitrogen fixation with zero annual ET
</entry>

<entry id="nfix_method" type="char*25" category="clm_nitrogen"
group="clm_nitrogen" valid_values="Houlton,Bytnerowicz" value="Houlton" >
Choice of nfix parameterization
</entry>

<entry id="use_undercanopy_stability" type="logical" category="clm_physics"
group="canopyfluxes_inparm" valid_values="" >
If TRUE use the undercanopy stability term used with CLM4.5 (Sakaguchi&amp;Zeng, 2008)
Expand Down
87 changes: 68 additions & 19 deletions src/biogeochem/CNFUNMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -26,7 +26,7 @@ module CNFUNMod
use pftconMod , only : pftcon, npcropmin
use decompMod , only : bounds_type
use clm_varctl , only : use_nitrif_denitrif,use_flexiblecn
use CNSharedParamsMod , only : use_matrixcn
use CNSharedParamsMod , only : use_matrixcn
use abortutils , only : endrun
use CNVegstateType , only : cnveg_state_type
use CNVegCarbonStateType , only : cnveg_carbonstate_type
Expand Down Expand Up @@ -206,7 +206,7 @@ subroutine CNFUN(bounds,num_soilc, filter_soilc,num_soilp&
use clm_time_manager, only : get_step_size_real, get_curr_date
use clm_varpar , only : nlevdecomp
use clm_varcon , only : secspday, smallValue, fun_period, tfrz, dzsoi_decomp, spval
use clm_varctl , only : use_nitrif_denitrif
use clm_varctl , only : use_nitrif_denitrif, nfix_method
use PatchType , only : patch
use subgridAveMod , only : p2c
use pftconMod , only : npcropmin
Expand Down Expand Up @@ -290,7 +290,7 @@ subroutine CNFUN(bounds,num_soilc, filter_soilc,num_soilp&
real(r8) :: litterfall_n(bounds%begp:bounds%endp) ! N loss based on the leafc to litter (gN/m2)
real(r8) :: litterfall_n_step(bounds%begp:bounds%endp,1:nstp) ! N loss based on the leafc to litter (gN/m2)
real(r8) :: litterfall_c_step(bounds%begp:bounds%endp,1:nstp) ! N loss based on the leafc to litter (gN/m2)
real(r8) :: tc_soisno(bounds%begc:bounds%endc,1:nlevdecomp) ! Soil temperature (degrees Celsius)
real(r8) :: tc_soisno(bounds%begc:bounds%endc,1:nlevdecomp) ! Soil temperature (degrees Celsius)
real(r8) :: npp_remaining(bounds%begp:bounds%endp,1:nstp) ! A temporary variable for npp_remaining(gC/m2)
real(r8) :: n_passive_step(bounds%begp:bounds%endp,1:nstp) ! N taken up by transpiration at substep(gN/m2)
real(r8) :: n_passive_acc(bounds%begp:bounds%endp) ! N acquired by passive uptake (gN/m2)
Expand Down Expand Up @@ -467,7 +467,6 @@ subroutine CNFUN(bounds,num_soilc, filter_soilc,num_soilp&
real(r8) :: total_c_spent_retrans
real(r8) :: total_c_accounted_retrans


!------end of not_use_nitrif_denitrif------!
!--------------------------------------------------------------------
!------------
Expand Down Expand Up @@ -495,9 +494,10 @@ subroutine CNFUN(bounds,num_soilc, filter_soilc,num_soilp&
! fixers, 2 for non fixers. This will become redundant with the
! 'fixer' parameter if it works.

character(len=32) :: subname = 'CNFUN'
!--------------------------------------------------------------------
!---------------------------------
associate(ivt => patch%itype , & ! Input: [integer (:) ] p
associate(ivt => patch%itype , & ! Input: [integer (:) ] p
leafcn => pftcon%leafcn , & ! Input: leaf C:N (gC/gN)
season_decid => pftcon%season_decid , & ! Input: binary flag for seasonal
! -deciduous leaf habit (0 or 1)
Expand All @@ -507,6 +507,9 @@ subroutine CNFUN(bounds,num_soilc, filter_soilc,num_soilp&
b_fix => pftcon%b_fix , & ! Input: A BNF parameter
c_fix => pftcon%c_fix , & ! Input: A BNF parameter
s_fix => pftcon%s_fix , & ! Input: A BNF parameter
nfix_tmin => pftcon%nfix_tmin , & ! Input: A BNF parameter
nfix_topt => pftcon%nfix_topt , & ! Input: A BNF parameter
nfix_tmax => pftcon%nfix_tmax , & ! Input: A BNF parameter
akc_active => pftcon%akc_active , & ! Input: A mycorrhizal uptake
! parameter
akn_active => pftcon%akn_active , & ! Input: A mycorrhizal uptake
Expand All @@ -522,10 +525,10 @@ subroutine CNFUN(bounds,num_soilc, filter_soilc,num_soilp&
perecm => pftcon%perecm , & ! Input: The fraction of ECM
! -associated PFT
grperc => pftcon%grperc , & ! Input: growth percentage
fun_cn_flex_a => pftcon%fun_cn_flex_a , & ! Parameter a of FUN-flexcn link code (def 5)
fun_cn_flex_b => pftcon%fun_cn_flex_b , & ! Parameter b of FUN-flexcn link code (def 200)
fun_cn_flex_c => pftcon%fun_cn_flex_c , & ! Parameter b of FUN-flexcn link code (def 80)
FUN_fracfixers => pftcon%FUN_fracfixers , & ! Fraction of C that can be used for fixation.
fun_cn_flex_a => pftcon%fun_cn_flex_a , & ! Parameter a of FUN-flexcn link code (def 5)
fun_cn_flex_b => pftcon%fun_cn_flex_b , & ! Parameter b of FUN-flexcn link code (def 200)
fun_cn_flex_c => pftcon%fun_cn_flex_c , & ! Parameter b of FUN-flexcn link code (def 80)
FUN_fracfixers => pftcon%FUN_fracfixers , & ! Fraction of C that can be used for fixation.
leafcn_offset => cnveg_state_inst%leafcn_offset_patch , & ! Output:
! [real(r8) (:)] Leaf C:N used by FUN
plantCN => cnveg_state_inst%plantCN_patch , & ! Output: [real(r8) (:)] Plant
Expand Down Expand Up @@ -1041,42 +1044,50 @@ subroutine CNFUN(bounds,num_soilc, filter_soilc,num_soilp&
npp_to_nonmyc_nh4(:) = 0.0_r8
npp_to_fixation(:) = 0.0_r8
npp_to_retrans(:) = 0.0_r8




unmetDemand = .TRUE.
plant_ndemand_pool_step(p,istp) = plant_ndemand_pool(p) * permyc(p,istp)
npp_remaining(p,istp) = availc_pool(p) * permyc(p,istp)


! if (plant_ndemand_pool_step(p,istp) .gt. 0._r8) then !
! plant_ndemand_pool_step > 0.0

do j = 1, nlevdecomp
tc_soisno(c,j) = t_soisno(c,j) - tfrz

if(pftcon%c3psn(patch%itype(p)).eq.1)then
fixer=1
else
fixer=0
endif
costNit(j,icostFix) = fun_cost_fix(fixer,a_fix(ivt(p)),b_fix(ivt(p))&
,c_fix(ivt(p)) ,big_cost,crootfr(p,j),s_fix(ivt(p)),tc_soisno(c,j))

select case (nfix_method)
case ('Houlton')
costNit(j,icostFix) = fun_cost_fix(fixer,&
a_fix(ivt(p)),b_fix(ivt(p)),c_fix(ivt(p)),&
big_cost,crootfr(p,j),s_fix(ivt(p)),tc_soisno(c,j))
case ('Bytnerowicz') ! no acclimation calculation
costNit(j,icostFix) = fun_cost_fix_Bytnerowicz_noAcc(fixer, &
nfix_tmin(ivt(p)),nfix_topt(ivt(p)),nfix_tmax(ivt(p)), &
big_cost,crootfr(p,j),s_fix(ivt(p)),tc_soisno(c,j))
case default
write(iulog,*) subname//' ERROR: unknown nfix_method value: ', nfix_method
call endrun(msg=errMsg(sourcefile, __LINE__))
end select

end do
cost_fix(p,1:nlevdecomp) = costNit(:,icostFix)


!--------------------------------------------------------------------
!------------
! If passive uptake is insufficient, consider fixation,
! mycorrhizal
! non-mycorrhizal, storage, and retranslocation.
!--------------------------------------------------------------------
!------------
!--------------------------------------------------------------------
!------------
! Costs of active uptake.
!--------------------------------------------------------------------
!------------
!------Mycorrhizal Uptake Cost-----------------!
do j = 1,nlevdecomp
rootc_dens_step = rootc_dens(p,j) * permyc(p,istp)
Expand Down Expand Up @@ -1608,6 +1619,44 @@ real(r8) function fun_cost_fix(fixer,a_fix,b_fix,c_fix,big_cost,crootfr,s_fix, t
end if ! ends up with the fixer or non-fixer decision

end function fun_cost_fix


!=========================================================================================
real(r8) function fun_cost_fix_Bytnerowicz_noAcc(fixer,nfix_tmin,nfix_topt,nfix_tmax,big_cost,crootfr,s_fix, tc_soisno)

! Description:
! Calculate the cost of fixing N by nodules.
! Code Description:
! This code is written to CTSM5.1 by Will Wieder 11/17/2022, modified for CLM6 11/01/2024

implicit none
!--------------------------------------------------------------------------
! Function result.
!--------------------------------------------------------------------------
! real(r8) , intent(out) :: cost_of_n !!! cost of fixing N (kgC/kgN)
!--------------------------------------------------------------------------
integer, intent(in) :: fixer ! flag indicating if plant is a fixer
! 1=yes, otherwise no.
real(r8), intent(in) :: nfix_tmin ! As in Bytnerowicz et al. (2022)
real(r8), intent(in) :: nfix_topt ! As in Bytnerowicz et al. (2022)
real(r8), intent(in) :: nfix_tmax ! As in Bytnerowicz et al. (2022)
real(r8), intent(in) :: big_cost ! an arbitrary large cost (gC/gN)
real(r8), intent(in) :: crootfr ! fraction of roots for carbon that are in this layer
real(r8), intent(in) :: s_fix ! Inverts the temperature function for a cost function
real(r8), intent(in) :: tc_soisno ! soil temperature (degrees Celsius)

if (fixer == 1 .and. crootfr > 1.e-6_r8 .and. tc_soisno > nfix_tmin .and. tc_soisno < nfix_tmax) then
fun_cost_fix_Bytnerowicz_noAcc = (-1*s_fix) * 1._r8 / ( ((nfix_tmax-tc_soisno)/(nfix_tmax-nfix_topt))*&
( ((tc_soisno-nfix_tmin)/(nfix_topt-nfix_tmin))**&
((nfix_topt- nfix_tmin)/(nfix_tmax-nfix_topt)) ) )
fun_cost_fix_Bytnerowicz_noAcc = min(fun_cost_fix_Bytnerowicz_noAcc,big_cost)
else
fun_cost_fix_Bytnerowicz_noAcc = big_cost
end if ! ends up with the fixer or non-fixer decision

end function fun_cost_fix_Bytnerowicz_noAcc
!=========================================================================================

!=========================================================================================
real(r8) function fun_cost_active(sminn_layer,big_cost,kc_active,kn_active,rootc_dens,crootfr,smallValue)

Expand Down
2 changes: 2 additions & 0 deletions src/main/clm_varctl.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -227,6 +227,8 @@ module clm_varctl
!
real(r8), public :: nfix_timeconst = -1.2345_r8

character(len=25), public :: nfix_method ! choice of nfix parameterization

!----------------------------------------------------------
! Physics
!----------------------------------------------------------
Expand Down
3 changes: 3 additions & 0 deletions src/main/controlMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -257,6 +257,8 @@ subroutine control_init(dtime)
CNratio_floating, lnc_opt, reduce_dayl_factor, vcmax_opt, &
CN_evergreen_phenology_opt, carbon_resp_opt

namelist /clm_nitrogen/ nfix_method

namelist /clm_inparm/ use_soil_moisture_streams

! excess ice flag
Expand Down Expand Up @@ -882,6 +884,7 @@ subroutine control_spmd()
call mpi_bcast (use_c13_timeseries, 1, MPI_LOGICAL, 0, mpicom, ier)
call mpi_bcast (atm_c13_filename, len(atm_c13_filename), MPI_CHARACTER, 0, mpicom, ier)
call mpi_bcast (use_fun, 1, MPI_LOGICAL, 0, mpicom, ier)
call mpi_bcast (nfix_method, len(nfix_method), MPI_CHARACTER, 0, mpicom, ier)
end if

call mpi_bcast (perchroot, 1, MPI_LOGICAL, 0, mpicom, ier)
Expand Down
18 changes: 18 additions & 0 deletions src/main/pftconMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -268,6 +268,9 @@ module pftconMod
real(r8), allocatable :: b_fix (:) ! A BNF parameter
real(r8), allocatable :: c_fix (:) ! A BNF parameter
real(r8), allocatable :: s_fix (:) ! A BNF parameter
real(r8), allocatable :: nfix_tmin (:) ! A BNF parameter
real(r8), allocatable :: nfix_topt (:) ! A BNF parameter
real(r8), allocatable :: nfix_tmax (:) ! A BNF parameter
real(r8), allocatable :: akc_active (:) ! A mycorrhizal uptake parameter
real(r8), allocatable :: akn_active (:) ! A mycorrhizal uptake parameter
real(r8), allocatable :: ekc_active (:) ! A mycorrhizal uptake parameter
Expand Down Expand Up @@ -485,6 +488,9 @@ subroutine InitAllocate (this)
allocate( this%b_fix (0:mxpft) )
allocate( this%c_fix (0:mxpft) )
allocate( this%s_fix (0:mxpft) )
allocate( this%nfix_tmin (0:mxpft) )
allocate( this%nfix_topt (0:mxpft) )
allocate( this%nfix_tmax (0:mxpft) )
allocate( this%akc_active (0:mxpft) )
allocate( this%akn_active (0:mxpft) )
allocate( this%ekc_active (0:mxpft) )
Expand Down Expand Up @@ -880,6 +886,15 @@ subroutine InitRead(this)
call ncd_io('s_fix', this%s_fix, 'read', ncid, readvar=readv, posNOTonfile=.true.)
if ( .not. readv ) call endrun(msg=' ERROR: error in reading in pft data'//errMsg(sourcefile, __LINE__))

call ncd_io('nfix_tmin', this%nfix_tmin, 'read', ncid, readvar=readv, posNOTonfile=.true.)
if ( .not. readv ) call endrun(msg=' ERROR: error in reading in pft data'//errMsg(sourcefile, __LINE__))

call ncd_io('nfix_topt', this%nfix_topt, 'read', ncid, readvar=readv, posNOTonfile=.true.)
if ( .not. readv ) call endrun(msg=' ERROR: error in reading in pft data'//errMsg(sourcefile, __LINE__))

call ncd_io('nfix_tmax', this%nfix_tmax, 'read', ncid, readvar=readv, posNOTonfile=.true.)
if ( .not. readv ) call endrun(msg=' ERROR: error in reading in pft data'//errMsg(sourcefile, __LINE__))

call ncd_io('akc_active', this%akc_active, 'read', ncid, readvar=readv, posNOTonfile=.true.)
if ( .not. readv ) call endrun(msg=' ERROR: error in reading in pft data'//errMsg(sourcefile, __LINE__))

Expand Down Expand Up @@ -1575,6 +1590,9 @@ subroutine Clean(this)
deallocate( this%b_fix)
deallocate( this%c_fix)
deallocate( this%s_fix)
deallocate( this%nfix_tmin)
deallocate( this%nfix_topt)
deallocate( this%nfix_tmax)
deallocate( this%akc_active)
deallocate( this%akn_active)
deallocate( this%ekc_active)
Expand Down

0 comments on commit 2db5c6a

Please sign in to comment.