Skip to content

Commit

Permalink
Merge pull request #2 from rgknox/xuchongang/HYDRO_UPDATE_v2
Browse files Browse the repository at this point in the history
minor fixes to the hydro updates
  • Loading branch information
rgknox authored Apr 29, 2019
2 parents 4bf50a4 + b4d643e commit ce61db1
Show file tree
Hide file tree
Showing 3 changed files with 56 additions and 65 deletions.
6 changes: 5 additions & 1 deletion biogeophys/FatesPlantHydraulicsMod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -41,7 +41,7 @@ module FatesPlantHydraulicsMod
use FatesConstantsMod, only : g_per_kg

use EDParamsMod , only : hydr_kmax_rsurf1
use EDParamsMod , only : hydr_kmax_rsurf2
! use EDParamsMod , only : hydr_kmax_rsurf2

use EDTypesMod , only : ed_site_type
use EDTypesMod , only : ed_patch_type
Expand Down Expand Up @@ -886,6 +886,10 @@ subroutine UpdateWaterDepTreeHydrCond(currentSite,ccohort,nlevsoi_hyd,bc_in)
! which is equiv to [kg m-1 s-1 MPa-1]
real(r8) :: kmax_root_surf ! maximum conducitivity for unit root surface (kg water/m2 root area/Mpa/s)

! (RGK 4-2019) THE FOLLOWING SHOULD BE ADDED TO THE PARAMETER FILE IN NEXT GROUPING
real(r8), parameter :: hydr_kmax_rsurf2 = 0.0001_r8 ! kg water/m2 root area/Mpa/s


ccohort_hydr => ccohort%co_hydr
csite_hydr => currentSite%si_hydr
k = 1 !only for the first soil shell
Expand Down
95 changes: 41 additions & 54 deletions biogeophys/FatesPlantRespPhotosynthMod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -45,6 +45,7 @@ module FATESPlantRespPhotosynthMod
use PRTGenericMod, only : store_organ
use PRTGenericMod, only : repro_organ
use PRTGenericMod, only : struct_organ
use EDParamsMod, only : ED_val_bbopt_c3, ED_val_bbopt_c4, ED_val_base_mr_20

! CIME Globals
use shr_log_mod , only : errMsg => shr_log_errMsg
Expand Down Expand Up @@ -96,7 +97,7 @@ subroutine FatesPlantRespPhotosynthDrive (nsites, sites,bc_in,bc_out,dtime)
use FatesConstantsMod, only : rgas => rgas_J_K_kmol
use FatesConstantsMod, only : tfrz => t_water_freeze_k_1atm
use FatesParameterDerivedMod, only : param_derived
use EDParamsMod, only : ED_val_bbopt_c3, ED_val_bbopt_c4, ED_val_base_mr_20

use FatesAllometryMod, only : bleaf
use FatesAllometryMod, only : storage_fraction_of_target
use FatesAllometryMod, only : set_root_fraction
Expand Down Expand Up @@ -394,10 +395,10 @@ subroutine FatesPlantRespPhotosynthDrive (nsites, sites,bc_in,bc_out,dtime)
(hlm_use_planthydro.eq.itrue) .or. &
(nleafage > 1) .or. &
(hlm_parteh_mode .ne. prt_carbon_allom_hyp ) ) then

if (hlm_use_planthydro.eq.itrue ) then

bbb = max (bbbopt(nint(c3psn(ft)))*currentCohort%co_hydr%btran(1), 1._r8)
if (hlm_use_planthydro.eq.itrue ) then
bbb = max( cf/rsmax0, bbbopt(nint(c3psn(ft)))*currentCohort%co_hydr%btran(1) )
btran_eff = currentCohort%co_hydr%btran(1)

! dinc_ed is the total vegetation area index of each "leaf" layer
Expand All @@ -414,7 +415,7 @@ subroutine FatesPlantRespPhotosynthDrive (nsites, sites,bc_in,bc_out,dtime)

else

