diff --git a/biogeochem/EDCanopyStructureMod.F90 b/biogeochem/EDCanopyStructureMod.F90 index 5ce6d6631f..f5419cedd7 100755 --- a/biogeochem/EDCanopyStructureMod.F90 +++ b/biogeochem/EDCanopyStructureMod.F90 @@ -742,7 +742,7 @@ subroutine canopy_summarization( nsites, sites, bc_in ) if(currentCohort%balive <= 0._r8)then write(fates_log(),*) 'ED: balive is zero in canopy_summarization',currentCohort%balive endif - + currentCohort => currentCohort%taller enddo ! ends 'do while(associated(currentCohort)) @@ -995,8 +995,6 @@ subroutine leaf_area_profile( currentSite , snow_depth_si, frac_sno_eff_si) layer_bottom_hite = currentCohort%hite-(((iv+1)/currentCohort%NV) * currentCohort%hite * & EDecophyscon%crown(currentCohort%pft)) ! pftcon%vertical_canopy_frac(ft)) - write(fates_log(), *) 'calc snow 2', snow_depth_si , frac_sno_eff_si - fraction_exposed =1.0_r8 currentPatch%tlai_profile(L,ft,iv) = currentPatch%tlai_profile(L,ft,iv)+ dinc_ed * fleaf * & @@ -1014,9 +1012,12 @@ subroutine leaf_area_profile( currentSite , snow_depth_si, frac_sno_eff_si) currentPatch%layer_height_profile(L,ft,iv) = currentPatch%layer_height_profile(L,ft,iv) + (dinc_ed * fleaf * & currentCohort%c_area/currentPatch%total_canopy_area *(layer_top_hite+layer_bottom_hite)/2.0_r8) !average height of layer. - write(fates_log(), *) 'LHP', currentPatch%layer_height_profile(L,ft,iv) - if ( DEBUG ) write(fates_log(), *) 'EDCLMLink 1246 ', currentPatch%elai_profile(1,ft,iv) - + if ( DEBUG ) then + write(fates_log(), *) 'calc snow 2', snow_depth_si , frac_sno_eff_si + write(fates_log(), *) 'LHP', currentPatch%layer_height_profile(L,ft,iv) + write(fates_log(), *) 'EDCLMLink 1246 ', currentPatch%elai_profile(1,ft,iv) + end if + end do !Bottom layer diff --git a/biogeochem/EDCohortDynamicsMod.F90 b/biogeochem/EDCohortDynamicsMod.F90 index 5fc649bb0d..2ccea2cad1 100755 --- a/biogeochem/EDCohortDynamicsMod.F90 +++ b/biogeochem/EDCohortDynamicsMod.F90 @@ -4,8 +4,10 @@ module EDCohortDynamicsMod ! Cohort stuctures in ED. ! ! !USES: - use shr_kind_mod , only : r8 => shr_kind_r8; - use clm_varctl , only : iulog + use abortutils , only : endrun + use FatesGlobals , only : fates_log + use FatesConstantsMod , only : r8 => fates_r8 + use shr_log_mod , only : errMsg => shr_log_errMsg use pftconMod , only : pftcon use EDEcophysContype , only : EDecophyscon use EDGrowthFunctionsMod , only : c_area, tree_lai @@ -32,6 +34,9 @@ module EDCohortDynamicsMod logical, parameter :: DEBUG = .false. ! local debug flag + character(len=*), parameter, private :: sourcefile = & + __FILE__ + ! 10/30/09: Created by Rosie Fisher !-------------------------------------------------------------------------------------! @@ -95,13 +100,18 @@ subroutine create_cohort(patchptr, pft, nn, hite, dbh, & call size_and_type_class_index(new_cohort%dbh,new_cohort%pft, & new_cohort%size_class,new_cohort%size_by_pft_class) - if ( DEBUG ) write(iulog,*) 'EDCohortDyn I ',bstore + if ( DEBUG ) write(fates_log(),*) 'EDCohortDyn I ',bstore - if (new_cohort%dbh <= 0.0_r8 .or. new_cohort%n == 0._r8 .or. new_cohort%pft == 0 & - .or. new_cohort%canopy_trim <= 0.0_r8 .or. new_cohort%balive <= 0._r8) then - write(iulog,*) 'ED: something is zero in create_cohort', & + ! This routine may be called during restarts, and at this point in the call sequence + ! the actual cohort data is unknown, as this is really only used for allocation + ! In these cases, testing if things like biomass are reasonable is pre-mature + ! However, in this part of the code, we will pass in nominal values for size, number and type + + if (new_cohort%dbh <= 0.0_r8 .or. new_cohort%n == 0._r8 .or. new_cohort%pft == 0 ) then + write(fates_log(),*) 'ED: something is zero in create_cohort', & new_cohort%indexnumber,new_cohort%dbh,new_cohort%n, & - new_cohort%pft,new_cohort%canopy_trim,new_cohort%balive + new_cohort%pft + call endrun(msg=errMsg(sourcefile, __LINE__)) endif if (new_cohort%siteptr%status==2 .and. pftcon%season_decid(pft) == 1) then @@ -274,12 +284,12 @@ subroutine allocate_live_biomass(cc_p,mode) endif if (abs(currentcohort%balive -currentcohort%bl- currentcohort%br - currentcohort%bsw)>1e-12) then - write(iulog,*) 'issue with carbon allocation in create_cohort,allocate_live_biomass',& + write(fates_log(),*) 'issue with carbon allocation in create_cohort,allocate_live_biomass',& currentcohort%balive -currentcohort%bl- currentcohort%br - currentcohort%bsw, & currentcohort%status_coh,currentcohort%balive - write(iulog,*) 'actual vs predicted balive',ideal_balive,currentcohort%balive ,ratio_balive,leaf_frac - write(iulog,*) 'leaf,root,stem',currentcohort%bl,currentcohort%br,currentcohort%bsw - write(iulog,*) 'pft',ft,pftcon%evergreen(ft),pftcon%season_decid(ft),leaves_off_switch + write(fates_log(),*) 'actual vs predicted balive',ideal_balive,currentcohort%balive ,ratio_balive,leaf_frac + write(fates_log(),*) 'leaf,root,stem',currentcohort%bl,currentcohort%br,currentcohort%bsw + write(fates_log(),*) 'pft',ft,pftcon%evergreen(ft),pftcon%season_decid(ft),leaves_off_switch endif currentCohort%b = currentCohort%bdead + currentCohort%balive @@ -347,16 +357,16 @@ subroutine nan_cohort(cc_p) currentCohort%treesai = nan ! stem area index of tree (total stem area (m2) / canopy area (m2) ! CARBON FLUXES - currentCohort%gpp = nan ! GPP: kgC/indiv/year - currentCohort%gpp_tstep = nan ! GPP: kgC/indiv/timestep + currentCohort%gpp_acc_hold = nan ! GPP: kgC/indiv/year + currentCohort%gpp_tstep = nan ! GPP: kgC/indiv/timestep currentCohort%gpp_acc = nan ! GPP: kgC/indiv/day - currentCohort%npp = nan ! NPP: kgC/indiv/year - currentCohort%npp_tstep = nan ! NPP: kGC/indiv/timestep + currentCohort%npp_acc_hold = nan ! NPP: kgC/indiv/year + currentCohort%npp_tstep = nan ! NPP: kGC/indiv/timestep currentCohort%npp_acc = nan ! NPP: kgC/indiv/day currentCohort%year_net_uptake(:) = nan ! Net uptake of individual leaf layers kgC/m2/year currentCohort%ts_net_uptake(:) = nan ! Net uptake of individual leaf layers kgC/m2/s - currentCohort%resp = nan ! RESP: kgC/indiv/year - currentCohort%resp_tstep = nan ! RESP: kgC/indiv/timestep + currentCohort%resp_acc_hold = nan ! RESP: kgC/indiv/year + currentCohort%resp_tstep = nan ! RESP: kgC/indiv/timestep currentCohort%resp_acc = nan ! RESP: kGC/cohort/day currentCohort%npp_leaf = nan @@ -435,10 +445,10 @@ subroutine zero_cohort(cc_p) currentcohort%npp_acc = 0._r8 currentcohort%gpp_acc = 0._r8 currentcohort%resp_acc = 0._r8 - currentcohort%npp_tstep = 0._r8 - currentcohort%gpp_tstep = 0._r8 - currentcohort%resp_tstep = 0._r8 - currentcohort%resp = 0._r8 + currentcohort%npp_tstep = 0._r8 + currentcohort%gpp_tstep = 0._r8 + currentcohort%resp_tstep = 0._r8 + currentcohort%resp_acc_hold = 0._r8 currentcohort%carbon_balance = 0._r8 currentcohort%leaf_litter = 0._r8 currentcohort%year_net_uptake(:) = 999 ! this needs to be 999, or trimming of new cohorts will break. @@ -448,8 +458,8 @@ subroutine zero_cohort(cc_p) currentcohort%md = 0._r8 currentcohort%root_md = 0._r8 currentcohort%leaf_md = 0._r8 - currentcohort%npp = 0._r8 - currentcohort%gpp = 0._r8 + currentcohort%npp_acc_hold = 0._r8 + currentcohort%gpp_acc_hold = 0._r8 currentcohort%storage_flux = 0._r8 currentcohort%dmort = 0._r8 currentcohort%gscan = 0._r8 @@ -496,7 +506,7 @@ subroutine terminate_cohorts( patchptr ) if (currentcohort%n < min_n_safemath) then terminate = 1 if ( DEBUG ) then - write(iulog,*) 'terminating cohorts 0',currentCohort%n/currentPatch%area,currentCohort%dbh + write(fates_log(),*) 'terminating cohorts 0',currentCohort%n/currentPatch%area,currentCohort%dbh endif endif @@ -510,7 +520,7 @@ subroutine terminate_cohorts( patchptr ) terminate = 1 if ( DEBUG ) then - write(iulog,*) 'terminating cohorts 1',currentCohort%n/currentPatch%area,currentCohort%dbh + write(fates_log(),*) 'terminating cohorts 1',currentCohort%n/currentPatch%area,currentCohort%dbh endif endif @@ -518,7 +528,7 @@ subroutine terminate_cohorts( patchptr ) if (currentCohort%canopy_layer > cp_nclmax ) then terminate = 1 if ( DEBUG ) then - write(iulog,*) 'terminating cohorts 2', currentCohort%canopy_layer + write(fates_log(),*) 'terminating cohorts 2', currentCohort%canopy_layer endif endif @@ -526,7 +536,7 @@ subroutine terminate_cohorts( patchptr ) if (currentCohort%balive < 1e-10_r8 .or. currentCohort%bstore < 1e-10_r8) then terminate = 1 if ( DEBUG ) then - write(iulog,*) 'terminating cohorts 3', currentCohort%balive,currentCohort%bstore + write(fates_log(),*) 'terminating cohorts 3', currentCohort%balive,currentCohort%bstore endif endif @@ -534,7 +544,7 @@ subroutine terminate_cohorts( patchptr ) if (currentCohort%balive+currentCohort%bdead+currentCohort%bstore < 0._r8) then terminate = 1 if ( DEBUG ) then - write(iulog,*) 'terminating cohorts 4', currentCohort%balive, & + write(fates_log(),*) 'terminating cohorts 4', currentCohort%balive, & currentCohort%bstore, currentCohort%bdead, & currentCohort%balive+currentCohort%bdead+& currentCohort%bstore, currentCohort%n @@ -667,11 +677,11 @@ subroutine fuse_cohorts(patchptr) currentCohort%balive = (currentCohort%n*currentCohort%balive + nextc%n*nextc%balive)/newn currentCohort%bdead = (currentCohort%n*currentCohort%bdead + nextc%n*nextc%bdead)/newn - if ( DEBUG ) write(iulog,*) 'EDcohortDyn I ',currentCohort%bstore + if ( DEBUG ) write(fates_log(),*) 'EDcohortDyn I ',currentCohort%bstore currentCohort%bstore = (currentCohort%n*currentCohort%bstore + nextc%n*nextc%bstore)/newn - if ( DEBUG ) write(iulog,*) 'EDcohortDyn II ',currentCohort%bstore + if ( DEBUG ) write(fates_log(),*) 'EDcohortDyn II ',currentCohort%bstore currentCohort%seed_prod = (currentCohort%n*currentCohort%seed_prod + nextc%n*nextc%seed_prod)/newn currentCohort%root_md = (currentCohort%n*currentCohort%root_md + nextc%n*nextc%root_md)/newn @@ -689,10 +699,10 @@ subroutine fuse_cohorts(patchptr) currentCohort%bsw = (currentCohort%n*currentCohort%bsw + nextc%n*nextc%bsw)/newn currentCohort%bl = (currentCohort%n*currentCohort%bl + nextc%n*nextc%bl)/newn - if ( DEBUG ) write(iulog,*) 'EDcohortDyn 569 ',currentCohort%br - if ( DEBUG ) write(iulog,*) 'EDcohortDyn 570 ',currentCohort%n - if ( DEBUG ) write(iulog,*) 'EDcohortDyn 571 ',nextc%br - if ( DEBUG ) write(iulog,*) 'EDcohortDyn 572 ',nextc%n + if ( DEBUG ) write(fates_log(),*) 'EDcohortDyn 569 ',currentCohort%br + if ( DEBUG ) write(fates_log(),*) 'EDcohortDyn 570 ',currentCohort%n + if ( DEBUG ) write(fates_log(),*) 'EDcohortDyn 571 ',nextc%br + if ( DEBUG ) write(fates_log(),*) 'EDcohortDyn 572 ',nextc%n currentCohort%br = (currentCohort%n*currentCohort%br + nextc%n*nextc%br)/newn currentCohort%hite = (currentCohort%n*currentCohort%hite + nextc%n*nextc%hite)/newn @@ -700,18 +710,25 @@ subroutine fuse_cohorts(patchptr) currentCohort%gpp_acc = (currentCohort%n*currentCohort%gpp_acc + nextc%n*nextc%gpp_acc)/newn - if ( DEBUG ) write(iulog,*) 'EDcohortDyn III ',currentCohort%npp_acc - if ( DEBUG ) write(iulog,*) 'EDcohortDyn IV ',currentCohort%resp_acc + if ( DEBUG ) write(fates_log(),*) 'EDcohortDyn III ',currentCohort%npp_acc + if ( DEBUG ) write(fates_log(),*) 'EDcohortDyn IV ',currentCohort%resp_acc currentCohort%npp_acc = (currentCohort%n*currentCohort%npp_acc + nextc%n*nextc%npp_acc)/newn currentCohort%resp_acc = (currentCohort%n*currentCohort%resp_acc + nextc%n*nextc%resp_acc)/newn - if ( DEBUG ) write(iulog,*) 'EDcohortDyn V ',currentCohort%npp_acc - if ( DEBUG ) write(iulog,*) 'EDcohortDyn VI ',currentCohort%resp_acc + if ( DEBUG ) write(fates_log(),*) 'EDcohortDyn V ',currentCohort%npp_acc + if ( DEBUG ) write(fates_log(),*) 'EDcohortDyn VI ',currentCohort%resp_acc + + currentCohort%resp_acc_hold = & + (currentCohort%n*currentCohort%resp_acc_hold + & + nextc%n*nextc%resp_acc_hold)/newn + currentCohort%npp_acc_hold = & + (currentCohort%n*currentCohort%npp_acc_hold + & + nextc%n*nextc%npp_acc_hold)/newn + currentCohort%gpp_acc_hold = & + (currentCohort%n*currentCohort%gpp_acc_hold + & + nextc%n*nextc%gpp_acc_hold)/newn - currentCohort%resp = (currentCohort%n*currentCohort%resp + nextc%n*nextc%resp)/newn - currentCohort%npp = (currentCohort%n*currentCohort%npp + nextc%n*nextc%npp)/newn - currentCohort%gpp = (currentCohort%n*currentCohort%gpp + nextc%n*nextc%gpp)/newn currentCohort%canopy_trim = (currentCohort%n*currentCohort%canopy_trim + nextc%n*nextc%canopy_trim)/newn currentCohort%dmort = (currentCohort%n*currentCohort%dmort + nextc%n*nextc%dmort)/newn currentCohort%fire_mort = (currentCohort%n*currentCohort%fire_mort + nextc%n*nextc%fire_mort)/newn @@ -791,7 +808,7 @@ subroutine fuse_cohorts(patchptr) !---------------------------------------------------------------------! dynamic_fusion_tolerance = dynamic_fusion_tolerance * 1.1_r8 - write(iulog,*) 'maxcohorts exceeded',dynamic_fusion_tolerance + write(fates_log(),*) 'maxcohorts exceeded',dynamic_fusion_tolerance else iterate = 0 @@ -1017,20 +1034,20 @@ subroutine copy_cohort( currentCohort,copyc ) n%excl_weight = o%excl_weight n%prom_weight = o%prom_weight - ! CARBON FLUXES - n%gpp = o%gpp + ! CARBON FLUXES + n%gpp_acc_hold = o%gpp_acc_hold n%gpp_acc = o%gpp_acc - n%gpp_tstep = o%gpp_tstep - n%npp = o%npp - n%npp_tstep = o%npp_tstep + n%gpp_tstep = o%gpp_tstep + n%npp_acc_hold = o%npp_acc_hold + n%npp_tstep = o%npp_tstep - if ( DEBUG ) write(iulog,*) 'EDcohortDyn Ia ',o%npp_acc - if ( DEBUG ) write(iulog,*) 'EDcohortDyn Ib ',o%resp_acc + if ( DEBUG ) write(fates_log(),*) 'EDcohortDyn Ia ',o%npp_acc + if ( DEBUG ) write(fates_log(),*) 'EDcohortDyn Ib ',o%resp_acc - n%npp_acc = o%npp_acc - n%resp_tstep = o%resp_tstep + n%npp_acc_hold = o%npp_acc_hold + n%resp_tstep = o%resp_tstep n%resp_acc = o%resp_acc - n%resp = o%resp + n%resp_acc_hold = o%resp_acc_hold n%year_net_uptake = o%year_net_uptake n%ts_net_uptake = o%ts_net_uptake @@ -1080,7 +1097,7 @@ subroutine copy_cohort( currentCohort,copyc ) n%dbdeaddt = o%dbdeaddt n%dbstoredt = o%dbstoredt - if ( DEBUG ) write(iulog,*) 'EDCohortDyn dpstoredt ',o%dbstoredt + if ( DEBUG ) write(fates_log(),*) 'EDCohortDyn dpstoredt ',o%dbstoredt n%storage_flux = o%storage_flux @@ -1129,7 +1146,7 @@ function count_cohorts( currentPatch ) result ( backcount ) enddo if (backcount /= currentPatch%countcohorts) then - write(iulog,*) 'problem with linked list, not symmetrical' + write(fates_log(),*) 'problem with linked list, not symmetrical' endif end function count_cohorts diff --git a/biogeochem/EDPhysiologyMod.F90 b/biogeochem/EDPhysiologyMod.F90 index c0d2b25dfb..fccd8c0843 100755 --- a/biogeochem/EDPhysiologyMod.F90 +++ b/biogeochem/EDPhysiologyMod.F90 @@ -789,9 +789,9 @@ subroutine Growth_Derivatives( currentSite, currentCohort) ! NPP if ( DEBUG ) write(iulog,*) 'EDphys 716 ',currentCohort%npp_acc - currentCohort%npp = currentCohort%npp_acc * udata%n_sub !Link to CLM. convert from kgC/indiv/day into kgC/indiv/year - currentCohort%gpp = currentCohort%gpp_acc * udata%n_sub !Link to CLM. convert from kgC/indiv/day into kgC/indiv/year - currentCohort%resp = currentCohort%resp_acc * udata%n_sub !Link to CLM. convert from kgC/indiv/day into kgC/indiv/year + currentCohort%npp_acc_hold = currentCohort%npp_acc * udata%n_sub ! convert from kgC/indiv/day into kgC/indiv/year + currentCohort%gpp_acc_hold = currentCohort%gpp_acc * udata%n_sub ! convert from kgC/indiv/day into kgC/indiv/year + currentCohort%resp_acc_hold = currentCohort%resp_acc * udata%n_sub ! convert from kgC/indiv/day into kgC/indiv/year currentSite%flux_in = currentSite%flux_in + currentCohort%npp_acc * currentCohort%n @@ -833,16 +833,17 @@ subroutine Growth_Derivatives( currentSite, currentCohort) ! Calculate carbon balance ! this is the fraction of maintenance demand we -have- to do... - if ( DEBUG ) write(iulog,*) 'EDphys 760 ',currentCohort%npp, currentCohort%md, & + if ( DEBUG ) write(iulog,*) 'EDphys 760 ',currentCohort%npp_acc_hold, currentCohort%md, & EDecophyscon%leaf_stor_priority(currentCohort%pft) - currentCohort%carbon_balance = currentCohort%npp - currentCohort%md * EDecophyscon%leaf_stor_priority(currentCohort%pft) + currentCohort%carbon_balance = currentCohort%npp_acc_hold - & + currentCohort%md * EDecophyscon%leaf_stor_priority(currentCohort%pft) ! Allowing only carbon from NPP pool to account for npp flux into the maintenance pools ! ie this does not include any use of storage carbon or balive to make up for missing carbon balance in the transfer - currentCohort%npp_leaf = min(currentCohort%npp*currentCohort%leaf_md/currentCohort%md, & + currentCohort%npp_leaf = min(currentCohort%npp_acc_hold*currentCohort%leaf_md/currentCohort%md, & currentCohort%leaf_md*EDecophyscon%leaf_stor_priority(currentCohort%pft)) - currentCohort%npp_froot = min(currentCohort%npp*currentCohort%root_md/currentCohort%md, & + currentCohort%npp_froot = min(currentCohort%npp_acc_hold*currentCohort%root_md/currentCohort%md, & currentCohort%root_md*EDecophyscon%leaf_stor_priority(currentCohort%pft)) @@ -944,12 +945,12 @@ subroutine Growth_Derivatives( currentSite, currentCohort) if ( DEBUG ) write(iulog,*) 'EDPhys dbstoredt I ',currentCohort%dbstoredt currentCohort%seed_prod = (1.0_r8 - gr_fract) * currentCohort%carbon_balance - if (abs(currentCohort%npp-(currentCohort%dbalivedt+currentCohort%dbdeaddt+currentCohort%dbstoredt+ & + if (abs(currentCohort%npp_acc_hold-(currentCohort%dbalivedt+currentCohort%dbdeaddt+currentCohort%dbstoredt+ & currentCohort%seed_prod+currentCohort%md)) > 0.0000000001_r8)then - write(iulog,*) 'error in carbon check growth derivs',currentCohort%npp- & + write(iulog,*) 'error in carbon check growth derivs',currentCohort%npp_acc_hold- & (currentCohort%dbalivedt+currentCohort%dbdeaddt+currentCohort%dbstoredt+currentCohort%seed_prod+currentCohort%md) write(iulog,*) 'cohort fluxes',currentCohort%pft,currentCohort%canopy_layer,currentCohort%n, & - currentCohort%npp,currentCohort%dbalivedt,balive_loss, & + currentCohort%npp_acc_hold,currentCohort%dbalivedt,balive_loss, & currentCohort%dbdeaddt,currentCohort%dbstoredt,currentCohort%seed_prod,currentCohort%md * & EDecophyscon%leaf_stor_priority(currentCohort%pft) write(iulog,*) 'proxies' ,target_balive,currentCohort%balive,currentCohort%dbh,va,vs,gr_fract diff --git a/biogeophys/EDSurfaceAlbedoMod.F90 b/biogeophys/EDSurfaceAlbedoMod.F90 index 5180c0f8b6..d76695916c 100644 --- a/biogeophys/EDSurfaceAlbedoMod.F90 +++ b/biogeophys/EDSurfaceAlbedoMod.F90 @@ -401,7 +401,8 @@ subroutine ED_Norman_Radiation (nsites, sites, bc_in, bc_out ) !How much diffuse light is intercepted and then reflected? refl_dif(L,ft,iv,ib) = (1._r8 - tr_dif_z(L,ft,iv)) * rhol(ft,ib) !How much diffuse light in this layer is transmitted? - tran_dif(L,ft,iv,ib) = (1._r8 - tr_dif_z(L,ft,iv)) * taul(ft,ib) + tr_dif_z(L,ft,iv) + tran_dif(L,ft,iv,ib) = (1._r8 - tr_dif_z(L,ft,iv)) * & + taul(ft,ib) + tr_dif_z(L,ft,iv) end do !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++! @@ -418,13 +419,17 @@ subroutine ED_Norman_Radiation (nsites, sites, bc_in, bc_out ) ! FIX(RF,032414) ray tracing eqution - need to find derivation of this... ! for each unit going down, there are x units going up. do iv = currentPatch%nrad(L,ft),1, -1 - dif_ratio(L,ft,iv,ib) = dif_ratio(L,ft,iv+1,ib) * tran_dif(L,ft,iv,ib)*tran_dif(L,ft,iv,ib) / & - (1._r8 - dif_ratio(L,ft,iv+1,ib) * refl_dif(L,ft,iv,ib)) + refl_dif(L,ft,iv,ib) - dif_ratio(L,ft,iv,ib) = dif_ratio(L,ft,iv,ib) * ftweight(L,ft,iv)/ftweight(L,ft,1) - dif_ratio(L,ft,iv,ib) = dif_ratio(L,ft,iv,ib) + dif_ratio(L,ft,iv+1,ib)* & + dif_ratio(L,ft,iv,ib) = dif_ratio(L,ft,iv+1,ib) * & + tran_dif(L,ft,iv,ib)*tran_dif(L,ft,iv,ib) / & + (1._r8 - dif_ratio(L,ft,iv+1,ib) * refl_dif(L,ft,iv,ib)) & + + refl_dif(L,ft,iv,ib) + dif_ratio(L,ft,iv,ib) = dif_ratio(L,ft,iv,ib) * & + ftweight(L,ft,iv)/ftweight(L,ft,1) + dif_ratio(L,ft,iv,ib) = dif_ratio(L,ft,iv,ib) + dif_ratio(L,ft,iv+1,ib) * & (ftweight(L,ft,1)-ftweight(L,ft,iv))/ftweight(L,ft,1) end do - weighted_dif_ratio(L,ib) = weighted_dif_ratio(L,ib) + dif_ratio(L,ft,1,ib) * ftweight(L,ft,1) + weighted_dif_ratio(L,ib) = weighted_dif_ratio(L,ib) + & + dif_ratio(L,ft,1,ib) * ftweight(L,ft,1) !instance where the first layer ftweight is used a proxy for the whole column. FTWA end do!cp_numSWb endif ! currentPatch%present diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 075f797b5c..be53100a71 100755 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -52,11 +52,8 @@ subroutine fire_model( currentSite, atm2lnd_inst, temperature_inst) type (ed_patch_type), pointer :: currentPatch - integer temporary_SF_switch - !zero fire things currentPatch => currentSite%youngest_patch - temporary_SF_switch = 0 do while(associated(currentPatch)) currentPatch%frac_burnt = 0.0_r8 currentPatch%AB = 0.0_r8 @@ -68,7 +65,7 @@ subroutine fire_model( currentSite, atm2lnd_inst, temperature_inst) write(iulog,*) 'use_ed_spit_fire',use_ed_spit_fire endif - if(use_ed_spit_fire.and.temporary_SF_switch==1)then + if(use_ed_spit_fire)then call fire_danger_index(currentSite, temperature_inst, atm2lnd_inst) call wind_effect(currentSite, atm2lnd_inst) call charecteristics_of_fuel(currentSite) @@ -222,7 +219,7 @@ subroutine charecteristics_of_fuel ( currentSite ) ! average water content !is this the correct metric? timeav_swc = sum(currentSite%water_memory(1:10)) / 10._r8 ! Equation B2 in Thonicke et al. 2010 - fuel_moisture(dg_sf) = max(0.0_r8, 10.0_r8/9._r8 * timeav_swc - 1.0_r8/9.0_r8) + fuel_moisture(lg_sf) = max(0.0_r8, 10.0_r8/9._r8 * timeav_swc - 1.0_r8/9.0_r8) ! Average properties over the first four litter pools (dead leaves, twigs, s branches, l branches) currentPatch%fuel_bulkd = sum(currentPatch%fuel_frac(dg_sf:lb_sf) * SF_val_FBD(dg_sf:lb_sf)) @@ -363,7 +360,7 @@ subroutine wind_effect ( currentSite, atm2lnd_inst) do while(associated(currentPatch)) currentPatch%total_tree_area = min(currentPatch%total_tree_area,currentPatch%area) - currentPatch%effect_wspeed = wind * (tree_fraction*0.6+grass_fraction*0.4+bare_fraction*1.0) + currentPatch%effect_wspeed = wind * (tree_fraction*0.4+(grass_fraction+bare_fraction)*0.6) currentPatch => currentPatch%younger enddo !end patch loop diff --git a/main/EDMainMod.F90 b/main/EDMainMod.F90 index 559c4b91ba..9499f93d02 100755 --- a/main/EDMainMod.F90 +++ b/main/EDMainMod.F90 @@ -202,7 +202,8 @@ subroutine ed_integrate_state_variables(currentSite, temperature_inst ) currentCohort%bstore+udata%deltat* & (currentCohort%md+currentCohort%seed_prod)-cohort_biomass_store-currentCohort%npp_acc) endif - !do we need these any more? + + ! THESE SHOULD BE MOVED TO A MORE "VISIBLE" LOCATION (RGK 10-2016) currentCohort%npp_acc = 0.0_r8 currentCohort%gpp_acc = 0.0_r8 currentCohort%resp_acc = 0.0_r8 diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index 35fdb699b8..3419386adf 100755 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -198,16 +198,36 @@ module EDTypesMod ! because we don't want to continually re-calculate the cohort's position when ! performing size diagnostics at high-frequency calls + ! CARBON FLUXES - real(r8) :: gpp ! GPP: kgC/indiv/year - real(r8) :: gpp_acc ! GPP: kgC/indiv/day - real(r8) :: gpp_tstep ! GPP: kgC/indiv/timestep - real(r8) :: npp ! NPP: kgC/indiv/year - real(r8) :: npp_acc ! NPP: kgC/indiv/day - real(r8) :: npp_tstep ! NPP: kgC/indiv/timestep - real(r8) :: resp ! Resp: kgC/indiv/year - real(r8) :: resp_acc ! Resp: kgC/indiv/day - real(r8) :: resp_tstep ! Resp: kgC/indiv/timestep + + ! ---------------------------------------------------------------------------------- + ! NPP, GPP and RESP: Instantaneous, accumulated and accumulated-hold types.* + ! + ! _tstep: The instantaneous estimate that is calculated at each rapid plant biophysics + ! time-step (ie photosynthesis, sub-hourly). (kgC/indiv/timestep) + ! _acc: The accumulation of the _tstep variable from the beginning to ending of + ! the dynamics time-scale. This variable is zero'd during initialization and + ! after the dynamics call-sequence is completed. (kgC/indiv/day) + ! _acc_hold: While _acc is zero'd after the dynamics call sequence and then integrated, + ! _acc_hold "holds" the integrated value until the next time dynamics is + ! called. This is necessary for restarts. This variable also has units + ! converted to a useful rate (kgC/indiv/yr) + ! ---------------------------------------------------------------------------------- + + real(r8) :: gpp_tstep ! Gross Primary Production (see above *) + real(r8) :: gpp_acc + real(r8) :: gpp_acc_hold + + real(r8) :: npp_tstep ! Net Primary Production (see above *) + real(r8) :: npp_acc + real(r8) :: npp_acc_hold + + real(r8) :: resp_tstep ! Autotrophic respiration (see above *) + real(r8) :: resp_acc + real(r8) :: resp_acc_hold + + ! Net Primary Production Partitions real(r8) :: npp_leaf ! NPP into leaves (includes replacement of turnover): KgC/indiv/day real(r8) :: npp_froot ! NPP into fine roots (includes replacement of turnover): KgC/indiv/day diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index 59721fc313..b2b090b24a 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -792,9 +792,9 @@ subroutine update_history_dyn(this,nc,nsites,sites) associate( scpf => ccohort%size_by_pft_class ) hio_gpp_si_scpf(io_si,scpf) = hio_gpp_si_scpf(io_si,scpf) + & - n_perm2*ccohort%gpp ! [kgC/m2/yr] + n_perm2*ccohort%gpp_acc_hold ! [kgC/m2/yr] hio_npp_totl_si_scpf(io_si,scpf) = hio_npp_totl_si_scpf(io_si,scpf) + & - ccohort%npp*n_perm2 + ccohort%npp_acc_hold *n_perm2 hio_npp_leaf_si_scpf(io_si,scpf) = hio_npp_leaf_si_scpf(io_si,scpf) + & ccohort%npp_leaf*n_perm2 hio_npp_fnrt_si_scpf(io_si,scpf) = hio_npp_fnrt_si_scpf(io_si,scpf) + & @@ -812,15 +812,15 @@ subroutine update_history_dyn(this,nc,nsites,sites) hio_npp_stor_si_scpf(io_si,scpf) = hio_npp_stor_si_scpf(io_si,scpf) + & ccohort%npp_store*n_perm2 - if( abs(ccohort%npp-(ccohort%npp_leaf+ccohort%npp_froot+ & + if( abs(ccohort%npp_acc_hold-(ccohort%npp_leaf+ccohort%npp_froot+ & ccohort%npp_bsw+ccohort%npp_bdead+ & ccohort%npp_bseed+ccohort%npp_store))>1.e-9) then write(fates_log(),*) 'NPP Partitions are not balancing' write(fates_log(),*) 'Fractional Error: ', & - abs(ccohort%npp-(ccohort%npp_leaf+ccohort%npp_froot+ & + abs(ccohort%npp_acc_hold-(ccohort%npp_leaf+ccohort%npp_froot+ & ccohort%npp_bsw+ccohort%npp_bdead+ & - ccohort%npp_bseed+ccohort%npp_store))/ccohort%npp - write(fates_log(),*) 'Terms: ',ccohort%npp,ccohort%npp_leaf,ccohort%npp_froot, & + ccohort%npp_bseed+ccohort%npp_store))/ccohort%npp_acc_hold + write(fates_log(),*) 'Terms: ',ccohort%npp_acc_hold,ccohort%npp_leaf,ccohort%npp_froot, & ccohort%npp_bsw,ccohort%npp_bdead, & ccohort%npp_bseed,ccohort%npp_store write(fates_log(),*) ' NPP components during FATES-HLM linking does not balance ' diff --git a/main/FatesRestartInterfaceMod.F90 b/main/FatesRestartInterfaceMod.F90 index 81d9246c2f..18b77bc6cf 100644 --- a/main/FatesRestartInterfaceMod.F90 +++ b/main/FatesRestartInterfaceMod.F90 @@ -84,8 +84,8 @@ module FatesRestartInterfaceMod integer, private :: ir_nplant_co integer, private :: ir_gpp_acc_co integer, private :: ir_npp_acc_co - integer, private :: ir_gpp_co ! IS THIS VARIABLE NECESSARY? ... (RGK 11-2016) - integer, private :: ir_npp_co ! IS THIS VARIABLE NECESSARY? ... (RGK 11-2016) + integer, private :: ir_gpp_acc_hold_co + integer, private :: ir_npp_acc_hold_co integer, private :: ir_npp_leaf_co integer, private :: ir_npp_froot_co integer, private :: ir_npp_sw_co @@ -661,15 +661,15 @@ subroutine define_restart_vars(this, initialize_variables) units='kgC/indiv', flushval = flushzero, & hlms='CLM:ALM', initialize=initialize_variables, ivar=ivar, index = ir_npp_acc_co ) - call this%set_restart_var(vname='fates_gpp', vtype=cohort_r8, & + call this%set_restart_var(vname='fates_gpp_acc_hold', vtype=cohort_r8, & long_name='ed cohort - current step gpp', & units='kgC/indiv/year', flushval = flushzero, & - hlms='CLM:ALM', initialize=initialize_variables, ivar=ivar, index = ir_gpp_co ) + hlms='CLM:ALM', initialize=initialize_variables, ivar=ivar, index = ir_gpp_acc_hold_co ) - call this%set_restart_var(vname='fates_npp', vtype=cohort_r8, & + call this%set_restart_var(vname='fates_npp_acc_hold', vtype=cohort_r8, & long_name='ed cohort - current step npp', & units='kgC/indiv/year', flushval = flushzero, & - hlms='CLM:ALM', initialize=initialize_variables, ivar=ivar, index = ir_npp_co ) + hlms='CLM:ALM', initialize=initialize_variables, ivar=ivar, index = ir_npp_acc_hold_co ) call this%set_restart_var(vname='fates_npp_leaf', vtype=cohort_r8, & long_name='ed cohort - npp sent to leaves', & @@ -994,8 +994,8 @@ subroutine set_restart_vectors(this,nc,nsites,sites) rio_nplant_co => this%rvars(ir_nplant_co)%r81d, & rio_gpp_acc_co => this%rvars(ir_gpp_acc_co)%r81d, & rio_npp_acc_co => this%rvars(ir_npp_acc_co)%r81d, & - rio_gpp_co => this%rvars(ir_gpp_co)%r81d, & - rio_npp_co => this%rvars(ir_npp_co)%r81d, & + rio_gpp_acc_hold_co => this%rvars(ir_gpp_acc_hold_co)%r81d, & + rio_npp_acc_hold_co => this%rvars(ir_npp_acc_hold_co)%r81d, & rio_npp_leaf_co => this%rvars(ir_npp_leaf_co)%r81d, & rio_npp_froot_co => this%rvars(ir_npp_froot_co)%r81d, & rio_npp_sw_co => this%rvars(ir_npp_sw_co)%r81d, & @@ -1102,8 +1102,8 @@ subroutine set_restart_vectors(this,nc,nsites,sites) rio_nplant_co(io_idx_co) = ccohort%n rio_gpp_acc_co(io_idx_co) = ccohort%gpp_acc rio_npp_acc_co(io_idx_co) = ccohort%npp_acc - rio_gpp_co(io_idx_co) = ccohort%gpp - rio_npp_co(io_idx_co) = ccohort%npp + rio_gpp_acc_hold_co(io_idx_co) = ccohort%gpp_acc_hold + rio_npp_acc_hold_co(io_idx_co) = ccohort%npp_acc_hold rio_npp_leaf_co(io_idx_co) = ccohort%npp_leaf rio_npp_froot_co(io_idx_co) = ccohort%npp_froot rio_npp_sw_co(io_idx_co) = ccohort%npp_bsw @@ -1556,8 +1556,8 @@ subroutine get_restart_vectors(this, nc, nsites, sites) rio_nplant_co => this%rvars(ir_nplant_co)%r81d, & rio_gpp_acc_co => this%rvars(ir_gpp_acc_co)%r81d, & rio_npp_acc_co => this%rvars(ir_npp_acc_co)%r81d, & - rio_gpp_co => this%rvars(ir_gpp_co)%r81d, & - rio_npp_co => this%rvars(ir_npp_co)%r81d, & + rio_gpp_acc_hold_co => this%rvars(ir_gpp_acc_hold_co)%r81d, & + rio_npp_acc_hold_co => this%rvars(ir_npp_acc_hold_co)%r81d, & rio_npp_leaf_co => this%rvars(ir_npp_leaf_co)%r81d, & rio_npp_froot_co => this%rvars(ir_npp_froot_co)%r81d, & rio_npp_sw_co => this%rvars(ir_npp_sw_co)%r81d, & @@ -1649,8 +1649,8 @@ subroutine get_restart_vectors(this, nc, nsites, sites) ccohort%n = rio_nplant_co(io_idx_co) ccohort%gpp_acc = rio_gpp_acc_co(io_idx_co) ccohort%npp_acc = rio_npp_acc_co(io_idx_co) - ccohort%gpp = rio_gpp_co(io_idx_co) - ccohort%npp = rio_npp_co(io_idx_co) + ccohort%gpp_acc_hold = rio_gpp_acc_hold_co(io_idx_co) + ccohort%npp_acc_hold = rio_npp_acc_hold_co(io_idx_co) ccohort%npp_leaf = rio_npp_leaf_co(io_idx_co) ccohort%npp_froot = rio_npp_froot_co(io_idx_co) ccohort%npp_bsw = rio_npp_sw_co(io_idx_co)