Skip to content

Commit

Permalink
Merge pull request #2 from lmbirch89/abz_update
Browse files Browse the repository at this point in the history 
Abz update
  • Loading branch information
@lmbirch89
lmbirch89 authored Mar 21, 2020
2 parents 95fc9e1 + 8a5252e commit ea984f5
Show file tree
Hide file tree
Showing 5 changed files with 230 additions and 44 deletions.
33 changes: 27 additions & 6 deletions src/biogeochem/CNPhenologyMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -299,7 +299,7 @@ subroutine CNPhenology (bounds, num_soilc, filter_soilc, num_soilp, &

call CNSeasonDecidPhenology(num_soilp, filter_soilp, &
temperature_inst, cnveg_state_inst, dgvs_inst, &
cnveg_carbonstate_inst, cnveg_nitrogenstate_inst, cnveg_carbonflux_inst, cnveg_nitrogenflux_inst)
cnveg_carbonstate_inst, cnveg_nitrogenstate_inst, cnveg_carbonflux_inst, cnveg_nitrogenflux_inst, waterstatebulk_inst)

call CNStressDecidPhenology(num_soilp, filter_soilp, &
soilstate_inst, temperature_inst, atm2lnd_inst, wateratm2lndbulk_inst, cnveg_state_inst, &
Expand Down Expand Up @@ -670,7 +670,8 @@ end subroutine CNEvergreenPhenology
!-----------------------------------------------------------------------
subroutine CNSeasonDecidPhenology (num_soilp, filter_soilp , &
temperature_inst, cnveg_state_inst, dgvs_inst , &
cnveg_carbonstate_inst, cnveg_nitrogenstate_inst, cnveg_carbonflux_inst, cnveg_nitrogenflux_inst)
cnveg_carbonstate_inst, cnveg_nitrogenstate_inst, cnveg_carbonflux_inst, cnveg_nitrogenflux_inst, &
waterstatebulk_inst)
!
! !DESCRIPTION:
! For coupled carbon-nitrogen code (CN).
Expand All @@ -692,12 +693,15 @@ subroutine CNSeasonDecidPhenology (num_soilp, filter_soilp , &
type(cnveg_nitrogenstate_type) , intent(inout) :: cnveg_nitrogenstate_inst
type(cnveg_carbonflux_type) , intent(inout) :: cnveg_carbonflux_inst
type(cnveg_nitrogenflux_type) , intent(inout) :: cnveg_nitrogenflux_inst
type(waterstatebulk_type) , intent(in) :: waterstatebulk_inst
!
! !LOCAL VARIABLES:
integer :: g,c,p !indices
integer :: fp !lake filter patch index
real(r8):: ws_flag !winter-summer solstice flag (0 or 1)
real(r8):: crit_onset_gdd !critical onset growing degree-day sum
real(r8):: crit_daylbirch !latitudinal gradient in arctic-boreal
real(r8):: onset_thresh !lbirch: flag onset threshold
real(r8):: soilt
!-----------------------------------------------------------------------

Expand Down Expand Up @@ -780,7 +784,10 @@ subroutine CNSeasonDecidPhenology (num_soilp, filter_soilp , &
livestemn_storage_to_xfer => cnveg_nitrogenflux_inst%livestemn_storage_to_xfer_patch , & ! Output: [real(r8) (:) ]
deadstemn_storage_to_xfer => cnveg_nitrogenflux_inst%deadstemn_storage_to_xfer_patch , & ! Output: [real(r8) (:) ]
livecrootn_storage_to_xfer => cnveg_nitrogenflux_inst%livecrootn_storage_to_xfer_patch , & ! Output: [real(r8) (:) ]
deadcrootn_storage_to_xfer => cnveg_nitrogenflux_inst%deadcrootn_storage_to_xfer_patch & ! Output: [real(r8) (:) ]
deadcrootn_storage_to_xfer => cnveg_nitrogenflux_inst%deadcrootn_storage_to_xfer_patch , & ! Output: [real(r8) (:) ]
t10 => temperature_inst%soila10_patch , & ! Output: [real(r8) (:) ]
tmin5 => temperature_inst%t_a5min_patch , & ! Output: [real(r8) (:) ]
snow_depth => waterstate_inst%snow_10day & ! Output: [real(r8) (:) ]
)

! start patch loop
Expand Down Expand Up @@ -907,13 +914,21 @@ subroutine CNSeasonDecidPhenology (num_soilp, filter_soilp , &
if (onset_gddflag(p) == 1.0_r8 .and. soilt > SHR_CONST_TKFRZ) then
onset_gdd(p) = onset_gdd(p) + (soilt-SHR_CONST_TKFRZ)*fracday
end if
!seperate into Arctic boreal and lower latitudes
if (onset_gdd(p) > crit_onset_gdd .and. abs(grc%latdeg(g))<45.0_r8) then
onset_thresh=1.0_r8
else if (onset_gddflag(p) == 1.0_r8 .and. t10(p) > SHR_CONST_TKFRZ .and. tmin5(p) > SHR_CONST_TKFRZ .and. ws_flag==1.0_r8 .and. dayl(g)>(crit_dayl/2.0_r8) .and. snow_depth(c)<0.1_r8) then
onset_thresh=1.0_r8
end if


! set onset_flag if critical growing degree-day sum is exceeded
if (onset_gdd(p) > crit_onset_gdd) then
if (onset_thresh == 1.0_r8) then
onset_flag(p) = 1.0_r8
dormant_flag(p) = 0.0_r8
onset_gddflag(p) = 0.0_r8
onset_gdd(p) = 0.0_r8
onset_thresh = 0.0_r8
onset_counter(p) = ndays_on * secspday

! move all the storage pools into transfer pools,
Expand Down Expand Up @@ -954,9 +969,15 @@ subroutine CNSeasonDecidPhenology (num_soilp, filter_soilp , &
days_active(p) = days_active(p) + fracday
if (days_active(p) > 355._r8) pftmayexist(p) = .false.
end if

! use 15 hr max to ~11hours in temperate regions
! only begin to test for offset daylength once past the summer sol
if (ws_flag == 0._r8 .and. dayl(g) < crit_dayl) then
crit_daylbirch=54000-360*(65-grc%latdeg(g))
if (crit_daylbirch < crit_dayl) then
crit_daylbirch = crit_dayl
end if

!print*,'lbirch',crit_daylbirch
if (ws_flag == 0._r8 .and. dayl(g) < crit_daylbirch) then
offset_flag(p) = 1._r8
offset_counter(p) = ndays_off * secspday
prev_leafc_to_litter(p) = 0._r8
Expand Down
36 changes: 21 additions & 15 deletions src/biogeophys/LunaMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -306,7 +306,9 @@ subroutine Update_Photosynthesis_Capacity(bounds, fn, filterp, &
vcmx25_z => photosyns_inst%vcmx25_z_patch , & ! Output: [real(r8) (:,:) ] patch leaf Vc,max25 (umol/m2 leaf/s) for canopy layer
jmx25_z => photosyns_inst%jmx25_z_patch , & ! Output: [real(r8) (:,:) ] patch leaf Jmax25 (umol electron/m**2/s) for canopy layer
pnlc_z => photosyns_inst%pnlc_z_patch , & ! Output: [real(r8) (:,:) ] patch proportion of leaf nitrogen allocated for light capture for canopy layer
enzs_z => photosyns_inst%enzs_z_patch & ! Output: [real(r8) (:,:) ] enzyme decay status 1.0-fully active; 0-all decayed during stress
enzs_z => photosyns_inst%enzs_z_patch , & ! Output: [real(r8) (:,:) ] enzyme decay status 1.0-fully active; 0-all decayed during stress
vcmx_prevyr => photosyns_inst%vcmx_prevyr , & ! Output: [real(r8) (:,:) ] patch leaf Vc,max25 from previous year avg
jmx_prevyr => photosyns_inst%jmx_prevyr & ! Output: [real(r8) (:,:) ] patch leaf Jmax25 from previous year avg
)
!----------------------------------------------------------------------------------------------------------------------------------------------------------
!set timestep
Expand All @@ -332,7 +334,7 @@ subroutine Update_Photosynthesis_Capacity(bounds, fn, filterp, &
hourpd = dayl(g) / 3600._r8
tleafd10 = t_veg10_day(p) - tfrz
tleafn10 = t_veg10_night(p) - tfrz
tleaf10 = (dayl(g)*tleafd10 +(86400._r8-dayl(g)) * tleafd10)/86400._r8
tleaf10 = (dayl(g)*tleafd10 +(86400._r8-dayl(g)) * tleafn10)/86400._r8
tair10 = t10(p)- tfrz
relh10 = min(1.0_r8, rh10_p(p))
rb10v = rb10_p(p)
Expand Down Expand Up @@ -409,11 +411,15 @@ subroutine Update_Photosynthesis_Capacity(bounds, fn, filterp, &

chg = vcmx25_opt-vcmx25_z(p, z)
chg_constrn = min(abs(chg),vcmx25_z(p, z)*max_daily_pchg)
vcmx_prevyr(p,z) = vcmx25_z(p,z)
vcmx25_z(p, z) = vcmx25_z(p, z)+sign(1.0_r8,chg)*chg_constrn

vcmx_prevyr(p,z) = (vcmx_prevyr(p,z)+vcmx25_z(p,z))/2.0_r8

chg = jmx25_opt-jmx25_z(p, z)
chg_constrn = min(abs(chg),jmx25_z(p, z)*max_daily_pchg)
jmx_prevyr(p,z) = jmx25_z(p,z)
jmx25_z(p, z) = jmx25_z(p, z)+sign(1.0_r8,chg)*chg_constrn
jmx_prevyr(p,z) = (jmx_prevyr(p,z)+jmx25_z(p,z))/2.0_r8

PNlc_z(p, z)= PNlcopt

Expand Down Expand Up @@ -472,8 +478,8 @@ subroutine Update_Photosynthesis_Capacity(bounds, fn, filterp, &
endif !if not C3 plants
else
do z = 1 , nrad(p)
jmx25_z(p, z) = 85._r8
vcmx25_z(p, z) = 50._r8
jmx25_z(p, z) = jmx_prevyr(p,z)
vcmx25_z(p, z) = vcmx_prevyr(p,z)
end do
endif !checking for LAI and LNC
endif !the first daycheck
Expand Down Expand Up @@ -792,7 +798,7 @@ end subroutine Clear24_Climate_LUNA
subroutine NitrogenAllocation(FNCa,forc_pbot10, relh10, CO2a10,O2a10, PARi10,PARimx10,rb10, hourpd, tair10, tleafd10, tleafn10, &
Jmaxb0, Jmaxb1, Wc2Wjb0, relhExp,&
PNlcold, PNetold, PNrespold, PNcbold, &
PNstoreopt, PNlcopt, PNetopt, PNrespopt, PNcbopt)
PNstoreopt, PNlcopt, PNetopt, PNrespopt, PNcbopt, dayl_factor)
implicit none
real(r8), intent (in) :: FNCa !Area based functional nitrogen content (g N/m2 leaf)
real(r8), intent (in) :: forc_pbot10 !10-day mean air pressure (Pa)
Expand All @@ -819,7 +825,7 @@ subroutine NitrogenAllocation(FNCa,forc_pbot10, relh10, CO2a10,O2a10, PARi10,PAR
real(r8), intent (out):: PNetopt !optimal proportion of nitrogen for electron transport
real(r8), intent (out):: PNrespopt !optimal proportion of nitrogen for respiration
real(r8), intent (out):: PNcbopt !optial proportion of nitrogen for carboxyaltion

real(r8), intent(in) :: dayl_factor !lbirch: added to scaled light:
!-------------------------------------------------------------------------------------------------------------------------------
!intermediate variables
real(r8) :: Carboncost1 !absolute amount of carbon cost associated with maintenance respiration due to deccrease in light capture nitrogen(g dry mass per day)
Expand Down Expand Up @@ -897,11 +903,11 @@ subroutine NitrogenAllocation(FNCa,forc_pbot10, relh10, CO2a10,O2a10, PARi10,PAR
tleafd10c = min(max(tleafd10, Trange1), Trange2) !constrain the physiological range
tleafn10c = min(max(tleafn10, Trange1), Trange2) !constrain the physiological range
ci = 0.7_r8 * CO2a10
JmaxCoef = Jmaxb1 * ((hourpd / 12.0_r8)**2.0_r8) * (1.0_r8 - exp(-relhExp * max(relh10 - minrelh, 0.0_r8) / &
JmaxCoef = Jmaxb1 * dayl_factor * (1.0_r8 - exp(-relhExp * max(relh10 - minrelh, 0.0_r8) / &
(1.0_r8 - minrelh)))
do while (PNlcoldi .NE. PNlc .and. jj < 100)
Fc = VcmxTKattge(tair10, tleafd10c) * Fc25
Fj = JmxTKattge(tair10, tleafd10c) * Fj25
Fc = VcmxTLeuning(tair10, tleafd10c) * Fc25
Fj = JmxTLeuning(tair10, tleafd10c) * Fj25
NUEr = Cv * NUEr25 * (RespTBernacchi(tleafd10c) * hourpd + RespTBernacchi(tleafn10c) * (24.0_r8 - hourpd)) !nitrogen use efficiency for respiration (g biomass/m2/day/g N)
!****************************************************
!Nitrogen Allocation Scheme: store the initial value
Expand Down Expand Up @@ -1054,7 +1060,7 @@ subroutine Nitrogen_investments (KcKjFlag, FNCa, Nlc, forc_pbot10, relh10, &
A = (1.0_r8 - theta_cj) * max(Wc, Wj) + theta_cj * min(Wc, Wj)
endif
PSN = Cv * A * hourpd
Vcmaxnight = VcmxTKattge(tair10, tleafn10) / VcmxTKattge(tair10, tleafd10) * Vcmax
Vcmaxnight = VcmxTLeuning(tair10, tleafn10) / VcmxTLeuning(tair10, tleafd10) * Vcmax
RESP = Cv * leaf_mr_vcm * (Vcmax * hourpd + Vcmaxnight * (24.0_r8 - hourpd))
Net = Jmax / Fj
Ncb = Vcmax / Fc
Expand Down Expand Up @@ -1209,8 +1215,8 @@ subroutine NUEref(NUEjref,NUEcref,Kj2Kcref)

tgrow = 25.0_r8
tleaf = 25.0_r8
Fc = VcmxTKattge(tgrow, tleaf) * Fc25
Fj = JmxTKattge(tgrow, tleaf) * Fj25
Fc = VcmxTLeuning(tgrow, tleaf) * Fc25
Fj = JmxTLeuning(tgrow, tleaf) * Fj25
CO2c = co2ref * forc_pbot_ref * 1.0e-6_r8 !pa
O2c = O2ref * forc_pbot_ref * 1.0e-6_r8 !pa
k_c = params_inst%kc25_coef * exp((79430.0_r8 / (rgas*1.e-3_r8 * (25.0_r8 + tfrz))) * (1.0_r8 - (tfrz + 25.0_r8) / (tfrz + tleaf)))
Expand Down Expand Up @@ -1250,8 +1256,8 @@ subroutine NUE(O2a, ci, tgrow, tleaf, NUEj,NUEc,Kj2Kc)
real(r8) :: awc !second deminator term for rubsico limited carboxylation rate based on Farquhar model
real(r8) :: c_p !CO2 compenstation point (Pa)

Fc = VcmxTKattge(tgrow, tleaf) * Fc25
Fj = JmxTKattge(tgrow, tleaf) * Fj25
Fc = VcmxTLenuning(tgrow, tleaf) * Fc25
Fj = JmxTLeuning(tgrow, tleaf) * Fj25
k_c = params_inst%kc25_coef * exp((79430.0_r8 / (rgas*1.e-3_r8 * (25.0_r8 + tfrz))) * (1.0_r8 - (tfrz + 25.0_r8) / (tfrz + tleaf)))
k_o = params_inst%ko25_coef * exp((36380.0_r8 / (rgas*1.e-3_r8 * (25.0_r8 + tfrz))) * (1.0_r8 - (tfrz + 25.0_r8) / (tfrz + tleaf)))
c_p = params_inst%cp25_yr2000 * exp((37830.0_r8 / (rgas*1.e-3_r8 * (25.0_r8 + tfrz))) * (1.0_r8 - (tfrz + 25.0_r8) / (tfrz + tleaf)))
Expand Down
31 changes: 22 additions & 9 deletions src/biogeophys/PhotosynthesisMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -183,6 +183,8 @@ module PhotosynthesisMod
! LUNA specific variables
real(r8), pointer, public :: vcmx25_z_patch (:,:) ! patch leaf Vc,max25 (umol CO2/m**2/s) for canopy layer
real(r8), pointer, public :: jmx25_z_patch (:,:) ! patch leaf Jmax25 (umol electron/m**2/s) for canopy layer
real(r8), pointer, public :: vcmx_prevyr (:,:) ! patch leaf Vc,max25 previous year avg
real(r8), pointer, public :: jmx2_prevyr (:,:) ! patch leaf Jmax25 previous year avg
real(r8), pointer, public :: pnlc_z_patch (:,:) ! patch proportion of leaf nitrogen allocated for light capture for canopy layer
real(r8), pointer, public :: enzs_z_patch (:,:) ! enzyme decay status 1.0-fully active; 0-all decayed during stress
real(r8), pointer, public :: fpsn24_patch (:) ! 24 hour mean patch photosynthesis (umol CO2/m**2 ground/day)
Expand Down Expand Up @@ -328,6 +330,8 @@ subroutine InitAllocate(this, bounds)
! statements.
allocate(this%vcmx25_z_patch (begp:endp,1:nlevcan)) ; this%vcmx25_z_patch (:,:) = 30._r8
allocate(this%jmx25_z_patch (begp:endp,1:nlevcan)) ; this%jmx25_z_patch (:,:) = 60._r8
allocate(this%vcmx_prevyr (begp:endp,1:nlevcan)) ; this%vcmx_prevyr (:,:) = 30._r8
allocate(this%jmx_prevyr (begp:endp,1:nlevcan)) ; this%jmx_prevyr (:,:) = 60._r8
allocate(this%pnlc_z_patch (begp:endp,1:nlevcan)) ; this%pnlc_z_patch (:,:) = 0.01_r8
allocate(this%fpsn24_patch (begp:endp)) ; this%fpsn24_patch (:) = nan
allocate(this%enzs_z_patch (begp:endp,1:nlevcan)) ; this%enzs_z_patch (:,:) = 1._r8
Expand Down Expand Up @@ -833,6 +837,14 @@ subroutine Restart(this, bounds, ncid, flag)
dim1name='pft', dim2name='levcan', switchdim=.true., &
long_name='Maximum carboxylation rate at 25 celcius for canopy layers', units='umol CO2/m**2/s', &
interpinic_flag='interp', readvar=readvar, data=this%jmx25_z_patch)
call restartvar(ncid=ncid, flag=flag, varname='vcmx_prevyr', xtype=ncd_double, &
dim1name='pft', dim2name='levcan', switchdim=.true., &
long_name='avg carboxylation rate at 25 celcius for canopy layers', units='umol CO2/m**2/s', &
interpinic_flag='interp', readvar=readvar, data=this%vcmx_prevyr)
call restartvar(ncid=ncid, flag=flag, varname='jmx_prevyr', xtype=ncd_double, &
dim1name='pft', dim2name='levcan', switchdim=.true., &
long_name='avg carboxylation rate at 25 celcius for canopy layers', units='umol CO2/m**2/s', &
interpinic_flag='interp', readvar=readvar, data=this%jmx_prevyr)
call restartvar(ncid=ncid, flag=flag, varname='pnlc_z', xtype=ncd_double, &
dim1name='pft', dim2name='levcan', switchdim=.true., &
long_name='proportion of leaf nitrogen allocated for light capture', units='unitless', &
Expand Down Expand Up @@ -2791,18 +2803,18 @@ subroutine PhotosynthesisHydraulicStress ( bounds, fn, filterp, &
kcha = 79430._r8
koha = 36380._r8
cpha = 37830._r8
vcmaxha = 72000._r8
jmaxha = 50000._r8
tpuha = 72000._r8
vcmaxha = 73637._r8
jmaxha = 50300._r8
tpuha = 73637._r8
lmrha = 46390._r8

! High temperature deactivation, from:
! Leuning (2002) Plant, Cell and Environment 25:1205-1210
! The factor "c" scales the deactivation to a value of 1.0 at 25C

vcmaxhd = 200000._r8
jmaxhd = 200000._r8
tpuhd = 200000._r8
vcmaxhd = 149252._r8
jmaxhd = 152044._r8
tpuhd = 149252._r8
lmrhd = 150650._r8
lmrse = 490._r8
lmrc = fth25 (lmrhd, lmrse)
Expand Down Expand Up @@ -3162,9 +3174,10 @@ subroutine PhotosynthesisHydraulicStress ( bounds, fn, filterp, &
kp25_sha = kp25top * nscaler_sha

! Adjust for temperature

vcmaxse = 668.39_r8 - 1.07_r8 * min(max((t10(p)-tfrz),11._r8),35._r8)
jmaxse = 659.70_r8 - 0.75_r8 * min(max((t10(p)-tfrz),11._r8),35._r8)
vcmaxse = 486.0_r8
jmaxse = 495.0_r8
!vcmaxse = 668.39_r8 - 1.07_r8 * min(max((t10(p)-tfrz),11._r8),35._r8)
!jmaxse = 659.70_r8 - 0.75_r8 * min(max((t10(p)-tfrz),11._r8),35._r8)
tpuse = vcmaxse
vcmaxc = fth25 (vcmaxhd, vcmaxse)
jmaxc = fth25 (jmaxhd, jmaxse)
Expand Down
Loading

0 comments on commit ea984f5

Please sign in to comment.