Skip to content

Commit

Permalink
Refactoring: factorize fast_fluxes in cohort
Browse files Browse the repository at this point in the history
  • Loading branch information
marcadella committed Jan 9, 2025
1 parent 3a39dc7 commit 08591f8
Show file tree
Hide file tree
Showing 4 changed files with 65 additions and 74 deletions.
60 changes: 33 additions & 27 deletions src/datatypes_biomee.mod.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -109,13 +109,10 @@ module datatypes_biomee
type(orgpool) :: plabl ! labile pool, temporary storage of N and C, kg tree-1

!===== Fast step fluxes, kg timestep-1 tree-1
real :: gpp = 0.0 ! gross primary productivity, kg C timestep-1 tree-1
real :: npp = 0.0 ! net primary productivity, kg C timestep-1 tree-1
real :: resp = 0.0 ! plant respiration, kg C timestep-1 tree-1
type(common_fluxes) :: fast_fluxes

real :: resl = 0.0 ! leaf respiration, kg C timestep-1 tree-1
real :: resr = 0.0 ! root respiration, kg C timestep-1 tree-1
real :: N_uptake = 0.0 ! N uptake at each step per tree, kg N timestep-1 tree-1
real :: fixedN = 0.0 ! fixed N at each step per tree, kg N timestep-1 tree-1

!===== Daily fluxes, kg day-1 tree-1
type(common_fluxes) :: daily_fluxes
Expand All @@ -142,7 +139,6 @@ module datatypes_biomee
real :: rootareaL(max_lev) = 0.0 ! Root length per layer, m of root/m
real :: WupL(max_lev) = 0.0 ! normalized vertical distribution of uptake
real :: W_supply = dummy ! potential water uptake rate per unit time per tree
real :: transp = dummy ! water transpired per tree per timestep, kg H2O timestep-1 tree-1

!===== Photosynthesis variables
real :: An_op = 0.0 ! mol C/(m2 of leaf per year)
Expand Down Expand Up @@ -207,7 +203,7 @@ module datatypes_biomee
real :: m_turnover = 0.0 ! C turnover due to mortality and tissue turnover (kg C m-2 yr-1)

!===== Leaf area index
real :: LAI ! leaf area index
real :: LAI ! leaf area index (surface of leaves per m2 of ground/tile)
real :: CAI ! crown area index

!===== Averaged quantities for PPA phenology
Expand All @@ -219,6 +215,7 @@ module datatypes_biomee
real :: totN = 0.0
real :: N_input = 0.0 ! annual N input (kg N m-2 yr-1)
real :: N_uptake = 0.0 ! N uptake at each time step, kg N m-2 timestep-1
real :: fixedN = 0.0 ! fixed N at each time step, kg N m-2 timestep-1
real :: annualN = 0.0 ! annual available N in a year
real :: Nloss_yr = 0.0 ! annual N loss
real :: N_P2S_yr = 0.0 ! N turnover (plant to soil) (kg N m-2 yr-1)
Expand Down Expand Up @@ -351,13 +348,12 @@ subroutine Zero_diagnostics(vegn)

cc%C_growth = 0.0
cc%N_growth = 0.0
cc%gpp = 0.0
cc%npp = 0.0
cc%resp = 0.0

cc%resl = 0.0
cc%resr = 0.0
cc%resg = 0.0
cc%transp = 0.0

cc%fast_fluxes = common_fluxes()

! Reset daily
cc%daily_fluxes = common_fluxes()
Expand Down Expand Up @@ -481,35 +477,45 @@ subroutine hourly_diagnostics(vegn, forcing)
vegn%age = vegn%age + myinterface%dt_fast_yr

! Tile summary
vegn%GPP = 0
vegn%NPP = 0
vegn%Resp = 0
vegn%transp = 0.0
vegn%GPP = 0.0
vegn%NPP = 0.0
vegn%Resp = 0.0
vegn%N_uptake = 0.0
vegn%fixedN = 0.0

do i = 1, vegn%n_cohorts
cc => vegn%cohorts(i)

! cohort daily
call update_fluxes(cc%daily_fluxes, common_fluxes(cc%transp, cc%gpp, cc%Npp, cc%Resp, cc%N_uptake, cc%fixedN))
call update_fluxes(cc%daily_fluxes, cc%fast_fluxes)

