diff --git a/biogeochem/EDPatchDynamicsMod.F90 b/biogeochem/EDPatchDynamicsMod.F90 index 6d7c84ebb1..df283d6c54 100755 --- a/biogeochem/EDPatchDynamicsMod.F90 +++ b/biogeochem/EDPatchDynamicsMod.F90 @@ -831,6 +831,7 @@ subroutine create_patch(currentSite, new_patch, age, areap, spread_local,cwd_ag_ new_patch%siteptr => currentSite new_patch%age = age + new_patch%age_class = 1 new_patch%area = areap new_patch%spread = spread_local new_patch%cwd_ag = cwd_ag_local @@ -899,6 +900,7 @@ subroutine zero_patch(cp_p) currentPatch%clm_pno = 999 currentPatch%age = nan + currentPatch%age_class = 1 currentPatch%area = nan currentPatch%canopy_layer_lai(:) = nan currentPatch%total_canopy_area = nan @@ -935,8 +937,6 @@ subroutine zero_patch(cp_p) currentPatch%lai = nan ! leaf area index of patch currentPatch%spread(:) = nan ! dynamic ratio of dbh to canopy area. currentPatch%pft_agb_profile(:,:) = nan - currentPatch%gpp = 0._r8 - currentPatch%npp = 0._r8 ! DISTURBANCE currentPatch%disturbance_rates = 0._r8 @@ -1150,6 +1150,7 @@ subroutine fuse_2_patches(dp, rp) ! associated with the secnd patch ! ! !USES: + use EDTypesMod, only: ageclass_ed ! ! !ARGUMENTS: type (ed_patch_type) , intent(inout), pointer :: dp ! Donor Patch @@ -1169,6 +1170,7 @@ subroutine fuse_2_patches(dp, rp) !area weighted average of ages & litter rp%age = (dp%age * dp%area + rp%age * rp%area)/(dp%area + rp%area) + rp%age_class = count(rp%age-ageclass_ed.ge.0.0_r8) do p = 1,numpft_ed rp%seeds_in(p) = (rp%seeds_in(p)*rp%area + dp%seeds_in(p)*dp%area)/(rp%area + dp%area) diff --git a/main/EDMainMod.F90 b/main/EDMainMod.F90 index 6882222f1e..2bddac71cf 100755 --- a/main/EDMainMod.F90 +++ b/main/EDMainMod.F90 @@ -143,6 +143,7 @@ subroutine ed_integrate_state_variables(currentSite, bc_in ) ! FIX(SPM,032414) refactor so everything goes through interface ! ! !USES: + use EDTypesMod, only : ageclass_ed ! ! !ARGUMENTS: type(ed_site_type) , intent(inout) :: currentSite @@ -177,6 +178,9 @@ subroutine ed_integrate_state_variables(currentSite, bc_in ) currentPatch%patchno,currentPatch%area endif + ! check to see if the patch has moved to the next age class + currentPatch%age_class = count(currentPatch%age-ageclass_ed.ge.0.0_r8) + ! Find the derivatives of the growth and litter processes. call canopy_derivs(currentSite, currentPatch, bc_in) diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index 1817de3e66..2470c8775e 100755 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -83,9 +83,12 @@ module EDTypesMod ! real(r8), parameter, dimension(16) :: sclass_ed = (/0.0_r8,1.0_r8,2.0_r8,3.0_r8,4.0_r8,5.0_r8,10.0_r8,20.0_r8,30.0_r8,40.0_r8, & ! 50.0_r8,60.0_r8,70.0_r8,80.0_r8,90.0_r8,100.0_r8/) - real(r8), parameter, dimension(13) :: sclass_ed = (/0.0_r8,5.0_r8,10.0_r8,15.0_r8,20.0_r8,30.0_r8,40.0_r8, & + real(r8), parameter, dimension(nlevsclass_ed) :: sclass_ed = (/0.0_r8,5.0_r8,10.0_r8,15.0_r8,20.0_r8,30.0_r8,40.0_r8, & 50.0_r8,60.0_r8,70.0_r8,80.0_r8,90.0_r8,100.0_r8/) + integer, parameter :: nlevage_ed = 7 ! Number of patch-age classes for age structured analyses + real(r8), parameter, dimension(nlevage_ed) :: ageclass_ed = (/0.0_r8,1.0_r8,2._r8,5.0_r8,10.0_r8,20.0_r8,50.0_r8/) + ! integer, parameter :: nlevsclass_ed = 17 ! real(r8), parameter, dimension(17) :: sclass_ed = (/0.1_r8, 5.0_r8,10.0_r8,15.0_r8,20.0_r8,25.0_r8, & @@ -99,13 +102,15 @@ module EDTypesMod character(len = 10), parameter,dimension(5) :: char_list = (/"background","hydraulic ","carbon ","impact ","fire "/) - ! These three vectors are used for history output mapping + ! These vectors are used for history output mapping real(r8) ,allocatable :: levsclass_ed(:) ! The lower bound on size classes for ED trees. This ! is used really for IO into the ! history tapes. It gets copied from ! the parameter array sclass_ed. integer , allocatable :: pft_levscpf_ed(:) integer , allocatable :: scls_levscpf_ed(:) + real(r8), allocatable :: levage_ed(:) + integer , allocatable :: levpft_ed(:) ! Control Parameters (cp_) @@ -325,6 +330,7 @@ module EDTypesMod ! PATCH INFO real(r8) :: age ! average patch age: years + integer :: age_class ! age class of the patch for history binning purposes real(r8) :: area ! patch area: m2 integer :: countcohorts ! Number of cohorts in patch integer :: ncl_p ! Number of occupied canopy layers @@ -388,8 +394,6 @@ module EDTypesMod ! PHOTOSYNTHESIS real(r8) :: psn_z(cp_nclmax,numpft_ed,cp_nlevcan) ! carbon assimilation in each canopy layer, pft, and leaf layer. umolC/m2/s - real(r8) :: gpp ! total patch gpp: KgC/m2/year - real(r8) :: npp ! total patch npp: KgC/m2/year ! ROOTS real(r8), allocatable :: rootfr_ft(:,:) ! root fraction of each PFT in each soil layer:- @@ -585,11 +589,20 @@ subroutine ed_hist_scpfmaps allocate( levsclass_ed(1:nlevsclass_ed )) allocate( pft_levscpf_ed(1:nlevsclass_ed*mxpft)) allocate(scls_levscpf_ed(1:nlevsclass_ed*mxpft)) + allocate( levpft_ed(1:mxpft )) + allocate( levage_ed(1:nlevage_ed )) ! Fill the IO array of plant size classes ! For some reason the history files did not like ! a hard allocation of sclass_ed levsclass_ed(:) = sclass_ed(:) + + levage_ed(:) = ageclass_ed(:) + + ! make pft array + do ipft=1,mxpft + levpft_ed(ipft) = ipft + end do ! Fill the IO arrays that match pft and size class to their combined array i=0 diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index b2b090b24a..568b9950d4 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -63,16 +63,7 @@ module FatesHistoryInterfaceMod integer, private :: ih_maint_resp_pa integer, private :: ih_growth_resp_pa - ! Indices to (patch x pft) variables (using nlevgrnd as surrogate) - - integer, private :: ih_biomass_pa_pft - integer, private :: ih_leafbiomass_pa_pft - integer, private :: ih_storebiomass_pa_pft - integer, private :: ih_nindivs_pa_pft - ! Indices to (site) variables - - integer, private :: ih_nep_si integer, private :: ih_nep_timeintegrated_si integer, private :: ih_npp_timeintegrated_si @@ -134,10 +125,28 @@ module FatesHistoryInterfaceMod integer, private :: ih_ar_crootm_si_scpf integer, private :: ih_ar_frootm_si_scpf + ! indices to (site x scls) variables + integer, private :: ih_ba_si_scls + + ! indices to (site x pft) variables + integer, private :: ih_biomass_si_pft + integer, private :: ih_leafbiomass_si_pft + integer, private :: ih_storebiomass_si_pft + integer, private :: ih_nindivs_si_pft + + + ! indices to (site x patch-age) variables + integer, private :: ih_area_si_age + integer, private :: ih_lai_si_age + integer, private :: ih_canopy_area_si_age + integer, private :: ih_gpp_si_age + integer, private :: ih_npp_si_age + integer, private :: ih_ncl_si_age + integer, private :: ih_npatches_si_age ! The number of variable dim/kind types we have defined (static) - integer, parameter :: fates_history_num_dimensions = 4 - integer, parameter :: fates_history_num_dim_kinds = 6 + integer, parameter :: fates_history_num_dimensions = 7 + integer, parameter :: fates_history_num_dim_kinds = 9 @@ -171,6 +180,7 @@ module FatesHistoryInterfaceMod type(iovar_map_type), pointer :: iovar_map(:) integer, private :: patch_index_, column_index_, levgrnd_index_, levscpf_index_ + integer, private :: levscls_index_, levpft_index_, levage_index_ contains procedure, public :: Init @@ -188,6 +198,9 @@ module FatesHistoryInterfaceMod procedure, public :: column_index procedure, public :: levgrnd_index procedure, public :: levscpf_index + procedure, public :: levscls_index + procedure, public :: levpft_index + procedure, public :: levage_index ! private work functions procedure, private :: define_history_vars @@ -200,6 +213,9 @@ module FatesHistoryInterfaceMod procedure, private :: set_column_index procedure, private :: set_levgrnd_index procedure, private :: set_levscpf_index + procedure, private :: set_levscls_index + procedure, private :: set_levpft_index + procedure, private :: set_levage_index end type fates_history_interface_type @@ -212,6 +228,7 @@ module FatesHistoryInterfaceMod subroutine Init(this, num_threads, fates_bounds) use FatesIODimensionsMod, only : patch, column, levgrnd, levscpf + use FatesIODimensionsMod, only : levscls, levpft, levage use FatesIODimensionsMod, only : fates_bounds_type implicit none @@ -240,7 +257,22 @@ subroutine Init(this, num_threads, fates_bounds) dim_count = dim_count + 1 call this%set_levscpf_index(dim_count) call this%dim_bounds(dim_count)%Init(levscpf, num_threads, & + fates_bounds%sizepft_class_begin, fates_bounds%sizepft_class_end) + + dim_count = dim_count + 1 + call this%set_levscls_index(dim_count) + call this%dim_bounds(dim_count)%Init(levscls, num_threads, & + fates_bounds%size_class_begin, fates_bounds%size_class_end) + + dim_count = dim_count + 1 + call this%set_levpft_index(dim_count) + call this%dim_bounds(dim_count)%Init(levpft, num_threads, & fates_bounds%pft_class_begin, fates_bounds%pft_class_end) + + dim_count = dim_count + 1 + call this%set_levage_index(dim_count) + call this%dim_bounds(dim_count)%Init(levage, num_threads, & + fates_bounds%age_class_begin, fates_bounds%age_class_end) ! FIXME(bja, 2016-10) assert(dim_count == FatesHistorydimensionmod::num_dimension_types) ! Allocate the mapping between FATES indices and the IO indices @@ -275,9 +307,21 @@ subroutine SetThreadBoundsEach(this, thread_index, thread_bounds) thread_bounds%ground_begin, thread_bounds%ground_end) index = this%levscpf_index() + call this%dim_bounds(index)%SetThreadBounds(thread_index, & + thread_bounds%sizepft_class_begin, thread_bounds%sizepft_class_end) + + index = this%levscls_index() + call this%dim_bounds(index)%SetThreadBounds(thread_index, & + thread_bounds%size_class_begin, thread_bounds%size_class_end) + + index = this%levpft_index() call this%dim_bounds(index)%SetThreadBounds(thread_index, & thread_bounds%pft_class_begin, thread_bounds%pft_class_end) + index = this%levage_index() + call this%dim_bounds(index)%SetThreadBounds(thread_index, & + thread_bounds%age_class_begin, thread_bounds%age_class_end) + end subroutine SetThreadBoundsEach ! =================================================================================== @@ -285,6 +329,7 @@ subroutine assemble_history_output_types(this) use FatesIOVariableKindMod, only : patch_r8, patch_ground_r8, patch_size_pft_r8 use FatesIOVariableKindMod, only : site_r8, site_ground_r8, site_size_pft_r8 + use FatesIOVariableKindMod, only : site_size_r8, site_pft_r8, site_age_r8 implicit none @@ -308,6 +353,15 @@ subroutine assemble_history_output_types(this) call this%set_dim_indices(site_size_pft_r8, 1, this%column_index()) call this%set_dim_indices(site_size_pft_r8, 2, this%levscpf_index()) + call this%set_dim_indices(site_size_r8, 1, this%column_index()) + call this%set_dim_indices(site_size_r8, 2, this%levscls_index()) + + call this%set_dim_indices(site_pft_r8, 1, this%column_index()) + call this%set_dim_indices(site_pft_r8, 2, this%levpft_index()) + + call this%set_dim_indices(site_age_r8, 1, this%column_index()) + call this%set_dim_indices(site_age_r8, 2, this%levage_index()) + end subroutine assemble_history_output_types ! =================================================================================== @@ -407,6 +461,48 @@ integer function levscpf_index(this) levscpf_index = this%levscpf_index_ end function levscpf_index + ! ======================================================================= + subroutine set_levscls_index(this, index) + implicit none + class(fates_history_interface_type), intent(inout) :: this + integer, intent(in) :: index + this%levscls_index_ = index + end subroutine set_levscls_index + + integer function levscls_index(this) + implicit none + class(fates_history_interface_type), intent(in) :: this + levscls_index = this%levscls_index_ + end function levscls_index + + ! ======================================================================= + subroutine set_levpft_index(this, index) + implicit none + class(fates_history_interface_type), intent(inout) :: this + integer, intent(in) :: index + this%levpft_index_ = index + end subroutine set_levpft_index + + integer function levpft_index(this) + implicit none + class(fates_history_interface_type), intent(in) :: this + levpft_index = this%levpft_index_ + end function levpft_index + + ! ======================================================================= + subroutine set_levage_index(this, index) + implicit none + class(fates_history_interface_type), intent(inout) :: this + integer, intent(in) :: index + this%levage_index_ = index + end subroutine set_levage_index + + integer function levage_index(this) + implicit none + class(fates_history_interface_type), intent(in) :: this + levage_index = this%levage_index_ + end function levage_index + ! ====================================================================================== subroutine flush_hvars(this,nc,upfreq_in) @@ -500,6 +596,7 @@ subroutine init_dim_kinds_maps(this) ! ---------------------------------------------------------------------------------- use FatesIOVariableKindMod, only : patch_r8, patch_ground_r8, patch_size_pft_r8 use FatesIOVariableKindMod, only : site_r8, site_ground_r8, site_size_pft_r8 + use FatesIOVariableKindMod, only : site_size_r8, site_pft_r8, site_age_r8 implicit none @@ -533,6 +630,18 @@ subroutine init_dim_kinds_maps(this) index = index + 1 call this%dim_kinds(index)%Init(site_size_pft_r8, 2) + ! site x size-class + index = index + 1 + call this%dim_kinds(index)%Init(site_size_r8, 2) + + ! site x pft + index = index + 1 + call this%dim_kinds(index)%Init(site_pft_r8, 2) + + ! site x patch-age clase + index = index + 1 + call this%dim_kinds(index)%Init(site_age_r8, 2) + ! FIXME(bja, 2016-10) assert(index == fates_history_num_dim_kinds) end subroutine init_dim_kinds_maps @@ -606,7 +715,10 @@ subroutine update_history_dyn(this,nc,nsites,sites) ed_patch_type, & AREA, & sclass_ed, & - nlevsclass_ed + nlevsclass_ed, & + levage_ed, & + nlevage_ed, & + levpft_ed use EDParamsMod , only : ED_val_ag_biomass ! Arguments @@ -618,7 +730,7 @@ subroutine update_history_dyn(this,nc,nsites,sites) ! Locals integer :: s ! The local site index integer :: io_si ! The site index of the IO array - integer :: ipa ! The local "I"ndex of "PA"tches + integer :: ipa, ipa2 ! The local "I"ndex of "PA"tches integer :: io_pa ! The patch index of the IO array integer :: io_pa1 ! The first patch index in the IO array for each site integer :: io_soipa @@ -637,6 +749,7 @@ subroutine update_history_dyn(this,nc,nsites,sites) real(r8), parameter :: daysecs = 86400.0_r8 ! What modeler doesn't recognize 86400? real(r8), parameter :: yeardays = 365.0_r8 ! ALM/CLM do not use leap-years + real(r8), parameter :: tiny = 1.e-5_r8 ! some small number associate( hio_npatches_si => this%hvars(ih_npatches_si)%r81d, & hio_ncohorts_si => this%hvars(ih_ncohorts_si)%r81d, & @@ -644,10 +757,10 @@ subroutine update_history_dyn(this,nc,nsites,sites) hio_area_plant_pa => this%hvars(ih_area_plant_pa)%r81d, & hio_area_treespread_pa => this%hvars(ih_area_treespread_pa)%r81d, & hio_canopy_spread_pa => this%hvars(ih_canopy_spread_pa)%r81d, & - hio_biomass_pa_pft => this%hvars(ih_biomass_pa_pft)%r82d, & - hio_leafbiomass_pa_pft => this%hvars(ih_leafbiomass_pa_pft)%r82d, & - hio_storebiomass_pa_pft => this%hvars(ih_storebiomass_pa_pft)%r82d, & - hio_nindivs_pa_pft => this%hvars(ih_nindivs_pa_pft)%r82d, & + hio_biomass_si_pft => this%hvars(ih_biomass_si_pft)%r82d, & + hio_leafbiomass_si_pft => this%hvars(ih_leafbiomass_si_pft)%r82d, & + hio_storebiomass_si_pft => this%hvars(ih_storebiomass_si_pft)%r82d, & + hio_nindivs_si_pft => this%hvars(ih_nindivs_si_pft)%r82d, & hio_nesterov_fire_danger_pa => this%hvars(ih_nesterov_fire_danger_pa)%r81d, & hio_spitfire_ros_pa => this%hvars(ih_spitfire_ROS_pa)%r81d, & hio_tfc_ros_pa => this%hvars(ih_TFC_ROS_pa)%r81d, & @@ -688,7 +801,14 @@ subroutine update_history_dyn(this,nc,nsites,sites) hio_m2_si_scpf => this%hvars(ih_m2_si_scpf)%r82d, & hio_m3_si_scpf => this%hvars(ih_m3_si_scpf)%r82d, & hio_m4_si_scpf => this%hvars(ih_m4_si_scpf)%r82d, & - hio_m5_si_scpf => this%hvars(ih_m5_si_scpf)%r82d ) + hio_m5_si_scpf => this%hvars(ih_m5_si_scpf)%r82d, & + hio_ba_si_scls => this%hvars(ih_ba_si_scls)%r82d, & + hio_area_si_age => this%hvars(ih_area_si_age)%r82d, & + hio_lai_si_age => this%hvars(ih_lai_si_age)%r82d, & + hio_canopy_area_si_age => this%hvars(ih_canopy_area_si_age)%r82d, & + hio_ncl_si_age => this%hvars(ih_ncl_si_age)%r82d, & + hio_npatches_si_age => this%hvars(ih_npatches_si_age)%r82d) + ! --------------------------------------------------------------------------------- ! Flush arrays to values defined by %flushval (see registry entry in @@ -720,6 +840,19 @@ subroutine update_history_dyn(this,nc,nsites,sites) ! Increment the number of patches per site hio_npatches_si(io_si) = hio_npatches_si(io_si) + 1._r8 + + ! Increment the fractional area in each age class bin + hio_area_si_age(io_si,cpatch%age_class) = hio_area_si_age(io_si,cpatch%age_class) & + + cpatch%area/AREA + + ! Increment some patch-age-resolved diagnostics + hio_lai_si_age(io_si,cpatch%age_class) = hio_lai_si_age(io_si,cpatch%age_class) & + + cpatch%lai * cpatch%area + hio_canopy_area_si_age(io_si,cpatch%age_class) = hio_canopy_area_si_age(io_si,cpatch%age_class) & + + cpatch%canopy_area/AREA + hio_ncl_si_age(io_si,cpatch%age_class) = hio_ncl_si_age(io_si,cpatch%age_class) & + + cpatch%ncl_p * cpatch%area + hio_npatches_si_age(io_si,cpatch%age_class) = hio_npatches_si_age(io_si,cpatch%age_class) + 1._r8 ccohort => cpatch%shortest do while(associated(ccohort)) @@ -768,17 +901,17 @@ subroutine update_history_dyn(this,nc,nsites,sites) hio_balive_pa(io_pa) = hio_balive_pa(io_pa) + n_density * ccohort%balive * 1.e3_r8 ! Update PFT partitioned biomass components - hio_biomass_pa_pft(io_pa,ft) = hio_biomass_pa_pft(io_pa,ft) + & - n_density * ccohort%b * 1.e3_r8 - - hio_leafbiomass_pa_pft(io_pa,ft) = hio_leafbiomass_pa_pft(io_pa,ft) + & - n_density * ccohort%bl * 1.e3_r8 + hio_leafbiomass_si_pft(io_si,ft) = hio_leafbiomass_si_pft(io_si,ft) + & + (ccohort%n / AREA) * ccohort%bl * 1.e3_r8 - hio_storebiomass_pa_pft(io_pa,ft) = hio_storebiomass_pa_pft(io_pa,ft) + & - n_density * ccohort%bstore * 1.e3_r8 + hio_storebiomass_si_pft(io_si,ft) = hio_storebiomass_si_pft(io_si,ft) + & + (ccohort%n / AREA) * ccohort%bstore * 1.e3_r8 - hio_nindivs_pa_pft(io_pa,ft) = hio_nindivs_pa_pft(io_pa,ft) + & - ccohort%n + hio_nindivs_si_pft(io_si,ft) = hio_nindivs_si_pft(io_si,ft) + & + ccohort%n / AREA + + hio_biomass_si_pft(io_si, ft) = hio_biomass_si_pft(io_si, ft) + & + (ccohort%n / AREA) * ccohort%b * 1.e3_r8 ! Site by Size-Class x PFT (SCPF) ! ------------------------------------------------------------------------ @@ -789,7 +922,8 @@ subroutine update_history_dyn(this,nc,nsites,sites) ! have any meaning, otherwise they are just inialization values if( .not.(ccohort%isnew) ) then - associate( scpf => ccohort%size_by_pft_class ) + associate( scpf => ccohort%size_by_pft_class, & + scls => ccohort%size_class ) hio_gpp_si_scpf(io_si,scpf) = hio_gpp_si_scpf(io_si,scpf) + & n_perm2*ccohort%gpp_acc_hold ! [kgC/m2/yr] @@ -840,6 +974,9 @@ subroutine update_history_dyn(this,nc,nsites,sites) ! basal area [m2/ha] hio_ba_si_scpf(io_si,scpf) = hio_ba_si_scpf(io_si,scpf) + & 0.25_r8*3.14159_r8*((dbh/100.0_r8)**2.0_r8)*n_perm2*AREA + ! also by size class only + hio_ba_si_scls(io_si,scls) = hio_ba_si_scls(io_si,scls) + & + 0.25_r8*3.14159_r8*((dbh/100.0_r8)**2.0_r8)*n_perm2*AREA ! number density [/ha] hio_nplant_si_scpf(io_si,scpf) = hio_nplant_si_scpf(io_si,scpf) + AREA*n_perm2 @@ -903,6 +1040,17 @@ subroutine update_history_dyn(this,nc,nsites,sites) ipa = ipa + 1 cpatch => cpatch%younger end do !patch loop + + ! divide so-far-just-summed but to-be-averaged patch-age-class variables by patch-age-class area to get mean values + do ipa2 = 1, nlevage_ed + if (hio_area_si_age(io_si, ipa2) .gt. tiny) then + hio_lai_si_age(io_si, ipa2) = hio_lai_si_age(io_si, ipa2) / (hio_area_si_age(io_si, ipa2)*AREA) + hio_ncl_si_age(io_si, ipa2) = hio_ncl_si_age(io_si, ipa2) / (hio_area_si_age(io_si, ipa2)*AREA) + else + hio_lai_si_age(io_si, ipa2) = 0._r8 + hio_ncl_si_age(io_si, ipa2) = 0._r8 + endif + end do enddo ! site loop @@ -924,6 +1072,7 @@ subroutine update_history_prod(this,nc,nsites,sites,dt_tstep) ed_cohort_type, & ed_patch_type, & AREA, & + nlevage_ed, & sclass_ed, & nlevsclass_ed ! Arguments @@ -945,6 +1094,9 @@ subroutine update_history_prod(this,nc,nsites,sites,dt_tstep) integer :: ft ! functional type index real(r8) :: n_density ! individual of cohort per m2. real(r8) :: n_perm2 ! individuals per m2 for the whole column + real(r8) :: patch_area_by_age(nlevage_ed) ! patch area in each bin for normalizing purposes + real(r8), parameter :: tiny = 1.e-5_r8 ! some small number + integer :: ipa2 ! patch incrementer type(fates_history_variable_type),pointer :: hvar type(ed_patch_type),pointer :: cpatch @@ -965,7 +1117,10 @@ subroutine update_history_prod(this,nc,nsites,sites,dt_tstep) hio_ar_agsapm_si_scpf => this%hvars(ih_ar_agsapm_si_scpf)%r82d, & hio_ar_darkm_si_scpf => this%hvars(ih_ar_darkm_si_scpf)%r82d, & hio_ar_crootm_si_scpf => this%hvars(ih_ar_crootm_si_scpf)%r82d, & - hio_ar_frootm_si_scpf => this%hvars(ih_ar_frootm_si_scpf)%r82d ) + hio_ar_frootm_si_scpf => this%hvars(ih_ar_frootm_si_scpf)%r82d, & + hio_gpp_si_age => this%hvars(ih_gpp_si_age)%r82d, & + hio_npp_si_age => this%hvars(ih_npp_si_age)%r82d & + ) ! Flush the relevant history variables @@ -979,10 +1134,15 @@ subroutine update_history_prod(this,nc,nsites,sites,dt_tstep) ipa = 0 cpatch => sites(s)%oldest_patch + + patch_area_by_age(:) = 0._r8 + do while(associated(cpatch)) io_pa = io_pa1 + ipa + patch_area_by_age(cpatch%age_class) = patch_area_by_age(cpatch%age_class) + cpatch%area + ccohort => cpatch%shortest do while(associated(ccohort)) @@ -1045,6 +1205,11 @@ subroutine update_history_prod(this,nc,nsites,sites,dt_tstep) hio_ar_frootm_si_scpf(io_si,scpf) = hio_ar_frootm_si_scpf(io_si,scpf) + & ccohort%froot_mr * n_perm2 * daysecs * yeardays + ! accumulate fluxes per patch age bin + hio_gpp_si_age(io_si,cpatch%age_class) = hio_gpp_si_age(io_si,cpatch%age_class) & + + ccohort%gpp_tstep * ccohort%n * 1.e3_r8 / dt_tstep + hio_npp_si_age(io_si,cpatch%age_class) = hio_npp_si_age(io_si,cpatch%age_class) & + + ccohort%npp_tstep * ccohort%n * 1.e3_r8 / dt_tstep end associate endif @@ -1053,6 +1218,16 @@ subroutine update_history_prod(this,nc,nsites,sites,dt_tstep) ipa = ipa + 1 cpatch => cpatch%younger end do !patch loop + + do ipa2 = 1, nlevage_ed + if (patch_area_by_age(ipa2) .gt. tiny) then + hio_gpp_si_age(io_si, ipa2) = hio_gpp_si_age(io_si, ipa2) / (patch_area_by_age(ipa2)) + hio_npp_si_age(io_si, ipa2) = hio_npp_si_age(io_si, ipa2) / (patch_area_by_age(ipa2)) + else + hio_gpp_si_age(io_si, ipa2) = 0._r8 + hio_npp_si_age(io_si, ipa2) = 0._r8 + endif + end do enddo ! site loop @@ -1124,6 +1299,7 @@ subroutine define_history_vars(this, initialize_variables) use FatesIOVariableKindMod, only : patch_r8, patch_ground_r8, patch_size_pft_r8 use FatesIOVariableKindMod, only : site_r8, site_ground_r8, site_size_pft_r8 + use FatesIOVariableKindMod, only : site_size_r8, site_pft_r8, site_age_r8 implicit none class(fates_history_interface_type), intent(inout) :: this @@ -1169,23 +1345,49 @@ subroutine define_history_vars(this, initialize_variables) call this%set_history_var(vname='PFTbiomass', units='gC/m2', & long='total PFT level biomass', use_default='active', & - avgflag='A', vtype=patch_ground_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & - ivar=ivar, initialize=initialize_variables, index = ih_biomass_pa_pft ) + avgflag='A', vtype=site_pft_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & + ivar=ivar, initialize=initialize_variables, index = ih_biomass_si_pft ) call this%set_history_var(vname='PFTleafbiomass', units='gC/m2', & long='total PFT level leaf biomass', use_default='active', & - avgflag='A', vtype=patch_ground_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & - ivar=ivar, initialize=initialize_variables, index = ih_leafbiomass_pa_pft ) + avgflag='A', vtype=site_pft_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & + ivar=ivar, initialize=initialize_variables, index = ih_leafbiomass_si_pft ) call this%set_history_var(vname='PFTstorebiomass', units='gC/m2', & long='total PFT level stored biomass', use_default='active', & - avgflag='A', vtype=patch_ground_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & - ivar=ivar, initialize=initialize_variables, index = ih_storebiomass_pa_pft ) + avgflag='A', vtype=site_pft_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & + ivar=ivar, initialize=initialize_variables, index = ih_storebiomass_si_pft ) call this%set_history_var(vname='PFTnindivs', units='indiv / m2', & long='total PFT level number of individuals', use_default='active', & - avgflag='A', vtype=patch_ground_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & - ivar=ivar, initialize=initialize_variables, index = ih_nindivs_pa_pft ) + avgflag='A', vtype=site_pft_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & + ivar=ivar, initialize=initialize_variables, index = ih_nindivs_si_pft ) + + ! patch age class variables + call this%set_history_var(vname='PATCH_AREA_BY_AGE', units='m2/m2', & + long='patch area by age bin', use_default='active', & + avgflag='A', vtype=site_age_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & + ivar=ivar, initialize=initialize_variables, index = ih_area_si_age ) + + call this%set_history_var(vname='LAI_BY_AGE', units='m2/m2', & + long='leaf area index by age bin', use_default='active', & + avgflag='A', vtype=site_age_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & + ivar=ivar, initialize=initialize_variables, index = ih_lai_si_age ) + + call this%set_history_var(vname='CANOPY_AREA_BY_AGE', units='m2/m2', & + long='canopy area by age bin', use_default='active', & + avgflag='A', vtype=site_age_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & + ivar=ivar, initialize=initialize_variables, index = ih_canopy_area_si_age ) + + call this%set_history_var(vname='NCL_BY_AGE', units='--', & + long='number of canopy levels by age bin', use_default='inactive', & + avgflag='A', vtype=site_age_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & + ivar=ivar, initialize=initialize_variables, index = ih_ncl_si_age ) + + call this%set_history_var(vname='NPATCH_BY_AGE', units='--', & + long='number of patches by age bin', use_default='active', & + avgflag='A', vtype=site_age_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & + ivar=ivar, initialize=initialize_variables, index = ih_npatches_si_age ) ! Fire Variables @@ -1341,6 +1543,18 @@ subroutine define_history_vars(this, initialize_variables) avgflag='A', vtype=patch_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=2, & ivar=ivar, initialize=initialize_variables, index = ih_maint_resp_pa ) + ! fast fluxes by age bin + call this%set_history_var(vname='NPP_BY_AGE', units='gC/m^2/s', & + long='net primary productivity by age bin', use_default='active', & + avgflag='A', vtype=site_age_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=2, & + ivar=ivar, initialize=initialize_variables, index = ih_npp_si_age ) + + call this%set_history_var(vname='GPP_BY_AGE', units='gC/m^2/s', & + long='gross primary productivity by age bin', use_default='active', & + avgflag='A', vtype=site_age_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=2, & + ivar=ivar, initialize=initialize_variables, index = ih_gpp_si_age ) + + ! Carbon Flux (grid dimension x scpf) (THESE ARE DEFAULT INACTIVE!!! ! (BECAUSE THEY TAKE UP SPACE!!! @@ -1474,6 +1688,11 @@ subroutine define_history_vars(this, initialize_variables) avgflag='A', vtype=site_size_pft_r8, hlms='CLM:ALM', flushval=0.0_r8, & upfreq=2, ivar=ivar, initialize=initialize_variables, index = ih_ar_frootm_si_scpf ) + ! size-class only variables + call this%set_history_var(vname='BA_SCLS', units = 'm2/ha', & + long='basal area by size class', use_default='active', & + avgflag='A', vtype=site_size_r8, hlms='CLM:ALM', flushval=0.0_r8, & + upfreq=1, ivar=ivar, initialize=initialize_variables, index = ih_ba_si_scls ) ! CARBON BALANCE VARIABLES THAT DEPEND ON HLM BGC INPUTS diff --git a/main/FatesHistoryVariableType.F90 b/main/FatesHistoryVariableType.F90 index 44d6458668..20abd41f89 100644 --- a/main/FatesHistoryVariableType.F90 +++ b/main/FatesHistoryVariableType.F90 @@ -45,6 +45,7 @@ subroutine Init(this, vname, units, long, use_default, & use FatesIOVariableKindMod, only : patch_r8, patch_ground_r8, patch_size_pft_r8 use FatesIOVariableKindMod, only : site_r8, site_ground_r8, site_size_pft_r8 + use FatesIOVariableKindMod, only : site_size_r8, site_pft_r8, site_age_r8 use FatesIOVariableKindMod, only : iotype_index implicit none @@ -118,6 +119,18 @@ subroutine Init(this, vname, units, long, use_default, & allocate(this%r82d(lb1:ub1, lb2:ub2)) this%r82d(:,:) = flushval + case(site_size_r8) + allocate(this%r82d(lb1:ub1, lb2:ub2)) + this%r82d(:,:) = flushval + + case(site_pft_r8) + allocate(this%r82d(lb1:ub1, lb2:ub2)) + this%r82d(:,:) = flushval + + case(site_age_r8) + allocate(this%r82d(lb1:ub1, lb2:ub2)) + this%r82d(:,:) = flushval + case default write(fates_log(),*) 'Incompatible vtype passed to set_history_var' write(fates_log(),*) 'vtype = ',trim(vtype),' ?' @@ -183,6 +196,7 @@ subroutine Flush(this, thread, dim_bounds, dim_kinds) use FatesIODimensionsMod, only : fates_io_dimension_type use FatesIOVariableKindMod, only : patch_r8, patch_ground_r8, patch_size_pft_r8 use FatesIOVariableKindMod, only : site_r8, site_ground_r8, site_size_pft_r8, patch_int + use FatesIOVariableKindMod, only : site_size_r8, site_pft_r8, site_age_r8 implicit none @@ -208,6 +222,12 @@ subroutine Flush(this, thread, dim_bounds, dim_kinds) this%r82d(lb1:ub1, lb2:ub2) = this%flushval case(site_size_pft_r8) this%r82d(lb1:ub1, lb2:ub2) = this%flushval + case(site_size_r8) + this%r82d(lb1:ub1, lb2:ub2) = this%flushval + case(site_pft_r8) + this%r82d(lb1:ub1, lb2:ub2) = this%flushval + case(site_age_r8) + this%r82d(lb1:ub1, lb2:ub2) = this%flushval case(patch_int) this%int1d(lb1:ub1) = nint(this%flushval) case default diff --git a/main/FatesIODimensionsMod.F90 b/main/FatesIODimensionsMod.F90 index 84c082e75c..83b2475aad 100644 --- a/main/FatesIODimensionsMod.F90 +++ b/main/FatesIODimensionsMod.F90 @@ -9,6 +9,9 @@ module FatesIODimensionsMod character(*), parameter :: column = 'column' character(*), parameter :: levgrnd = 'levgrnd' character(*), parameter :: levscpf = 'levscpf' + character(*), parameter :: levscls = 'levscls' + character(*), parameter :: levpft = 'levpft' + character(*), parameter :: levage = 'levage' ! patch = This is a structure that records where FATES patch boundaries ! on each thread point to in the host IO array, this structure @@ -24,6 +27,15 @@ module FatesIODimensionsMod ! levscpf = This is a structure that records the boundaries for the ! number of size-class x pft dimension + ! levscls = This is a structure that records the boundaries for the + ! number of size-class dimension + + ! levpft = This is a structure that records the boundaries for the + ! number of pft dimension + + ! levage = This is a structure that records the boundaries for the + ! number of patch-age-class dimension + type, public :: fates_bounds_type integer :: patch_begin @@ -34,8 +46,14 @@ module FatesIODimensionsMod integer :: column_end ! we call this a "site" (rgk 11-2016) integer :: ground_begin integer :: ground_end + integer :: sizepft_class_begin + integer :: sizepft_class_end + integer :: size_class_begin + integer :: size_class_end integer :: pft_class_begin integer :: pft_class_end + integer :: age_class_begin + integer :: age_class_end end type fates_bounds_type diff --git a/main/FatesIOVariableKindMod.F90 b/main/FatesIOVariableKindMod.F90 index 343d3b4364..2c8eb98216 100644 --- a/main/FatesIOVariableKindMod.F90 +++ b/main/FatesIOVariableKindMod.F90 @@ -16,9 +16,13 @@ module FatesIOVariableKindMod character(*), parameter :: site_int = 'SI_INT' character(*), parameter :: site_ground_r8 = 'SI_GRND_R8' character(*), parameter :: site_size_pft_r8 = 'SI_SCPF_R8' + character(*), parameter :: site_size_r8 = 'SI_SCLS_R8' character(*), parameter :: patch_int = 'PA_INT' character(*), parameter :: cohort_r8 = 'CO_R8' character(*), parameter :: cohort_int = 'CO_INT' + character(*), parameter :: site_pft_r8 = 'SI_PFT_R8' + character(*), parameter :: site_age_r8 = 'SI_AGE_R8' + ! NOTE(RGK, 2016) %active is not used yet. Was intended as a check on the HLM->FATES ! control parameter passing to ensure all active dimension types received all