bbb = max (bbbopt(nint(c3psn(ft)))*currentPatch%btran_ft(ft), 1._r8)
bbb = max( cf/rsmax0, bbbopt(nint(c3psn(ft)))*currentPatch%btran_ft(ft) )
btran_eff = currentPatch%btran_ft(ft)
! For consistency sake, we use total LAI here, and not exposed
! if the plant is under-snow, it will be effectively dormant for
Expand Down Expand Up @@ -844,7 +845,7 @@ subroutine LeafLayerPhotosynthesis(f_sun_lsl, & ! in
! ------------------------------------------------------------------------------------

use EDPftvarcon , only : EDPftvarcon_inst
use EDParamsMod, only : ED_val_bbopt_c3, ED_val_bbopt_c4


! Arguments
! ------------------------------------------------------------------------------------
Expand Down Expand Up @@ -919,11 +920,19 @@ subroutine LeafLayerPhotosynthesis(f_sun_lsl, & ! in
real(r8) :: leaf_co2_ppress ! CO2 partial pressure at leaf surface (Pa)
real(r8) :: init_co2_inter_c ! First guess intercellular co2 specific to C path

real(r8), dimension(0:1) :: bbbopt ! Cuticular conductance at full water potential (umol H2O /m2/s)



! Parameters
! ------------------------------------------------------------------------
! Fraction of light absorbed by non-photosynthetic pigments
real(r8),parameter :: fnps = 0.15_r8

! For plants with no leaves, a miniscule amount of conductance
! can happen through the stems, at a partial rate of cuticular conductance
real(r8),parameter :: stem_cuticle_loss_frac = 0.1_r8

! empirical curvature parameter for electron transport rate
real(r8),parameter :: theta_psii = 0.7_r8

Expand All @@ -943,18 +952,15 @@ subroutine LeafLayerPhotosynthesis(f_sun_lsl, & ! in
! empirical curvature parameter for ap photosynthesis co-limitation
real(r8),parameter :: theta_ip = 0.999_r8

real(r8), dimension(0:1) :: bbbopt !cuticular conductance

associate( bb_slope => EDPftvarcon_inst%BB_slope) ! slope of BB relationship



bbbopt(0) = ED_val_bbopt_c4
bbbopt(1) = ED_val_bbopt_c3

! photosynthetic pathway: 0. = c4, 1. = c3
c3c4_path_index = nint(EDPftvarcon_inst%c3psn(ft))

bbbopt(0) = ED_val_bbopt_c4
bbbopt(1) = ED_val_bbopt_c3

if (c3c4_path_index == 1) then
init_co2_inter_c = init_a2l_co2_c3 * can_co2_ppress
else
Expand All @@ -968,25 +974,19 @@ subroutine LeafLayerPhotosynthesis(f_sun_lsl, & ! in

anet_av_out = -lmr
psn_out = 0._r8
if(btran>0._r8) then
rstoma_out = min(rsmax0, cf*1._r8/(bbbopt(c3c4_path_index)*btran))
else
rstoma_out = rsmax0
endif

! The cuticular conductance already factored in maximum resistance as a bound
! no need to re-bound it

rstoma_out = cf/bbb

c13disc_z = 0.0_r8 !carbon 13 discrimination in night time carbon flux, note value of 1.0 is used in CLM

else ! day time (a little bit more complicated ...)

! if ( debug ) write(fates_log(),*) 'EDphot 594 ',laisun_lsl
! if ( debug ) write(fates_log(),*) 'EDphot 595 ',laisha_lsl

!is there leaf area? - (NV can be larger than 0 with only stem area if deciduous)
if ( laisun_lsl + laisha_lsl > 0._r8 ) then
if(bbb > 1.0_r8)then !only if stomata open larger than cuticular conductance
!is there leaf area? - (NV can be larger than 0 with only stem area if deciduous)
if ( laisun_lsl + laisha_lsl > 0._r8 ) then

! if ( debug ) write(fates_log(),*) '600 in laisun, laisha loop '

!Loop aroun shaded and unshaded leaves
psn_out = 0._r8 ! psn is accumulated across sun and shaded leaves.
rstoma_out = 0._r8 ! 1/rs is accumulated across sun and shaded leaves.
Expand Down Expand Up @@ -1136,19 +1136,16 @@ subroutine LeafLayerPhotosynthesis(f_sun_lsl, & ! in
! Convert gs_mol (umol /m**2/s) to gs (m/s) and then to rs (s/m)
gs = gs_mol / cf

! estimate carbon 13 discrimination in leaf level carbon flux Liang WEI and Hang ZHOU 2018, based on
! estimate carbon 13 discrimination in leaf level carbon
! flux Liang WEI and Hang ZHOU 2018, based on
! Ubierna and Farquhar, 2014 doi:10.1111/pce.12346, using the simplified model:
! $\Delta ^{13} C = \alpha_s + (b - \alpha_s) \cdot \frac{C_i}{C_a}$
! just hard code b and \alpha_s for now, might move to parameter set in future
! b = 27.0 alpha_s = 4.4
! TODO, not considering C4 or CAM right now, may need to address this
! note co2_inter_c is intracelluar CO2, not intercelluar
c13disc_z = 4.4_r8 + (27.0_r8 - 4.4_r8) * min (can_co2_ppress, max (co2_inter_c, 0._r8)) / can_co2_ppress


! if ( debug ) write(fates_log(),*) 'EDPhoto 737 ', psn_out
! if ( debug ) write(fates_log(),*) 'EDPhoto 738 ', agross
! if ( debug ) write(fates_log(),*) 'EDPhoto 739 ', f_sun_lsl
! note co2_inter_c is intracelluar CO2, not intercelluar
c13disc_z = 4.4_r8 + (27.0_r8 - 4.4_r8) * &
min (can_co2_ppress, max (co2_inter_c, 0._r8)) / can_co2_ppress

! Accumulate total photosynthesis umol/m2 ground/s-1.
! weight per unit sun and sha leaves.
Expand All @@ -1162,10 +1159,6 @@ subroutine LeafLayerPhotosynthesis(f_sun_lsl, & ! in
gstoma = gstoma + &
1._r8/(min(1._r8/gs, rsmax0)) * (1.0_r8-f_sun_lsl)
end if

! if ( debug ) write(fates_log(),*) 'EDPhoto 758 ', psn_out
! if ( debug ) write(fates_log(),*) 'EDPhoto 759 ', agross
! if ( debug ) write(fates_log(),*) 'EDPhoto 760 ', f_sun_lsl

! Make sure iterative solution is correct
if (gs_mol < 0._r8) then
Expand All @@ -1188,25 +1181,19 @@ subroutine LeafLayerPhotosynthesis(f_sun_lsl, & ! in
! This is the stomatal resistance of the leaf layer
rstoma_out = 1._r8/gstoma

else !!set the situations with only cuticular conductance and no photosynthesis

psn_out = 0._r8
anet_av_out = -lmr
if(btran>0._r8) then
rstoma_out = min(rsmax0, cf*1._r8/(bbbopt(c3c4_path_index)*btran))
else
rstoma_out = rsmax0
endif
endif

else
!No leaf area. This layer is present only because of stems.

! No leaf area. This layer is present only because of stems.
! Net assimilation is zero, not negative because there are
! no leaves to even respire
! (leaves are off, or have reduced to 0)
psn_out = 0._r8
rstoma_out = min(rsmax0, cf*1._r8/(0.1_r8*bbbopt(c3c4_path_index))) !assume stem loss is only 10% of leaf culticular
c13disc_z = 0.0_r8

psn_out = 0._r8
anet_av_out = 0._r8
rstoma_out = min(rsmax0, cf/(stem_cuticle_loss_frac*bbbopt(c3c4_path_index)))
c13disc_z = 0.0_r8

end if !is there leaf area?
end if !is there leaf area?


end if ! night or day
Expand Down
20 changes: 10 additions & 10 deletions main/EDParamsMod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -80,10 +80,10 @@ module EDParamsMod
! Hydraulics Control Parameters (ONLY RELEVANT WHEN USE_FATES_HYDR = TRUE)
! ----------------------------------------------------------------------------------------------
real(r8),protected :: hydr_kmax_rsurf1 ! maximum conducitivity for unit root surface for water uptake(kg water/m2 root area/Mpa/s)
real(r8),protected :: hydr_kmax_rsurf2 ! maximum conducitivity for unit root surface for water loss (kg water/m2 root area/Mpa/s)
! real(r8),protected :: hydr_kmax_rsurf2 ! maximum conducitivity for unit root surface for water loss (kg water/m2 root area/Mpa/s)

character(len=param_string_length),parameter :: hydr_name_kmax_rsurf1 = "fates_hydr_kmax_rsurf1"
character(len=param_string_length),parameter :: hydr_name_kmax_rsurf2 = "fates_hydr_kmax_rsurf2"
character(len=param_string_length),parameter :: hydr_name_kmax_rsurf1 = "fates_hydr_kmax_rsurf"
! character(len=param_string_length),parameter :: hydr_name_kmax_rsurf2 = "fates_hydr_kmax_rsurf2"

real(r8),protected :: hydr_psi0 ! sapwood water potential at saturation (MPa)
character(len=param_string_length),parameter :: hydr_name_psi0 = "fates_hydr_psi0"
Expand Down Expand Up @@ -161,7 +161,7 @@ subroutine FatesParamsInit()
ED_val_canopy_closure_thresh = nan

hydr_kmax_rsurf1 = nan
hydr_kmax_rsurf2 = nan
! hydr_kmax_rsurf2 = nan
hydr_psi0 = nan
hydr_psicap = nan

Expand Down Expand Up @@ -266,8 +266,8 @@ subroutine FatesRegisterParams(fates_params)
call fates_params%RegisterParameter(name=hydr_name_kmax_rsurf1, dimension_shape=dimension_shape_1d, &
dimension_names=dim_names)

call fates_params%RegisterParameter(name=hydr_name_kmax_rsurf2, dimension_shape=dimension_shape_1d, &
dimension_names=dim_names)
! call fates_params%RegisterParameter(name=hydr_name_kmax_rsurf2, dimension_shape=dimension_shape_1d, &
! dimension_names=dim_names)

call fates_params%RegisterParameter(name=hydr_name_psi0, dimension_shape=dimension_shape_1d, &
dimension_names=dim_names)
Expand Down Expand Up @@ -390,8 +390,8 @@ subroutine FatesReceiveParams(fates_params)
call fates_params%RetreiveParameter(name=hydr_name_kmax_rsurf1, &
data=hydr_kmax_rsurf1)

call fates_params%RetreiveParameter(name=hydr_name_kmax_rsurf2, &
data=hydr_kmax_rsurf2)
! call fates_params%RetreiveParameter(name=hydr_name_kmax_rsurf2, &
! data=hydr_kmax_rsurf2)

call fates_params%RetreiveParameter(name=hydr_name_psi0, &
data=hydr_psi0)
Expand Down Expand Up @@ -470,8 +470,8 @@ subroutine FatesReportParams(is_master)
write(fates_log(),fmt0) 'ED_val_cohort_fusion_tol = ',ED_val_cohort_fusion_tol
write(fates_log(),fmt0) 'ED_val_patch_fusion_tol = ',ED_val_patch_fusion_tol
write(fates_log(),fmt0) 'ED_val_canopy_closure_thresh = ',ED_val_canopy_closure_thresh
write(fates_log(),fmt0) 'hydr_kmax_rsurf1 = ',hydr_kmax_rsurf1
write(fates_log(),fmt0) 'hydr_kmax_rsurf2 = ',hydr_kmax_rsurf2
write(fates_log(),fmt0) 'hydr_kmax_rsurf1 = ',hydr_kmax_rsurf1
! write(fates_log(),fmt0) 'hydr_kmax_rsurf2 = ',hydr_kmax_rsurf2
write(fates_log(),fmt0) 'hydr_psi0 = ',hydr_psi0
write(fates_log(),fmt0) 'hydr_psicap = ',hydr_psicap
write(fates_log(),fmt0) 'bgc_soil_salinity = ', bgc_soil_salinity
Expand Down

0 comments on commit ce61db1

Please sign in to comment.