! Tile hourly
vegn%GPP = vegn%GPP + cc%gpp * cc%nindivs
vegn%NPP = vegn%NPP + cc%Npp * cc%nindivs
vegn%Resp = vegn%Resp + cc%Resp * cc%nindivs
vegn%transp = vegn%transp + cc%fast_fluxes%trsp * cc%nindivs
vegn%GPP = vegn%GPP + cc%fast_fluxes%gpp * cc%nindivs
vegn%NPP = vegn%NPP + cc%fast_fluxes%Npp * cc%nindivs
vegn%Resp = vegn%Resp + cc%fast_fluxes%Resp * cc%nindivs
vegn%N_uptake = vegn%N_uptake + cc%fast_fluxes%Nup * cc%nindivs
vegn%fixedN = vegn%fixedN + cc%fast_fluxes%fixedN * cc%nindivs

! Reset fast fluxes
cc%fast_fluxes = common_fluxes()
enddo

! NEP is equal to NNP minus soil respiration
vegn%nep = vegn%npp - vegn%rh

! Daily summary:
vegn%dailyNup = vegn%dailyNup + vegn%N_uptake
vegn%dailyGPP = vegn%dailyGPP + vegn%gpp
vegn%dailyNPP = vegn%dailyNPP + vegn%npp
vegn%dailyResp = vegn%dailyResp + vegn%resp
vegn%dailyRh = vegn%dailyRh + vegn%rh
vegn%dailyTrsp = vegn%dailyTrsp + vegn%transp
vegn%dailyEvap = vegn%dailyEvap + vegn%evap
vegn%dailyRoff = vegn%dailyRoff + vegn%runoff
vegn%dailyPrcp = vegn%dailyPrcp + forcing%rain * myinterface%step_seconds
vegn%dailyNup = vegn%dailyNup + vegn%N_uptake
vegn%dailyfixedN = vegn%dailyfixedN + vegn%fixedN
vegn%dailyGPP = vegn%dailyGPP + vegn%gpp
vegn%dailyNPP = vegn%dailyNPP + vegn%npp
vegn%dailyResp = vegn%dailyResp + vegn%resp
vegn%dailyRh = vegn%dailyRh + vegn%rh
vegn%dailyTrsp = vegn%dailyTrsp + vegn%transp
vegn%dailyEvap = vegn%dailyEvap + vegn%evap
vegn%dailyRoff = vegn%dailyRoff + vegn%runoff
vegn%dailyPrcp = vegn%dailyPrcp + forcing%rain * myinterface%step_seconds

end subroutine hourly_diagnostics

