| Version 31 (modified by dgoll, 10 years ago) (diff) |
|---|
Daniel's page
CNP-Dev version based on MERGE-OCN revision 2698
Bugs fixed
1. Mass conservation issue: stomate_growth_fun_all.f90
The carbon being allocated to biomass pools must no be substracted before nutrient limitation of allocation is computed. This can be fixed by following:
!DSGdebug_01
! do NOT this now, do it after nutrient limitation on allocation is considered in bm_alloc_tot(ipts,j)
! ! Update the labile carbon pool
! biomass(ipts,j,ilabile,icarbon) = biomass(ipts,j,ilabile,icarbon) - &
! bm_alloc_tot(ipts,j)
!DSGdebug_01
!! 3.10 Maintenance respiration
! The calculation of ::resp_maint is solely based on the demand i.e.
! given the biomass and the condition of the plant, how much should be
! respired. It is not sure that this demand can be satisfied i.e. the
! calculated maintenance respiration may exceed the available carbon
!DSGdebug_01
! DEFAULT CASE: There is no deficit which must be subtracted from labile
deficit = zero
!DSGdebug_01
IF ( bm_alloc_tot(ipts,j) - resp_maint(ipts,j) .LT. zero ) THEN
[...]
! Not enough carbon to pay the deficit, the individual
! is going to die at the end of this day
bm_alloc_tot(ipts,j) = bm_alloc_tot(ipts,j) + &
biomass(ipts,j,icarbres,icarbon)
biomass(ipts,j,icarbres,icarbon) = zero
! Truncate the maintenance respiration to the available carbon
resp_maint(ipts,j) = bm_alloc_tot(ipts,j)
!DSGdebug_01
! There is no deficit which must be subtracted from labile
deficit = zero
!DSGdebug_01
ENDIF
[...]
! Final ::resp_maint is know
bm_alloc_tot(ipts,j) = bm_alloc_tot(ipts,j) - resp_maint(ipts,j)
!DSGdebug_01
! Subtracted the deficit from labile pool
biomass(ipts,j,ilabile,icarbon) = biomass(ipts,j,ilabile,icarbon) - &
(resp_maint(ipts,j) + deficit)
!DSGdebug_01
!! 3.11 Growth respiration
! Calculate total growth respiration and update allocatable carbon
! Growth respiration is a tax on productivity, not actual allocation
! Total growth respiration has be calculated before the allocation
! takes place because the allocation itself is not linear. After
! the allocation has been calculated, growth respiration can be
! calculated for each biomass component separatly. The unit of
! resp_growth is gC m-2 dt-1
resp_growth(ipts,j) = frac_growthresp(j) * MAX(zero, bm_alloc_tot(ipts,j))
bm_alloc_tot(ipts,j) = bm_alloc_tot(ipts,j) - resp_growth(ipts,j)
!DSGdebug_01
biomass(ipts,j,ilabile,icarbon) = biomass(ipts,j,ilabile,icarbon) - &
resp_growth(ipts,j)
!DSGdebug_01
[...]
!=======================================================
! Block from OCN but I did not find it in DOFOCO ???
!
! 5.1 retrieve allocated biomass from labile pool (nitrogen, or new allocation)
!
!DSGdebug_01
! This is the right spot to remove bm_alloc_tot:
biomass(:,:,ilabile,icarbon) = biomass(:,:,ilabile,icarbon) - bm_alloc_tot(:,:)
!DSGdebug_01
biomass(:,:,ilabile,initrogen) = biomass(:,:,ilabile,initrogen) - n_alloc_tot(:,:)
2. Mass conservation issue: stomate_phenology.f90
The nutrient demand must be calculated AFTER the final Cl_init and Cr_init are known. This can be fixed by:
! The biomass available to use is set to be the minimum of the biomass of
! the labile pool (if carbon not taken from the atmosphere), and
! the wanted biomass.
bm_use(i) = MIN( biomass(i,j,ilabile,icarbon) + biomass(i,j,icarbres,icarbon), &
bm_wanted(i) )
!DSGdebug_02 ! the nutrients need to support the biomass:
!DSGdebug_02 bm_wanted_n(i) = (Cl_init + Cr_init*fcn_root(j))/cn_leaf_prescribed(j)
[...]
! In case nitrogen or phosphorus is not sufficiently available
! downregulate new leaf biomass to respect leaf stoichiometry;
! DSG: this violates the ratio used to calculate the
! leave-root-sapwood relationships: is this OK?
!DSGdebug_02: moved after Cl_init and Cr_init are updated
! the nutrients need to support the biomass:
bm_wanted_n(i) = (Cl_init + Cr_init*fcn_root(j))/cn_leaf_prescribed(j)
!DSGdebug_02
3. Parameter value issue: constantes_mtc.f90
The parameter k_latosa_max and k_latosa_min were initially designed for tree PFT only. However these variables are also used for grass PFT (stomate_growth_fun_all.f90). Therefore these parameter cannot be set to undef. Parameter set to value of tree PFT.
!DSGdebug_03 REAL(r_std), PARAMETER, DIMENSION(nvmc) :: k_latosa_max_mtc = & !! Maximum leaf-to-sapwood area ratio as defined in McDowell et al & (/ undef, 5000., 5000., 5000., 3000., 5000., 5000., & !! 2002, Oecologia and compiled in Hickler et al 2006, Appendix S2 ! & 5000., 5000., undef, undef, undef, undef /) !! (unitless) & 5000., 5000., 5000., 5000., 5000., 5000. /) !! (unitless) REAL(r_std), PARAMETER, DIMENSION(nvmc) :: k_latosa_min_mtc = & !! Minimum leaf-to-sapwood area ratio as defined in McDowell et al & (/ undef, 1500., 1500., 1500., 1000., 1500., 1500., & !! 2002, Oecologia and compiled in Hickler et al 2006, Appendix S2 ! & 1500., 1500., undef, undef, undef, undef /) !! (unitless) & 1500., 1500., 1500., 1500., 1500., 1500. /) !! (unitless) !DSGdebug_03
4. Stoichiometry issue: stomate_turnover.f90
There must be Nitrogen losses from abovesap wood when there are Carbon losses to ensure CN ratio. Here I assume that sapwood nitrogen can be recycled by plants according to leaves (this can be discussed).
dturnover(:) = biomass(:,ivm,iroot,initrogen) * leaf_frac(:,ivm,ilage) * turnover_rate(:)
biomass(:,ivm,ilabile,initrogen) = biomass(:,ivm,ilabile,initrogen) + recycle_root(ivm) * dturnover(:)
biomass(:,ivm,iroot, initrogen) = biomass(:,ivm,iroot, initrogen) - dturnover(:)
turnover(:,ivm,iroot, initrogen) = turnover(:,ivm,iroot,initrogen) + ( un - recycle_root(ivm) ) * dturnover(:)
dturnover(:) = biomass(:,ivm,ifruit,initrogen) * leaf_frac(:,ivm,ilage) * turnover_rate(:)
biomass(:,ivm,ifruit, initrogen) = biomass(:,ivm,ifruit,initrogen) - dturnover(:)
turnover(:,ivm,ifruit,initrogen) = turnover(:,ivm,ifruit,initrogen) + dturnover(:)
!DSGdebug_04
dturnover(:) = biomass(:,ivm,isapabove,initrogen) * leaf_frac(:,ivm,ilage) * turnover_rate(:)
biomass(:,ivm,ilabile,initrogen) = biomass(:,ivm,ilabile,initrogen) + recycle_leaf(ivm) * dturnover(:)
biomass(:,ivm,isapabove,initrogen) = biomass(:,ivm,isapabove,initrogen) - dturnover(:)
turnover(:,ivm,isapabove,initrogen) = turnover(:,ivm,isapabove,initrogen) + ( un - recycle_leaf(ivm) ) * dturnover(:)
!DSGdebug_04
5. Negative biomass: stomate_growth_fun_all.f90
In case the variable "residual" is not zero and nutrient availability reduces bm_alloc_tot to a value smaller than the value of residual, the allocated biomass gets negative. This can be fixed by correcting the biomass to be allocated by the residual before calculating nutrient limitation.
! Move the unallocated carbon back into the labile pool
biomass(ipts,j,ilabile,icarbon) = &
biomass(ipts,j,ilabile,icarbon) + residual(ipts,j)
ENDIF
!DSGdebug_05
! correct the biomass to be allocated by the residual.
! this has to be done HERE as we need to known the actual biomass
! being allocated for the calculation of nutrient limitation
bm_alloc_tot(ipts,j) = bm_alloc_tot(ipts,j) - residual(ipts,j)
!DSGdebug_05
[...]
!DSGdebug_05
! bm_alloc(:,j,k,icarbon) = f_alloc(:,j,k) * (bm_alloc_tot(:,j) - residual(:,j))
bm_alloc(:,j,k,icarbon) = f_alloc(:,j,k) * bm_alloc_tot(:,j) ! residual was already removed from bm_alloc_tot
!DSGdebug_05
Further changes needed regarding residual to ensure mass conservation:
- The residual has to be initialized with zero:
! If npp is not initialized, bare soil value is never set. npp(:,:) = zero !DSGdebug_05a ! Not having this results in an unitilized error !DSGdebug_05a ! with valgrid, but I can't figure out why. It always !DSGdebug_05a ! seems to be set before being used. ! DSG: Of course you get errors, when a variable was never set but is used ! (for example in case first_call=true) ! I set it to zero, as it should be zero in the aformentioned case !DSGdebug_05a residual(:,:) = val_exp !DSGdebug_05a residual(:,:) = zero !DSGdebug_05a - In case of nutrient down regulating bm_alloc the labile pool has to be corrected for not allocated residual
!DSGdebug_01 biomass(:,:,ilabile,:) = biomass(:,:,ilabile,:) - alloc_tot(:,:,:) !DSGdebug_01 !DSGdebug_05a WHERE ((residual(:,:).GT. zero).AND.(residual(:,:) .GT. bm_alloc_tot(:,:))) biomass(:,:,ilabile,icarbon) = biomass(:,:,ilabile,icarbon) - (residual(:,:) - bm_alloc_tot(:,:)) END WHERE !DSGdebug_05a
6. Mass conservation issue: stomate_phenology.f90
When calculating the share of new biomass to leaves (Cl_init) and roots (Cr_init) of each respective circ_class, the biomass has to be divided by the number of tree per class (circ_class_n); which wasn't done.
! Distribute the biomass over the leaves and roots (gC tree-1)
! Use Cl_init to estimate the share for each circumference class
! leaf biomass = FK * Cs / height (allometric relationship)
! root biomass = KF / LF * Cs / height
! Convert from gC m-2 to gC tree-1
! +++CHECK+++
! Cl_init + Cr_init can exceed bm_use. bm_use should be used in these equations
! DSG: I cannot confirm the statement in the line before
! but the equation were wrong, the biomass has to be divided by numbers of tree (circ_class_n)
circ_class_biomass(i,j,l,ileaf,icarbon) = circ_class_biomass(i,j,l,ileaf,icarbon) + &
Cl_init * ( KF(i,j) * Cs_tree(l) / height(i,j,l) * circ_class_n(i,j,l) ) / &
!DSGdebug_06 SUM( KF(i,j) * Cs_tree(:) / height(i,j,:) * circ_class_n(i,j,:) )
SUM( KF(i,j) * Cs_tree(:) / height(i,j,:) * circ_class_n(i,j,:) ) / circ_class_n(i,j,l)
!DSGdebug_06
circ_class_biomass(i,j,l,iroot,icarbon) = circ_class_biomass(i,j,l,iroot,icarbon) + &
Cr_init * ( KF(i,j) * Cs_tree(l) / height(i,j,l) * circ_class_n(i,j,l) ) / &
!DSGdebug_06 SUM( KF(i,j) * Cs_tree(:) / height(i,j,:) * circ_class_n(i,j,:) )
SUM( KF(i,j) * Cs_tree(:) / height(i,j,:) * circ_class_n(i,j,:) ) / circ_class_n(i,j,l)
!DSGdebug_06
!++++++++++++
7. Mass conservation issue: stomate_phenology.f90
see in code comment for explanation
! Calculate the available biomass in roots, sapwood and leaves (gC ind-1)
IF ( biomass(i,j,ileaf,icarbon) .LT. min_stomate .AND. biomass(i,j,iroot,icarbon) .LT. min_stomate) THEN
Cs_grass = biomass(i,j,isapabove,icarbon)
ELSE
WRITE(6,*) 'There is leaf or root carbon that should not be here, something could have gone wrong in senescence'
!DSGdebug_07: that doesn't justify to violate mass conservation
Cs_grass = biomass(i,j,isapabove,icarbon)
ENDIF
! Cl_init and Cr_init were recalculated to properly account for bm_use (the available C)
IF (j==test_pft) WRITE (6,*) 'Cr_init (end)=',Cr_init
!DSGdebug_07 biomass(i,j,ileaf,icarbon) = Cl_init
!DSGdebug_07 biomass(i,j,iroot,icarbon) = Cr_init
biomass(i,j,ileaf,icarbon) = biomass(i,j,ileaf,icarbon) + Cl_init
biomass(i,j,iroot,icarbon) = biomass(i,j,iroot,icarbon) + Cr_init
!DSGdebug_07
Conceptual problems
1. hydraulic architecture model for grass PFTs
The use of the hydraulic architecture model of Hickler et al (2006) designed for tree for grass seems problematic. Technically, the parameter "k_sap" is used in function "calculate_c0_allo"c but not defined for non-Tree PFTs resulting in "NaN" for calc_c0_alloc. I did a dirty fix by setting "c0_alloc= 2.0" to avoid NaN for non-Tree PFTs.
Mass conservation checks
Mass closure given by: mass_before + mass_change = mass_after
stomate_lpj.f90
!DSG mass conservation ========================================
mass_before(:,:,:) = SUM(biomass(:,:,:,:),DIM=3)
!! 5. Grow new biomass - respiration, npp and allocation
! Call the allometry based allocation (based on Sitch et al 2003 and Zaehle et al 2010)
CALL growth_fun_all (npts, dt_days, veget_max, PFTpresent, &
senescence, when_growthinit, moiavail_growingseason, t2m_week, &
gpp_daily, gpp_week, resp_maint_part, resp_maint, &
resp_growth, npp_daily, biomass, age, &
leaf_age, leaf_frac, use_reserve, t_photo_stress, &
lab_fac, lai_target, ind, rue_longterm, &
circ_class_n, circ_class_biomass, c0_alloc, cn_leaf_season, np_leaf_season, &
KF, n_uptake_daily, p_uptake_daily)
!DSG mass conservation ============================================
mass_change(:,:,icarbon) = npp_daily(:,:)*dt_days
mass_change(:,:,initrogen) = SUM(n_uptake_daily(:,:,:),DIM=3)
mass_change(:,:,iphosphorus) = p_uptake_daily(:,:)
!DSG mass conservation ========================================
mass_before(:,:,:) = SUM(biomass(:,:,:,:),DIM=3)
CALL phenology_prognostic (npts, dt_days, PFTpresent, veget_max, &
tlong_ref, t2m_month, t2m_week, gpp_daily, &
maxmoiavail_lastyear, minmoiavail_lastyear, moiavail_month, moiavail_week, &
gdd_m5_dormance, gdd_midwinter, ncd_dormance, ngd_minus5, &
senescence, time_hum_min, biomass, leaf_frac, &
leaf_age, when_growthinit, co2_to_bm, circ_class_n, &
circ_class_biomass, ind, c0_alloc, KF)
!DSG mass conservation ============================================
mass_change(:,:,icarbon) = zero
mass_change(:,:,initrogen) = zero
mass_change(:,:,iphosphorus) = zero
