Skip to content

Commit

Permalink
Merge pull request ESCOMP#434 from xuchongang/xuchongang/fix_hydrauli…
Browse files Browse the repository at this point in the history 
…c_bug_ready_to_merge

Bug fixes and parameter updates for the HYDRO codes
  • Loading branch information
@rgknox
rgknox authored Oct 30, 2018
2 parents 614a9e3 + 82594ee commit af76d76
Show file tree
Hide file tree
Showing 9 changed files with 40 additions and 37 deletions.
12 changes: 6 additions & 6 deletions biogeochem/EDCanopyStructureMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -504,10 +504,10 @@ subroutine DemoteFromLayer(currentSite,currentPatch,i_lyr)

allocate(copyc)
call InitPRTCohort(copyc)
call copy_cohort(currentCohort, copyc)
if( hlm_use_planthydro.eq.itrue ) then
call InitHydrCohort(currentSite,copyc)
endif
call copy_cohort(currentCohort, copyc)

newarea = currentCohort%c_area - cc_loss
copyc%n = currentCohort%n*newarea/currentCohort%c_area
Expand Down Expand Up @@ -851,11 +851,11 @@ subroutine PromoteIntoLayer(currentSite,currentPatch,i_lyr)

allocate(copyc)
call InitPRTCohort(copyc)
call copy_cohort(currentCohort, copyc) !makes an identical copy...
if( hlm_use_planthydro.eq.itrue ) then
call InitHydrCohort(CurrentSite,copyc)
endif

call copy_cohort(currentCohort, copyc) !makes an identical copy...

newarea = currentCohort%c_area - cc_gain !new area of existing cohort