Expand Down
18 changes: 9 additions & 9 deletions src/gpp_biomee.mod.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -172,17 +172,17 @@ subroutine gpp( forcing, vegn, first_simu_step )
cc%An_op = psyn ! molC s-1 m-2 of leaves ! net photosynthesis, mol C/(m2 of leaves s)
cc%An_cl = -resp ! molC s-1 m-2 of leaves
cc%w_scale = w_scale2
cc%transp = transp * h2o_molmass * 1e-3 * cc%leafarea * myinterface%step_seconds ! Transpiration (kgH2O/(tree step), Weng, 2017-10-16
cc%resl = -resp * c_molmass * 1e-3 * cc%leafarea * myinterface%step_seconds ! kgC tree-1 step-1
cc%gpp = (psyn - resp) * c_molmass * 1e-3 * cc%leafarea * myinterface%step_seconds ! kgC step-1 tree-1
cc%fast_fluxes%trsp = transp * h2o_molmass * 1e-3 * cc%leafarea * myinterface%step_seconds ! Transpiration (kgH2O/(tree step), Weng, 2017-10-16
cc%resl = -resp * c_molmass * 1e-3 * cc%leafarea * myinterface%step_seconds ! kgC tree-1 step-1
cc%fast_fluxes%gpp = (psyn - resp) * c_molmass * 1e-3 * cc%leafarea * myinterface%step_seconds ! kgC step-1 tree-1

else

! no leaves means no photosynthesis and no stomatal conductance either
cc%An_op = 0.0
cc%An_cl = 0.0
cc%gpp = 0.0
cc%transp = 0.0
cc%fast_fluxes%gpp = 0.0
cc%fast_fluxes%trsp = 0.0
cc%w_scale = dummy

endif
Expand Down Expand Up @@ -264,11 +264,11 @@ subroutine gpp( forcing, vegn, first_simu_step )
! irrelevant variables for this setup
cc%An_op = 0.0
cc%An_cl = 0.0
cc%transp = 0.0
cc%fast_fluxes%trsp = 0.0
cc%w_scale = dummy

! quantities per tree and cumulated over seconds in time step (kgC step-1 tree-1 )
cc%gpp = par * fapar_tree(i) * out_pmodel%lue * cc%crownarea * myinterface%step_seconds * 1.0e-3
cc%fast_fluxes%gpp = par * fapar_tree(i) * out_pmodel%lue * cc%crownarea * myinterface%step_seconds * 1.0e-3
cc%resl = fapar_tree(i) * out_pmodel%vcmax25 * params_gpp%rd_to_vcmax * calc_ftemp_inst_rd( forcing%Tair - kTkelvin ) &
* cc%crownarea * myinterface%step_seconds * c_molmass * 1.0e-3

Expand All @@ -277,10 +277,10 @@ subroutine gpp( forcing, vegn, first_simu_step )
! no leaves means no photosynthesis and no stomatal conductance either
cc%An_op = 0.0
cc%An_cl = 0.0
cc%transp = 0.0
cc%w_scale = dummy
cc%resl = 0.0
cc%gpp = 0.0
cc%fast_fluxes%gpp = 0.0
cc%fast_fluxes%trsp = 0.0

endif

Expand Down
12 changes: 1 addition & 11 deletions src/soil_biomee.mod.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -105,7 +105,6 @@ subroutine SoilWaterDynamicsLayer(forcing,vegn) !outputs

! Water uptaken by roots, per timestep
WaterBudgetL = 0.0
vegn%transp = 0.0
do j = 1, vegn%n_cohorts
cc => vegn%cohorts(j)
! Compare with soil water
Expand All @@ -115,10 +114,9 @@ subroutine SoilWaterDynamicsLayer(forcing,vegn) !outputs
if(cc%W_supply>0.0)then
do i=1,max_lev
fsupply = cc%WupL(i)/cc%W_supply
transp = fsupply * cc%transp * cc%nindivs
transp = fsupply * cc%fast_fluxes%trsp * cc%nindivs
!vegn%wcl(i) = vegn%wcl(i) - transp/(thksl(i)*1000.0)
WaterBudgetL(i) = WaterBudgetL(i) - transp
vegn%transp = vegn%transp + transp
enddo
endif
enddo ! all cohorts
Expand All @@ -128,16 +126,8 @@ subroutine SoilWaterDynamicsLayer(forcing,vegn) !outputs
kappa = 0.75
! thermodynamic parameters for air

! print*, forcing%radiation, kappa, vegn%LAI ! xxx debug

Rsoilabs = forcing%radiation * exp(-kappa*vegn%LAI)

! print*,'forcing%radiation', forcing%radiation
! print*,'kappa', kappa
! print*,'vegn%LAI', vegn%LAI
! print*,'Rsoilabs', Rsoilabs
! print*, 'transp', transp

Hgrownd = 0.0
TairK = forcing%Tair
Tair = forcing%Tair - 273.16
Expand Down
49 changes: 22 additions & 27 deletions src/vegetation_biomee.mod.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -61,13 +61,13 @@ subroutine vegn_CNW_budget( vegn, forcing, first_simu_step, tsoil )
call plant_respiration( cc, forcing%tair ) ! get resp per tree per time step

! We add the growth respiration scaled from daily to timestap
cc%resp = cc%resp + (cc%resg * myinterface%step_seconds) / secs_per_day
cc%resp = cc%resp * myinterface%params_tile%tf_base ! scaling for calibration
cc%npp = cc%gpp - cc%resp ! kgC tree-1 step-1
cc%fast_fluxes%resp = cc%fast_fluxes%resp + (cc%resg * myinterface%step_seconds) / secs_per_day
cc%fast_fluxes%resp = cc%fast_fluxes%resp * myinterface%params_tile%tf_base ! scaling for calibration
cc%fast_fluxes%npp = cc%fast_fluxes%gpp - cc%fast_fluxes%resp ! kgC tree-1 step-1

! detach photosynthesis model from plant growth
cc%plabl%c%c12 = cc%plabl%c%c12 + cc%npp
cc%plabl%n%n14 = cc%plabl%n%n14 + cc%fixedN
cc%plabl%c%c12 = cc%plabl%c%c12 + cc%fast_fluxes%npp
cc%plabl%n%n14 = cc%plabl%n%n14 + cc%fast_fluxes%fixedN

enddo ! all cohorts

Expand Down Expand Up @@ -122,13 +122,13 @@ subroutine plant_respiration( cc, tairK )
!endif

! Obligate Nitrogen Fixation
cc%fixedN = fnsc*spdata(sp)%NfixRate0 * cc%proot%c%c12 * tf * myinterface%dt_fast_yr ! kgN tree-1 step-1
r_Nfix = spdata(sp)%NfixCost0 * cc%fixedN ! + 0.25*spdata(sp)%NfixCost0 * cc%N_uptake ! tree-1 step-1
cc%fast_fluxes%fixedN = fnsc*spdata(sp)%NfixRate0 * cc%proot%c%c12 * tf * myinterface%dt_fast_yr ! kgN tree-1 step-1
r_Nfix = spdata(sp)%NfixCost0 * cc%fast_fluxes%fixedN ! + 0.25*spdata(sp)%NfixCost0 * cc%N_uptake ! tree-1 step-1

! LeafN = spdata(sp)%LNA * cc%leafarea ! gamma_SW is sapwood respiration rate (kgC m-2 Acambium yr-1)
r_stem = fnsc*spdata(sp)%gamma_SW * Acambium * tf * myinterface%dt_fast_yr ! kgC tree-1 step-1
r_root = fnsc*spdata(sp)%gamma_FR * cc%proot%n%n14 * tf * myinterface%dt_fast_yr ! root respiration ~ root N
cc%resp = cc%resl + r_stem + r_root + r_Nfix !kgC tree-1 step-1
cc%fast_fluxes%resp = cc%resl + r_stem + r_root + r_Nfix !kgC tree-1 step-1
cc%resr = r_root + r_Nfix ! tree-1 step-1

end associate
Expand Down Expand Up @@ -182,7 +182,7 @@ end subroutine fetch_CN_for_growth

subroutine vegn_growth_EW( vegn )
!////////////////////////////////////////////////////////////////
! updates cohort biomass pools, LAI, and height using accumulated
! updates cohort biomass pools, LAI, and height using accumulated
! C_growth and bHW_gain
! Code from BiomeE-Allocation
!---------------------------------------------------------------
Expand All @@ -209,15 +209,13 @@ subroutine vegn_growth_EW( vegn )
real :: LF_deficit, FR_deficit
real :: N_demand,Nsupplyratio,extraN
real :: r_N_SD
logical :: do_editor_scheme = .False.
integer :: i
do_editor_scheme = .False. ! .True.

! Turnover of leaves and fine roots
call vegn_tissue_turnover( vegn )

! Allocate C_gain to tissues
vegn%LAI = 0.0 ! added here, otherwise LAI is the sum of current and new, Beni 24 Aug 2022
vegn%LAI = 0.0
do i = 1, vegn%n_cohorts
cc => vegn%cohorts(i)

Expand Down Expand Up @@ -1394,20 +1392,18 @@ subroutine vegn_N_uptake(vegn, tsoil)
! It considers competition here. How much N one can absorp depends on
! how many roots it has and how many roots other individuals have.
N_Roots = 0.0
vegn%N_uptake = 0.0

do i = 1, vegn%n_cohorts
cc => vegn%cohorts(i)
associate (sp => myinterface%params_species(cc%species))

cc%N_uptake = 0.0 ! We MUST ALWAYS reset it (whether or not vegn%ninorg%n14 > 0.0)
cc%NSNmax = sp%fNSNmax*(cc%bl_max/(sp%CNleaf0*sp%leafLS)+cc%br_max/sp%CNroot0)
if (cc%plabl%n%n14 < cc%NSNmax) N_Roots = N_Roots + cc%proot%c%c12 * cc%nindivs
if (vegn%ninorg%n14 > 0.0) then

end associate
enddo
do i = 1, vegn%n_cohorts
cc => vegn%cohorts(i)
associate (sp => myinterface%params_species(cc%species))

if (vegn%ninorg%n14 > 0.0) then
cc%NSNmax = sp%fNSNmax*(cc%bl_max/(sp%CNleaf0*sp%leafLS)+cc%br_max/sp%CNroot0)
if (cc%plabl%n%n14 < cc%NSNmax) N_Roots = N_Roots + cc%proot%c%c12 * cc%nindivs

end associate
enddo

! M-M equation for Nitrogen absoption, McMurtrie et al. 2012, Ecology & Evolution
! rate at given root biomass and period of time
Expand All @@ -1426,12 +1422,11 @@ subroutine vegn_N_uptake(vegn, tsoil)
cc => vegn%cohorts(i)
if (cc%plabl%n%n14 < cc%NSNmax) then

cc%N_uptake = cc%proot%c%c12 * avgNup ! min(cc%proot%c%c12*avgNup, cc%NSNmax-cc%plabl%n%n14)
cc%plabl%n%n14 = cc%plabl%n%n14 + cc%N_uptake
cc%fast_fluxes%Nup = cc%proot%c%c12 * avgNup ! min(cc%proot%c%c12*avgNup, cc%NSNmax-cc%plabl%n%n14)
cc%plabl%n%n14 = cc%plabl%n%n14 + cc%fast_fluxes%Nup

! subtract N from mineral N
vegn%ninorg%n14 = vegn%ninorg%n14 - cc%N_uptake * cc%nindivs
vegn%N_uptake = vegn%N_uptake + cc%N_uptake * cc%nindivs
vegn%ninorg%n14 = vegn%ninorg%n14 - cc%fast_fluxes%Nup * cc%nindivs
endif
enddo
cc => null()
Expand Down

0 comments on commit 08591f8

Please sign in to comment.