call carea_allom(currentCohort%dbh,currentCohort%n,currentSite%spread, &
Expand Down Expand Up @@ -1422,11 +1422,11 @@ subroutine leaf_area_profile( currentSite , snow_depth_si, frac_sno_eff_si)
! is obscured by snow.

layer_top_hite = currentCohort%hite - &
( dble(iv-1.0)/currentCohort%NV * currentCohort%hite * &
( real(iv-1,r8)/currentCohort%NV * currentCohort%hite * &
EDPftvarcon_inst%crown(currentCohort%pft) )

layer_bottom_hite = currentCohort%hite - &
( dble(iv)/currentCohort%NV * currentCohort%hite * &
( real(iv,r8)/currentCohort%NV * currentCohort%hite * &
EDPftvarcon_inst%crown(currentCohort%pft) )

fraction_exposed = 1.0_r8
Expand All @@ -1449,7 +1449,7 @@ subroutine leaf_area_profile( currentSite , snow_depth_si, frac_sno_eff_si)

if(iv==currentCohort%NV) then
remainder = (currentCohort%treelai + currentCohort%treesai) - &
(dinc_ed*dble(currentCohort%nv-1.0_r8))
(dinc_ed*real(currentCohort%nv-1,r8))
if(remainder > dinc_ed )then
write(fates_log(), *)'ED: issue with remainder', &
currentCohort%treelai,currentCohort%treesai,dinc_ed, &
Expand Down
2 changes: 1 addition & 1 deletion biogeochem/EDPhysiologyMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -589,7 +589,7 @@ subroutine phenology( currentSite, bc_in )
if ((t >= currentSite%dleafondate - 30.and.t <= currentSite%dleafondate + 30).or.(t > 360 - 15.and. &
currentSite%dleafondate < 15))then ! are we in the window?
! TODO: CHANGE THIS MATH, MOVE THE DENOMENATOR OUTSIDE OF THE SUM (rgk 01-2017)
if (sum(currentSite%water_memory(1:numWaterMem)/dble(numWaterMem)) &
if (sum(currentSite%water_memory(1:numWaterMem)/real(numWaterMem,r8)) &
>= ED_val_phen_drought_threshold.and.currentSite%dstatus == 1.and.t >= 10)then
! leave some minimum time between leaf off and leaf on to prevent 'flickering'.
if (timesincedleafoff > ED_val_phen_doff_time)then
Expand Down
2 changes: 1 addition & 1 deletion biogeophys/EDBtranMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -197,7 +197,7 @@ subroutine btran_ed( nsites, sites, bc_in, bc_out)
cpatch%rootr_ft(ft,j) * pftgs(ft)/sum_pftgs
else
bc_out(s)%rootr_pasl(ifp,j) = bc_out(s)%rootr_pasl(ifp,j) + &
cpatch%rootr_ft(ft,j) * 1._r8/dble(numpft)
cpatch%rootr_ft(ft,j) * 1._r8/real(numpft,r8)
end if
enddo
enddo
Expand Down
25 changes: 13 additions & 12 deletions biogeophys/FatesPlantHydraulicsMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -402,15 +402,16 @@ subroutine updateSizeDepTreeHydProps(currentSite,cc_p,bc_in)

! a_sapwood = a_leaf_tot / EDPftvarcon_inst%allom_latosa_int(FT)*1.e-4_r8
! m2 sapwood = m2 leaf * cm2 sapwood/m2 leaf *1.0e-4m2


! or ...
! a_sapwood = a_leaf_tot / ( 0.001_r8 + 0.025_r8 * cCohort%hite ) * 1.e-4_r8

v_sapwood = a_sapwood * z_stem
ccohort_hydr%v_ag(n_hypool_leaf+1:n_hypool_ag) = v_sapwood / n_hypool_stem

! TRANSPORTING ROOT DEPTH & VOLUME
!in special case where n_hypool_troot = 1, the node depth of the single troot pool
!is the depth at which 50% total root distribution is attained
dcumul_rf = 1._r8/dble(n_hypool_troot)
dcumul_rf = 1._r8/real(n_hypool_troot,r8)

do k=1,n_hypool_troot
cumul_rf = dcumul_rf*k
Expand Down Expand Up @@ -1026,11 +1027,11 @@ subroutine updateSizeDepRhizHydProps(currentSite, bc_in )
kmax_soil_total = 2._r8*pi_const*csite_hydr%l_aroot_layer(j) / &
log(csite_hydr%r_node_shell(j,k)/csite_hydr%rs1(j))*hksat_s
csite_hydr%kmax_upper_shell(j,k) = (1._r8/kmax_root_surf_total + &
1._r8/kmax_soil_total)**-1._r8
1._r8/kmax_soil_total)**(-1._r8)
csite_hydr%kmax_bound_shell(j,k) = (1._r8/kmax_root_surf_total + &
1._r8/kmax_soil_total)**-1._r8
1._r8/kmax_soil_total)**(-1._r8)
csite_hydr%kmax_lower_shell(j,k) = (1._r8/kmax_root_surf_total + &
1._r8/kmax_soil_total)**-1._r8
1._r8/kmax_soil_total)**(-1._r8)
end if
if(j == 1) then
if(csite_hydr%r_node_shell(j,k) <= csite_hydr%rs1(j)) then
Expand All @@ -1041,11 +1042,11 @@ subroutine updateSizeDepRhizHydProps(currentSite, bc_in )
kmax_soil_total = 2._r8*pi_const*csite_hydr%l_aroot_1D / &
log(csite_hydr%r_node_shell_1D(k)/csite_hydr%rs1(j))*hksat_s
csite_hydr%kmax_upper_shell_1D(k) = (1._r8/kmax_root_surf_total + &
1._r8/kmax_soil_total)**-1._r8
1._r8/kmax_soil_total)**(-1._r8)
csite_hydr%kmax_bound_shell_1D(k) = (1._r8/kmax_root_surf_total + &
1._r8/kmax_soil_total)**-1._r8
1._r8/kmax_soil_total)**(-1._r8)
csite_hydr%kmax_lower_shell_1D(k) = (1._r8/kmax_root_surf_total + &
1._r8/kmax_soil_total)**-1._r8
1._r8/kmax_soil_total)**(-1._r8)
end if
end if
else
Expand Down Expand Up @@ -2554,7 +2555,7 @@ subroutine Hydraulics_1DSolve(cc_p, ft, z_node, v_node, ths_node, thr_node, kmax
iterh1 = 0
do while( iterh1 == 0 .or. ((abs(we_local) > thresh .or. supsub_flag /= 0) .and. iterh1 < maxiter) )
dt_fac = max(imult*iterh1,1)
dt_fac2 = DBLE(dt_fac)
dt_fac2 = real(dt_fac,r8)
dt_new = dtime/dt_fac2

!! restore initial states for a fresh attempt using new sub-timesteps
Expand Down Expand Up @@ -4100,7 +4101,7 @@ subroutine swcCampbell_satfrac_from_psi(psi, psisat, B, satfrac)
! !LOCAL VARIABLES:
!------------------------------------------------------------------------------

satfrac = (psi/psisat)**(-1/B)
satfrac = (psi/psisat)**(-1.0_r8/B)

end subroutine swcCampbell_satfrac_from_psi

Expand Down Expand Up @@ -4438,7 +4439,7 @@ subroutine shellGeom(l_aroot, rs1, area, dz, r_out_shell, r_node_shell, v_shell)
r_out_shell(nshell) = (pi_const*l_aroot/(area*dz))**(-0.5_r8) ! eqn(8) S98
if(nshell > 1) then
do k = 1,nshell-1
r_out_shell(k) = rs1*(r_out_shell(nshell)/rs1)**((k+0._r8)/nshell) ! eqn(7) S98
r_out_shell(k) = rs1*(r_out_shell(nshell)/rs1)**((real(k,r8))/real(nshell,r8)) ! eqn(7) S98
enddo
end if

Expand Down
8 changes: 5 additions & 3 deletions fire/SFMainMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -245,7 +245,7 @@ subroutine charecteristics_of_fuel ( currentSite )
endif
! FIX(RF,032414): needs refactoring.
! average water content !is this the correct metric?
timeav_swc = sum(currentSite%water_memory(1:numWaterMem)) / dble(numWaterMem)
timeav_swc = sum(currentSite%water_memory(1:numWaterMem)) / real(numWaterMem,r8)
! Equation B2 in Thonicke et al. 2010
! live grass moisture content depends on upper soil layer
fuel_moisture(lg_sf) = max(0.0_r8, 10.0_r8/9._r8 * timeav_swc - 1.0_r8/9.0_r8)
Expand Down Expand Up @@ -769,7 +769,8 @@ subroutine area_burnt ( currentSite )
! THIS SHOULD HAVE THE COLUMN AND LU AREA WEIGHT ALSO, NO?

gridarea = km2_to_m2 ! 1M m2 in a km2
!NF = number of lighting strikes per day per km2

! NF = number of lighting strikes per day per km2
currentPatch%NF = ED_val_nignitions * currentPatch%area/area /365

! If there are 15 lightening strickes per year, per km2. (approx from NASA product)
Expand All @@ -784,9 +785,10 @@ subroutine area_burnt ( currentSite )
size_of_fire = ((3.1416_r8/(4.0_r8*lb))*((df+db)**2.0_r8))

!AB = daily area burnt = size fires in m2 * num ignitions * prob ignition starts fire
! m2 per km2 per day
currentPatch%AB = size_of_fire * currentPatch%NF * currentSite%FDI

patch_area_in_m2 = gridarea*currentPatch%area/area
patch_area_in_m2 = gridarea *currentPatch%area/area

currentPatch%frac_burnt = currentPatch%AB / patch_area_in_m2
if(write_SF == itrue)then
Expand Down
12 changes: 6 additions & 6 deletions main/EDMainMod.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -330,9 +330,9 @@ subroutine ed_integrate_state_variables(currentSite, bc_in )
currentCohort%npp_acc = currentCohort%npp_acc_hold / hlm_days_per_year
endif
else
currentCohort%npp_acc_hold = currentCohort%npp_acc * dble(hlm_days_per_year)
currentCohort%gpp_acc_hold = currentCohort%gpp_acc * dble(hlm_days_per_year)
currentCohort%resp_acc_hold = currentCohort%resp_acc * dble(hlm_days_per_year)
currentCohort%npp_acc_hold = currentCohort%npp_acc * real(hlm_days_per_year,r8)
currentCohort%gpp_acc_hold = currentCohort%gpp_acc * real(hlm_days_per_year,r8)
currentCohort%resp_acc_hold = currentCohort%resp_acc * real(hlm_days_per_year,r8)
endif

currentSite%flux_in = currentSite%flux_in + currentCohort%npp_acc * currentCohort%n
Expand Down Expand Up @@ -657,9 +657,9 @@ subroutine bypass_dynamics(currentSite)

currentCohort%isnew=.false.

currentCohort%npp_acc_hold = currentCohort%npp_acc * dble(hlm_days_per_year)
currentCohort%gpp_acc_hold = currentCohort%gpp_acc * dble(hlm_days_per_year)
currentCohort%resp_acc_hold = currentCohort%resp_acc * dble(hlm_days_per_year)
currentCohort%npp_acc_hold = currentCohort%npp_acc * real(hlm_days_per_year,r8)
currentCohort%gpp_acc_hold = currentCohort%gpp_acc * real(hlm_days_per_year,r8)
currentCohort%resp_acc_hold = currentCohort%resp_acc * real(hlm_days_per_year,r8)

currentCohort%npp_acc = 0.0_r8
currentCohort%gpp_acc = 0.0_r8
Expand Down
2 changes: 1 addition & 1 deletion main/EDPftvarcon.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -203,7 +203,7 @@ module EDPftvarcon

! PFT Dimension
real(r8), allocatable :: hydr_p_taper(:) ! xylem taper exponent
real(r8), allocatable :: hydr_rs2(:) ! absorbing root radius (mm)
real(r8), allocatable :: hydr_rs2(:) ! absorbing root radius (m)
real(r8), allocatable :: hydr_srl(:) ! specific root length (m g-1)
real(r8), allocatable :: hydr_rfrac_stem(:) ! fraction of total tree resistance from troot to canopy
real(r8), allocatable :: hydr_avuln_gs(:) ! shape parameter for stomatal control of water vapor exiting leaf
Expand Down
2 changes: 1 addition & 1 deletion parameter_files/fates_params_14pfts.cdl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -315,7 +315,7 @@ variables:
fates_hydr_rfrac_stem:units = "fraction" ;
fates_hydr_rfrac_stem:long_name = "fraction of total tree resistance from troot to canopy" ;
float fates_hydr_rs2(fates_pft) ;
fates_hydr_rs2:units = "mm" ;
fates_hydr_rs2:units = "m" ;
fates_hydr_rs2:long_name = "absorbing root radius" ;
float fates_hydr_srl(fates_pft) ;
fates_hydr_srl:units = "m g-1" ;
Expand Down
12 changes: 6 additions & 6 deletions parameter_files/fates_params_default.cdl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -258,7 +258,7 @@ variables:
fates_hydr_rfrac_stem:units = "fraction" ;
fates_hydr_rfrac_stem:long_name = "fraction of total tree resistance from troot to canopy" ;
float fates_hydr_rs2(fates_pft) ;
fates_hydr_rs2:units = "mm" ;
fates_hydr_rs2:units = "m" ;
fates_hydr_rs2:long_name = "absorbing root radius" ;
float fates_hydr_srl(fates_pft) ;
fates_hydr_srl:units = "m g-1" ;
Expand Down Expand Up @@ -657,10 +657,10 @@ data:
-2.25, -2.25 ;

fates_hydr_pinot_node =
-999, -999,
-999, -999,
-999, -999,
-999, -999 ;
-1.465984, -1.465984,
-1.228070, -1.228070,
-1.228070, -1.228070,
-1.043478, -1.043478 ;

fates_hydr_pitlp_node =
-1.67, -1.67,
Expand Down Expand Up @@ -1004,7 +1004,7 @@ data:

fates_fire_nignitions = 15 ;

fates_hydr_kmax_rsurf = 0.001;
fates_hydr_kmax_rsurf = 20;

fates_hydr_psi0 = 0 ;

Expand Down

0 comments on commit af76d76

Please sign in to comment.