diff --git a/biogeochem/EDCanopyStructureMod.F90 b/biogeochem/EDCanopyStructureMod.F90 index 240c373965..2877d093ba 100644 --- a/biogeochem/EDCanopyStructureMod.F90 +++ b/biogeochem/EDCanopyStructureMod.F90 @@ -939,10 +939,25 @@ subroutine UpdateFatesAvgSnowDepth(sites,bc_in) type(ed_site_type) , intent(inout), target :: sites(:) type(bc_in_type) , intent(in) :: bc_in(:) + ! Local + type (fates_patch_type) , pointer :: currentPatch + integer :: s + integer :: ifp do s = 1, size(sites,dim=1) - sites(s)%snow_depth = bc_in(s)%snow_depth_si * bc_in(s)%frac_sno_eff_si + + currentPatch => sites(s)%oldest_patch + + do while(associated(currentPatch)) + + ifp = currentPatch%patchno + currentPatch%snow_depth = sites(s)%bc_in(ifp)%snow_depth * sites(s)%bc_in(ifp)%frac_snow_eff + + currentPatch => currentPatch%younger + + end do + end do return @@ -1124,15 +1139,15 @@ subroutine leaf_area_profile( currentSite ) ( real(iv,r8)/currentCohort%NV * crown_depth ) fraction_exposed = 1.0_r8 - if(currentSite%snow_depth > layer_top_height)then + if(cpatch%snow_depth > layer_top_height)then fraction_exposed = 0._r8 endif - if(currentSite%snow_depth < layer_bottom_height)then + if(cpatch%snow_depth < layer_bottom_height)then fraction_exposed = 1._r8 endif - if(currentSite%snow_depth >= layer_bottom_height .and. & - currentSite%snow_depth <= layer_top_height) then !only partly hidden... - fraction_exposed = 1._r8 - max(0._r8,(min(1.0_r8,(currentSite%snow_depth -layer_bottom_height)/ & + if(cpatch%snow_depth >= layer_bottom_height .and. & + cpatch%snow_depth <= layer_top_height) then !only partly hidden... + fraction_exposed = 1._r8 - max(0._r8,(min(1.0_r8,(cpatch%snow_depth -layer_bottom_height)/ & (layer_top_height-layer_bottom_height )))) endif @@ -1170,7 +1185,7 @@ subroutine leaf_area_profile( currentSite ) currentCohort%treesai, & currentCohort%height, & iv,currentCohort%nv,currentCohort%pft, & - currentSite%snow_depth, & + cpatch%snow_depth, & vai_top,vai_bot, & elai_layer,esai_layer,tlai_layer,tsai_layer) diff --git a/biogeochem/EDCohortDynamicsMod.F90 b/biogeochem/EDCohortDynamicsMod.F90 index ee1e2d4a7c..81c303c579 100644 --- a/biogeochem/EDCohortDynamicsMod.F90 +++ b/biogeochem/EDCohortDynamicsMod.F90 @@ -486,7 +486,7 @@ subroutine terminate_cohort(currentSite, currentPatch, currentCohort, bc_in, ter if (currentCohort%n.gt.0.0_r8) then call SendCohortToLitter(currentSite,currentPatch, & - currentCohort,currentCohort%n,bc_in) + currentCohort,currentCohort%n) end if ! Set pointers and deallocate the current cohort from the list @@ -513,7 +513,7 @@ end subroutine terminate_cohort ! ===================================================================================== - subroutine SendCohortToLitter(csite,cpatch,ccohort,nplant,bc_in) + subroutine SendCohortToLitter(csite,cpatch,ccohort,nplant) ! ----------------------------------------------------------------------------------- ! This routine transfers the existing mass in all pools and all elements @@ -537,7 +537,6 @@ subroutine SendCohortToLitter(csite,cpatch,ccohort,nplant,bc_in) type (fates_cohort_type) , target :: ccohort real(r8) :: nplant ! Number (absolute) ! of plants to transfer - type(bc_in_type), intent(in) :: bc_in type(litter_type), pointer :: litt ! Litter object for each element type(elem_diag_type),pointer :: elflux_diags @@ -564,7 +563,7 @@ subroutine SendCohortToLitter(csite,cpatch,ccohort,nplant,bc_in) plant_dens = nplant/cpatch%area call set_root_fraction(csite%rootfrac_scr, pft, csite%zi_soil, & - bc_in%max_rooting_depth_index_col) + csite%bc_in(cpatch%patchno)%max_rooting_depth_index_col) do el=1,num_elements diff --git a/biogeochem/EDLoggingMortalityMod.F90 b/biogeochem/EDLoggingMortalityMod.F90 index a134257e10..d9da188f15 100644 --- a/biogeochem/EDLoggingMortalityMod.F90 +++ b/biogeochem/EDLoggingMortalityMod.F90 @@ -770,7 +770,7 @@ end subroutine get_harvest_rate_carbon ! ============================================================================ - subroutine logging_litter_fluxes(currentSite, currentPatch, newPatch, patch_site_areadis, bc_in) + subroutine logging_litter_fluxes(currentSite, currentPatch, newPatch, patch_site_areadis) ! ------------------------------------------------------------------------------------------- ! @@ -816,7 +816,6 @@ subroutine logging_litter_fluxes(currentSite, currentPatch, newPatch, patch_site type(fates_patch_type) , intent(inout), target :: currentPatch type(fates_patch_type) , intent(inout), target :: newPatch real(r8) , intent(in) :: patch_site_areadis - type(bc_in_type) , intent(in) :: bc_in !LOCAL VARIABLES: @@ -953,7 +952,7 @@ subroutine logging_litter_fluxes(currentSite, currentPatch, newPatch, patch_site call set_root_fraction(currentSite%rootfrac_scr, pft, & currentSite%zi_soil, & - bc_in%max_rooting_depth_index_col) + currentSite%bc_in(currentPatch%patchno)%max_rooting_depth_index_col) ag_wood = (direct_dead+indirect_dead) * (struct_m + sapw_m ) * & prt_params%allom_agb_frac(currentCohort%pft) diff --git a/biogeochem/EDPatchDynamicsMod.F90 b/biogeochem/EDPatchDynamicsMod.F90 index 30a41e8790..be00a49713 100644 --- a/biogeochem/EDPatchDynamicsMod.F90 +++ b/biogeochem/EDPatchDynamicsMod.F90 @@ -758,7 +758,7 @@ subroutine spawn_patches( currentSite, bc_in ) select case(i_disturbance_type) case (dtype_ilog) call logging_litter_fluxes(currentSite, currentPatch, & - newPatch, patch_site_areadis,bc_in) + newPatch, patch_site_areadis) ! if transitioning from primary to secondary, then may need to change nocomp pft, ! so tag as having transitioned LU @@ -767,13 +767,13 @@ subroutine spawn_patches( currentSite, bc_in ) end if case (dtype_ifire) call fire_litter_fluxes(currentSite, currentPatch, & - newPatch, patch_site_areadis,bc_in) + newPatch, patch_site_areadis) case (dtype_ifall) call mortality_litter_fluxes(currentSite, currentPatch, & - newPatch, patch_site_areadis,bc_in) + newPatch, patch_site_areadis) case (dtype_ilandusechange) call landusechange_litter_fluxes(currentSite, currentPatch, & - newPatch, patch_site_areadis,bc_in, & + newPatch, patch_site_areadis, & clearing_matrix(i_donorpatch_landuse_type,i_landusechange_receiverpatchlabel)) ! if land use change, then may need to change nocomp pft, so tag as having transitioned LU @@ -2142,7 +2142,7 @@ end subroutine TransLitterNewPatch ! ============================================================================ subroutine fire_litter_fluxes(currentSite, currentPatch, & - newPatch, patch_site_areadis, bc_in) + newPatch, patch_site_areadis) ! ! !DESCRIPTION: ! CWD pool burned by a fire. @@ -2161,7 +2161,6 @@ subroutine fire_litter_fluxes(currentSite, currentPatch, & type(fates_patch_type) , intent(inout), target :: currentPatch ! Donor Patch type(fates_patch_type) , intent(inout), target :: newPatch ! New Patch real(r8) , intent(in) :: patch_site_areadis ! Area being donated - type(bc_in_type) , intent(in) :: bc_in ! ! !LOCAL VARIABLES: @@ -2304,7 +2303,7 @@ subroutine fire_litter_fluxes(currentSite, currentPatch, & site_mass%burn_flux_to_atm = site_mass%burn_flux_to_atm + burned_mass call set_root_fraction(currentSite%rootfrac_scr, pft, currentSite%zi_soil, & - bc_in%max_rooting_depth_index_col) + currentSite%bc_in(currentPatch%patchno)%max_rooting_depth_index_col) ! Contribution of dead trees to root litter (no root burn flux to atm) do dcmpy=1,ndcmpy @@ -2381,7 +2380,7 @@ end subroutine fire_litter_fluxes ! ============================================================================ subroutine mortality_litter_fluxes(currentSite, currentPatch, & - newPatch, patch_site_areadis,bc_in) + newPatch, patch_site_areadis) ! ! !DESCRIPTION: ! Carbon going from mortality associated with disturbance into CWD pools. @@ -2403,7 +2402,6 @@ subroutine mortality_litter_fluxes(currentSite, currentPatch, & type(fates_patch_type) , intent(inout), target :: currentPatch type(fates_patch_type) , intent(inout), target :: newPatch real(r8) , intent(in) :: patch_site_areadis - type(bc_in_type) , intent(in) :: bc_in ! ! !LOCAL VARIABLES: type(fates_cohort_type), pointer :: currentCohort @@ -2531,7 +2529,7 @@ subroutine mortality_litter_fluxes(currentSite, currentPatch, & bg_wood = num_dead * (struct_m + sapw_m) * (1.0_r8-prt_params%allom_agb_frac(pft)) call set_root_fraction(currentSite%rootfrac_scr, pft, currentSite%zi_soil, & - bc_in%max_rooting_depth_index_col) + currentSite%bc_in(currentPatch%patchno)%max_rooting_depth_index_col) ! Adjust how wood is partitioned between the cwd classes based on cohort dbh call adjust_SF_CWD_frac(currentCohort%dbh,ncwd,SF_val_CWD_frac,SF_val_CWD_frac_adj) @@ -2614,8 +2612,8 @@ end subroutine mortality_litter_fluxes ! ============================================================================ subroutine landusechange_litter_fluxes(currentSite, currentPatch, & - newPatch, patch_site_areadis, bc_in, & - clearing_matrix_element) + newPatch, patch_site_areadis, clearing_matrix_element) + ! ! !DESCRIPTION: ! CWD pool from land use change. @@ -2630,7 +2628,6 @@ subroutine landusechange_litter_fluxes(currentSite, currentPatch, & type(fates_patch_type) , intent(inout), target :: currentPatch ! Donor Patch type(fates_patch_type) , intent(inout), target :: newPatch ! New Patch real(r8) , intent(in) :: patch_site_areadis ! Area being donated - type(bc_in_type) , intent(in) :: bc_in logical , intent(in) :: clearing_matrix_element ! whether or not to clear vegetation ! @@ -2771,7 +2768,7 @@ subroutine landusechange_litter_fluxes(currentSite, currentPatch, & site_mass%burn_flux_to_atm = site_mass%burn_flux_to_atm + burned_mass call set_root_fraction(currentSite%rootfrac_scr, pft, currentSite%zi_soil, & - bc_in%max_rooting_depth_index_col) + currentSite%bc_in(currentPatch%patchno)%max_rooting_depth_index_col) ! Contribution of dead trees to root litter (no root burn flux to atm) do dcmpy=1,ndcmpy diff --git a/biogeochem/EDPhysiologyMod.F90 b/biogeochem/EDPhysiologyMod.F90 index 601bb48fe8..84dbd48960 100644 --- a/biogeochem/EDPhysiologyMod.F90 +++ b/biogeochem/EDPhysiologyMod.F90 @@ -425,7 +425,7 @@ end subroutine GenerateDamageAndLitterFluxes ! ============================================================================ - subroutine PreDisturbanceLitterFluxes( currentSite, currentPatch, bc_in ) + subroutine PreDisturbanceLitterFluxes( currentSite, currentPatch) ! ----------------------------------------------------------------------------------- ! @@ -446,7 +446,6 @@ subroutine PreDisturbanceLitterFluxes( currentSite, currentPatch, bc_in ) ! !ARGUMENTS type(ed_site_type), intent(inout) :: currentSite type(fates_patch_type), intent(inout) :: currentPatch - type(bc_in_type), intent(in) :: bc_in ! ! !LOCAL VARIABLES: @@ -476,12 +475,12 @@ subroutine PreDisturbanceLitterFluxes( currentSite, currentPatch, bc_in ) ! Send fluxes from newly created litter into the litter pools ! This litter flux is from non-disturbance inducing mortality, as well ! as litter fluxes from live trees - call CWDInput(currentSite, currentPatch, litt,bc_in) + call CWDInput(currentSite, currentPatch, litt) ! Only calculate fragmentation flux over layers that are active ! (RGK-Mar2019) SHOULD WE MAX THIS AT 1? DONT HAVE TO - nlev_eff_decomp = max(bc_in%max_rooting_depth_index_col, 1) + nlev_eff_decomp = max(currentSite%bc_in(currentPatch%patchno)%max_rooting_depth_index_col, 1) call CWDOut(litt,currentPatch%fragmentation_scaler,nlev_eff_decomp) ! Fragmentation flux to soil decomposition model [kg/site/day] @@ -924,6 +923,7 @@ subroutine phenology( currentSite, bc_in ) integer :: i_wmem ! Loop counter for water mem days integer :: i_tmem ! Loop counter for veg temp mem days integer :: ipft ! plant functional type index + integer :: ifp ! fates patch index integer :: j ! Soil layer index real(r8) :: mean_10day_liqvol ! mean soil liquid volume over last 10 days [m3/m3] real(r8) :: mean_10day_smp ! mean soil matric potential over last 10 days [mm] @@ -1187,10 +1187,14 @@ subroutine phenology( currentSite, bc_in ) currentSite%smp_memory (i_wmem,ipft) = currentSite%smp_memory (i_wmem-1,ipft) end do + ! Temporarily set the bc index to one for multi-column fates refactor + ifp = 1 + ! Find the rooting depth distribution for PFT call set_root_fraction( currentSite%rootfrac_scr, ipft, currentSite%zi_soil, & - bc_in%max_rooting_depth_index_col ) - nlevroot = max(2,min(ubound(currentSite%zi_soil,1),bc_in%max_rooting_depth_index_col)) + currentSite%bc_in(ifp)%max_rooting_depth_index_col ) + nlevroot = max(2,min(ubound(currentSite%zi_soil,1), & + currentSite%bc_in(ifp)%max_rooting_depth_index_col)) ! The top most layer is typically very thin (~ 2cm) and dries rather quickly. Despite ! being thin, it can have a non-negligible rooting fraction (e.g., using @@ -2799,7 +2803,7 @@ end subroutine recruitment ! ====================================================================================== - subroutine CWDInput( currentSite, currentPatch, litt, bc_in) + subroutine CWDInput( currentSite, currentPatch, litt ) ! ! !DESCRIPTION: @@ -2818,7 +2822,6 @@ subroutine CWDInput( currentSite, currentPatch, litt, bc_in) type(ed_site_type), intent(inout), target :: currentSite type(fates_patch_type),intent(inout), target :: currentPatch type(litter_type),intent(inout),target :: litt - type(bc_in_type),intent(in) :: bc_in ! ! !LOCAL VARIABLES: @@ -2893,7 +2896,7 @@ subroutine CWDInput( currentSite, currentPatch, litt, bc_in) pft = currentCohort%pft call set_root_fraction(currentSite%rootfrac_scr, pft, currentSite%zi_soil, & - bc_in%max_rooting_depth_index_col) + currentSite%bc_in(currentPatch%patchno)%max_rooting_depth_index_col) store_m_turnover = currentCohort%prt%GetTurnover(store_organ,element_id) fnrt_m_turnover = currentCohort%prt%GetTurnover(fnrt_organ,element_id) diff --git a/biogeochem/FatesPatchMod.F90 b/biogeochem/FatesPatchMod.F90 index e673110118..20e5568426 100644 --- a/biogeochem/FatesPatchMod.F90 +++ b/biogeochem/FatesPatchMod.F90 @@ -128,6 +128,7 @@ module FatesPatchMod real(r8) :: total_tree_area ! area that is covered by woody vegetation [m2] real(r8) :: total_grass_area ! area that is covered by non-woody vegetation [m2] real(r8) :: zstar ! height of smallest canopy tree, only meaningful in "strict PPA" mode [m] + real(r8) :: snow_depth ! patch-level snow depth (used for ELAI/TLAI calcs) ! exposed leaf area in each canopy layer, pft, and leaf layer [m2 leaf/m2 contributing crown area] real(r8), allocatable :: elai_profile(:,:,:) ! nclmax,maxpft,nlevleaf) diff --git a/biogeochem/FatesSoilBGCFluxMod.F90 b/biogeochem/FatesSoilBGCFluxMod.F90 index 805edff6cf..3c0c235d26 100644 --- a/biogeochem/FatesSoilBGCFluxMod.F90 +++ b/biogeochem/FatesSoilBGCFluxMod.F90 @@ -331,7 +331,7 @@ subroutine PrepCH4BCs(csite,bc_in,bc_out) pft = ccohort%pft call set_root_fraction(csite%rootfrac_scr, pft, csite%zi_soil, & - bc_in%max_rooting_depth_index_col ) + csite%bc_in(cpatch%patchno)%max_rooting_depth_index_col ) fnrt_c = ccohort%prt%GetState(fnrt_organ, carbon12_element) @@ -488,7 +488,7 @@ subroutine PrepNutrientAquisitionBCs(csite, bc_in, bc_out) bc_out%ft_index(icomp) = pft call set_root_fraction(csite%rootfrac_scr, pft, csite%zi_soil, & - bc_in%max_rooting_depth_index_col ) + csite%bc_in(cpatch%patchno)%max_rooting_depth_index_col ) fnrt_c = ccohort%prt%GetState(fnrt_organ, carbon12_element) @@ -541,7 +541,7 @@ end subroutine PrepNutrientAquisitionBCs ! ===================================================================================== - subroutine EffluxIntoLitterPools(csite, cpatch, ccohort, bc_in ) + subroutine EffluxIntoLitterPools(csite, cpatch, ccohort) ! ----------------------------------------------------------------------------------- ! This subroutine just handles the transfer of exudation/efflux from plants @@ -554,7 +554,6 @@ subroutine EffluxIntoLitterPools(csite, cpatch, ccohort, bc_in ) type(ed_site_type), intent(inout) :: csite type(fates_patch_type), intent(inout) :: cpatch type(fates_cohort_type), intent(inout),target :: ccohort - type(bc_in_type), intent(in) :: bc_in ! locals integer :: el ! element loop index @@ -564,7 +563,7 @@ subroutine EffluxIntoLitterPools(csite, cpatch, ccohort, bc_in ) call set_root_fraction(csite%rootfrac_scr, & ccohort%pft, csite%zi_soil, & - bc_in%max_rooting_depth_index_col ) + csite%bc_in(cpatch%patchno)%max_rooting_depth_index_col ) ! Loop over the different elements. do el = 1, num_elements diff --git a/biogeophys/EDBtranMod.F90 b/biogeophys/EDBtranMod.F90 index 177557aedf..f9ab26db4f 100644 --- a/biogeophys/EDBtranMod.F90 +++ b/biogeophys/EDBtranMod.F90 @@ -124,6 +124,8 @@ subroutine btran_ed( nsites, sites, bc_in, bc_out) real(r8) :: temprootr real(r8) :: sum_pftgs ! sum of weighted conductances (for normalization) real(r8), allocatable :: root_resis(:,:) ! Root resistance in each pft x layer + real(r8) :: effective_porosity ! Effective porosity in each soil layer + real(r8) :: water_saturation ! Water saturation in each soil layer !------------------------------------------------------------------------------ associate( & @@ -149,7 +151,7 @@ subroutine btran_ed( nsites, sites, bc_in, bc_out) do ft = 1,numpft call set_root_fraction(sites(s)%rootfrac_scr, ft, sites(s)%zi_soil, & - bc_in(s)%max_rooting_depth_index_col ) + sites(s)%bc_in(cpatch%patchno)%max_rooting_depth_index_col ) cpatch%btran_ft(ft) = 0.0_r8 do j = 1,bc_in(s)%nlevsoil @@ -160,8 +162,18 @@ subroutine btran_ed( nsites, sites, bc_in, bc_out) if ( check_layer_water(bc_in(s)%h2o_liqvol_sl(j),bc_in(s)%tempk_sl(j)) ) then smp_node = max(smpsc(ft), bc_in(s)%smp_sl(j)) - - rresis = min( (bc_in(s)%eff_porosity_sl(j)/bc_in(s)%watsat_sl(j))* & + + ! Check that the patch has exposed vegetation + ! If it does, locally override the inbound effective porosity values from the host + ! Note that filter_btran will need to be converted to handle MCF + effective_porosity = sites(s)%bc_in(ifp)%eff_porosity_sl(j) + water_saturation = sites(s)%bc_in(ifp)%watsat_sl(j) + if (.not. bc_in(s)%filter_btran) then + effective_porosity = -999._r8 + water_saturation = -999._r8 + end if + + rresis = min( (effective_porosity/water_saturation) * & (smp_node - smpsc(ft)) / (smpso(ft) - smpsc(ft)), 1._r8) root_resis(ft,j) = sites(s)%rootfrac_scr(j)*rresis @@ -247,6 +259,7 @@ subroutine btran_ed( nsites, sites, bc_in, bc_out) enddo end if + endif if_bare cpatch => cpatch%younger end do diff --git a/biogeophys/FatesBstressMod.F90 b/biogeophys/FatesBstressMod.F90 index f37ab8ccb1..d51a9c497b 100644 --- a/biogeophys/FatesBstressMod.F90 +++ b/biogeophys/FatesBstressMod.F90 @@ -68,7 +68,7 @@ subroutine btran_sal_stress_fates( nsites, sites, bc_in) call set_root_fraction(sites(s)%rootfrac_scr, ft, & sites(s)%zi_soil, & - bc_in(s)%max_rooting_depth_index_col ) + sites(s)%bc_in(cpatch%patchno)%max_rooting_depth_index_col ) do j = 1,bc_in(s)%nlevsoil diff --git a/biogeophys/FatesPlantHydraulicsMod.F90 b/biogeophys/FatesPlantHydraulicsMod.F90 index 1672ec5ff8..31f8651520 100644 --- a/biogeophys/FatesPlantHydraulicsMod.F90 +++ b/biogeophys/FatesPlantHydraulicsMod.F90 @@ -1530,6 +1530,7 @@ subroutine HydrSiteColdStart(sites, bc_in ) integer :: j,j_t,j_b integer :: nsites integer :: nlevrhiz + integer :: ifp class(wrf_type_vg), pointer :: wrf_vg class(wkf_type_vg), pointer :: wkf_vg class(wrf_type_cch), pointer :: wrf_cch @@ -1541,13 +1542,15 @@ subroutine HydrSiteColdStart(sites, bc_in ) do s = 1,nsites + ifp = 1 + csite_hydr => sites(s)%si_hydr nlevrhiz = csite_hydr%nlevrhiz do j = 1,nlevrhiz j_t = csite_hydr%map_r2s(j,1) ! top soil layer matching rhiz layer j_b = csite_hydr%map_r2s(j,2) ! bottom soil layer matching rhiz layer - eff_por = csite_hydr%AggBCToRhiz(bc_in(s)%eff_porosity_sl,j,bc_in(s)%dz_sisl) + eff_por = csite_hydr%AggBCToRhiz(sites(s)%bc_in(ifp)%eff_porosity_sl,j,bc_in(s)%dz_sisl) ! [kg/m2] / ([m] * [kg/m3]) = [m3/m3] h2osoi_liqvol = min(eff_por, & @@ -2756,12 +2759,14 @@ subroutine hydraulics_bc ( nsites, sites, bc_in, bc_out, dtime) sumweight = 0._r8 do j_bc = j_t,j_b + ! Temporarily set ifp index for mcf update + ifp = 1 if(rootflux_disagg == soilk_disagg)then if(qflx_soil2root_rhiz>0._r8 )then ! Weight disaggregation by K*dz, but only for flux ! into the root, othersize weight by root length ! h2osoi_liqvol: [kg/m2] / [m] / [kg/m3] = [m3/m3] - eff_por = bc_in(s)%eff_porosity_sl(j_bc) + eff_por = sites(s)%bc_in(ifp)%eff_porosity_sl(j_bc) h2osoi_liqvol = min(eff_por, bc_in(s)%h2o_liq_sisl(j_bc)/(bc_in(s)%dz_sisl(j_bc)*denh2o)) psi_layer = csite_hydr%wrf_soil(j)%p%psi_from_th(h2osoi_liqvol) ftc_layer = csite_hydr%wkf_soil(j)%p%ftc_from_psi(psi_layer) diff --git a/biogeophys/FatesPlantRespPhotosynthMod.F90 b/biogeophys/FatesPlantRespPhotosynthMod.F90 index 60cf234a64..4921162eb3 100644 --- a/biogeophys/FatesPlantRespPhotosynthMod.F90 +++ b/biogeophys/FatesPlantRespPhotosynthMod.F90 @@ -325,18 +325,17 @@ subroutine FatesPlantRespPhotosynthDrive (nsites, sites,bc_in,bc_out,dtime) allocate(rootfr_ft(numpft, bc_in(s)%nlevsoil)) - do ft = 1,numpft - call set_root_fraction(rootfr_ft(ft,:), ft, & - bc_in(s)%zi_sisl, & - bc_in(s)%max_rooting_depth_index_col) - end do - - currentpatch => sites(s)%oldest_patch do while (associated(currentpatch)) ifp = currentPatch%patchno + do ft = 1,numpft + call set_root_fraction(rootfr_ft(ft,:), ft, & + bc_in(s)%zi_sisl, & + sites(s)%bc_in(ifp)%max_rooting_depth_index_col ) + end do + if_notbare: if(currentpatch%nocomp_pft_label.ne.nocomp_bareground)then NCL_p = currentPatch%NCL_p @@ -437,7 +436,7 @@ subroutine FatesPlantRespPhotosynthDrive (nsites, sites,bc_in,bc_out,dtime) iv, & currentCohort%nv, & currentCohort%pft, & - sites(s)%snow_depth, & + currentPatch%snow_depth, & cohort_vaitop(iv), & cohort_vaibot(iv), & cohort_layer_elai(iv), & diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 7e6de3323b..d8ba644eab 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -220,6 +220,8 @@ subroutine CalculateIgnitionsandFDI(currentSite, bc_in) ! initialize site parameters to zero currentSite%NF_successful = 0.0_r8 + currentSite%rxfire_area_fi = 0.0_r8 + currentSite%rxfire_area_fuel = 0.0_r8 ! Equation 7 from Venevsky et al GCB 2002 (modification of equation 8 in Thonicke et al. 2010) ! FDI 0.1 = low, 0.3 moderate, 0.75 high, and 1 = extreme ignition potential for alpha 0.000337 diff --git a/main/EDInitMod.F90 b/main/EDInitMod.F90 index 260e1cf63f..1246de5e3c 100644 --- a/main/EDInitMod.F90 +++ b/main/EDInitMod.F90 @@ -297,7 +297,6 @@ subroutine zero_site( site_in ) site_in%cstatus = fates_unset_int ! are leaves in this pixel on or off? site_in%dstatus(:) = fates_unset_int site_in%grow_deg_days = nan ! growing degree days - site_in%snow_depth = nan site_in%nchilldays = fates_unset_int site_in%ncolddays = fates_unset_int site_in%cleafondate = fates_unset_int ! doy of leaf on (cold) diff --git a/main/EDMainMod.F90 b/main/EDMainMod.F90 index d3914b017c..171157ae80 100644 --- a/main/EDMainMod.F90 +++ b/main/EDMainMod.F90 @@ -648,7 +648,7 @@ subroutine ed_integrate_state_variables(currentSite, bc_in, bc_out ) ! Then zero out the daily uptakes, they have been used - call EffluxIntoLitterPools(currentSite, currentPatch, currentCohort, bc_in ) + call EffluxIntoLitterPools(currentSite, currentPatch, currentCohort) if(element_pos(nitrogen_element)>0) then ! Mass balance for N uptake @@ -775,7 +775,7 @@ subroutine ed_integrate_state_variables(currentSite, bc_in, bc_out ) call GenerateDamageAndLitterFluxes( currentSite, currentPatch) - call PreDisturbanceLitterFluxes( currentSite, currentPatch, bc_in) + call PreDisturbanceLitterFluxes( currentSite, currentPatch ) call PreDisturbanceIntegrateLitter(currentPatch ) diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index bd50db4d8d..8e6fe3db03 100644 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -422,7 +422,6 @@ module EDTypesMod ! PHENOLOGY real(r8) :: grow_deg_days ! Phenology growing degree days - real(r8) :: snow_depth ! site-level snow depth (used for ELAI/TLAI calcs) integer :: cstatus ! are leaves in this pixel on or off for cold decid ! 0 = this site has not experienced a cold period over at least diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index ee8cc1ed3c..106de2687b 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -2366,7 +2366,7 @@ subroutine update_history_dyn(this,nc,nsites,sites,bc_in) if(hlm_hist_level_dynam>0) then call update_history_dyn_sitelevel(this,nc,nsites,sites) if(hlm_hist_level_dynam>1) then - call update_history_dyn_subsite(this,nc,nsites,sites,bc_in) + call update_history_dyn_subsite(this,nc,nsites,sites) call update_history_dyn_subsite_ageclass(this,nc,nsites,sites) call reset_history_dyn_subsite(this, nsites, sites) end if @@ -3063,7 +3063,7 @@ end subroutine update_history_dyn_sitelevel ! ========================================================================================= - subroutine update_history_dyn_subsite(this,nc,nsites,sites,bc_in) + subroutine update_history_dyn_subsite(this,nc,nsites,sites) ! --------------------------------------------------------------------------------- ! This subroutine is intended to update all history variables with upfreq == @@ -3077,7 +3077,6 @@ subroutine update_history_dyn_subsite(this,nc,nsites,sites,bc_in) integer , intent(in) :: nc ! clump index integer , intent(in) :: nsites type(ed_site_type) , intent(inout), target :: sites(nsites) - type(bc_in_type) , intent(in) :: bc_in(nsites) type(fates_cohort_type), pointer :: ccohort type(fates_patch_type), pointer :: cpatch @@ -3534,7 +3533,7 @@ subroutine update_history_dyn_subsite(this,nc,nsites,sites,bc_in) endif call set_root_fraction(sites(s)%rootfrac_scr, ccohort%pft, sites(s)%zi_soil, & - bc_in(s)%max_rooting_depth_index_col ) + sites(s)%bc_in(cpatch%patchno)%max_rooting_depth_index_col ) ! Update biomass components ! Mass pools [kg] @@ -5150,8 +5149,9 @@ subroutine update_history_hifrq_sitelevel(this,nc,nsites,sites,bc_in,dt_tstep) real(r8) , intent(in) :: dt_tstep ! Locals - integer :: s ! The local site index - integer :: io_si ! The site index of the IO array + integer :: s ! The local site index + integer :: io_si ! The site index of the IO array + integer :: ifp ! fates patch index integer :: age_class ! class age index real(r8) :: site_area_veg_inv ! inverse canopy area of the site (1/m2) real(r8) :: site_area_rad_inv ! inverse canopy area of site for only @@ -5201,9 +5201,12 @@ subroutine update_history_hifrq_sitelevel(this,nc,nsites,sites,bc_in,dt_tstep) call this%zero_site_hvars(sites(s), upfreq_in=group_hifr_simple) io_si = sites(s)%h_gid + + ! Temporarily hard code ifp to 1. This will be updated for future multicolumn fates + ifp = 1 - hio_nep_si(io_si) = -bc_in(s)%tot_het_resp * kg_per_g - hio_hr_si(io_si) = bc_in(s)%tot_het_resp * kg_per_g + hio_nep_si(io_si) = -sites(s)%bc_in(ifp)%tot_het_resp * kg_per_g + hio_hr_si(io_si) = sites(s)%bc_in(ifp)%tot_het_resp * kg_per_g ! Diagnostics that are only incremented if we called the radiation solver ! We do not call the radiation solver if @@ -5888,6 +5891,7 @@ subroutine update_history_hydraulics(this,nc,nsites,sites,bc_in,dt_tstep) integer :: j_t,j_b ! top and bottom soil layer matching current rhiz layer integer :: nlevrhiz ! number of rhizosphere layers integer :: nlevsoil ! number of soil layers + integer :: ifp ! index of the fates patch number real(r8) :: mean_soil_vwc ! mean soil volumetric water content [m3/m3] real(r8) :: mean_soil_vwcsat ! mean soil saturated volumetric water content [m3/m3] real(r8) :: mean_soil_matpot ! mean soil water potential [MPa] @@ -5926,6 +5930,10 @@ subroutine update_history_hydraulics(this,nc,nsites,sites,bc_in,dt_tstep) do s = 1,nsites + ! Prior to multi-column fates conversion, we can use the fates patch number + ! associated with the first patch on the site for the boundary condition indexing. + ifp = 1 + call this%zero_site_hvars(sites(s),upfreq_in=group_hydr_simple) site_hydr => sites(s)%si_hydr @@ -5960,9 +5968,14 @@ subroutine update_history_hydraulics(this,nc,nsites,sites,bc_in,dt_tstep) ! with model internals and physics. Should ! implement caps inside the functions ! if desired. (RGK 12-2021) - vwc_sat = bc_in(s)%watsat_sl(j_bc) + vwc_sat = sites(s)%bc_in(ifp)%watsat_sl(j_bc) depth_frac = bc_in(s)%dz_sisl(j_bc)/site_hydr%dz_rhiz(j) + ! Override the watsat_sl if patch doesn't have exposed vegetation + if (.not. bc_in(s)%filter_btran) then + vwc_sat = -999._r8 + end if + ! If there are any roots, we use root weighting if(sum(site_hydr%l_aroot_layer(:),dim=1) > nearzero) then layer_areaweight = site_hydr%l_aroot_layer(j)*depth_frac*pi_const*site_hydr%rs1(j)**2.0 @@ -6034,6 +6047,10 @@ subroutine update_history_hydraulics(this,nc,nsites,sites,bc_in,dt_tstep) call this%flush_hvars(nc,upfreq_in=group_hydr_complx) do s = 1,nsites + + ! Prior to multi-column fates conversion, we can use the fates patch number + ! associated with the first patch on the site for the boundary condition indexing. + ifp = 1 site_hydr => sites(s)%si_hydr nlevrhiz = site_hydr%nlevrhiz @@ -6067,9 +6084,14 @@ subroutine update_history_hydraulics(this,nc,nsites,sites,bc_in,dt_tstep) ! with model internals and physics. Should ! implement caps inside the functions ! if desired. (RGK 12-2021) - vwc_sat = bc_in(s)%watsat_sl(j_bc) + vwc_sat = sites(s)%bc_in(ifp)%watsat_sl(j_bc) depth_frac = bc_in(s)%dz_sisl(j_bc)/site_hydr%dz_rhiz(j) + ! Override the watsat_sl if patch doesn't have exposed vegetation + if (.not. bc_in(s)%filter_btran) then + vwc_sat = -999._r8 + end if + ! If there are any roots, we use root weighting if(sum(site_hydr%l_aroot_layer(:),dim=1) > nearzero) then layer_areaweight = site_hydr%l_aroot_layer(j)*depth_frac*pi_const*site_hydr%rs1(j)**2.0 diff --git a/main/FatesInterfaceMod.F90 b/main/FatesInterfaceMod.F90 index 264f1c82e1..3d210b5b04 100644 --- a/main/FatesInterfaceMod.F90 +++ b/main/FatesInterfaceMod.F90 @@ -190,6 +190,7 @@ module FatesInterfaceMod procedure :: SetRegistryActiveState procedure :: SetRegistryLastState procedure, public :: UpdateInterfaceVariables + procedure, public :: UpdateInterfaceVariablesTimestep procedure, public :: UpdateLitterFluxes end type fates_interface_type @@ -318,19 +319,11 @@ subroutine zero_bcs(fates,s) fates%bc_in(s)%solad_parb(:,:) = 0.0_r8 fates%bc_in(s)%solai_parb(:,:) = 0.0_r8 fates%bc_in(s)%smp_sl(:) = 0.0_r8 - fates%bc_in(s)%eff_porosity_sl(:) = 0.0_r8 - fates%bc_in(s)%watsat_sl(:) = 0.0_r8 fates%bc_in(s)%tempk_sl(:) = 0.0_r8 fates%bc_in(s)%h2o_liqvol_sl(:) = 0.0_r8 fates%bc_in(s)%fcansno_pa(:) = 0.0_r8 fates%bc_in(s)%albgr_dir_rb(:) = 0.0_r8 fates%bc_in(s)%albgr_dif_rb(:) = 0.0_r8 - fates%bc_in(s)%max_rooting_depth_index_col = 0 - fates%bc_in(s)%tot_het_resp = 0.0_r8 - fates%bc_in(s)%tot_somc = 0.0_r8 - fates%bc_in(s)%tot_litc = 0.0_r8 - fates%bc_in(s)%snow_depth_si = 0.0_r8 - fates%bc_in(s)%frac_sno_eff_si = 0.0_r8 if(do_fates_salinity)then fates%bc_in(s)%salinity_sl(:) = 0.0_r8 @@ -542,8 +535,6 @@ subroutine allocate_bcin(bc_in, nlevsoil_in, nlevdecomp_in, num_lu_harvest_cats, ! Hydrology allocate(bc_in%smp_sl(nlevsoil_in)) - allocate(bc_in%eff_porosity_sl(nlevsoil_in)) - allocate(bc_in%watsat_sl(nlevsoil_in)) allocate(bc_in%tempk_sl(nlevsoil_in)) allocate(bc_in%h2o_liqvol_sl(nlevsoil_in)) @@ -2982,6 +2973,7 @@ subroutine InitializeBoundaryConditions(this, patches_per_site) r = this%sites(s)%GetRegistryIndex(ifp) ! Register the boundary conditions that are necessary for allocating other boundary conditions first + call this%registry(r)%Register(key=hlm_fates_nlevground, data=bc_in%nlevgrnd, hlm_flag=.false.) call this%registry(r)%Register(key=hlm_fates_decomp_max, data=bc_in%nlevdecomp_full, hlm_flag=.false.) call this%registry(r)%Register(key=hlm_fates_decomp, data=bc_in%nlevdecomp, hlm_flag=.false.) call this%registry(r)%Register(key=hlm_fates_soil_level, data=bc_in%nlevsoil, hlm_flag=.false.) @@ -3004,7 +2996,17 @@ subroutine InitializeBoundaryConditions(this, patches_per_site) data=bc_in%w_scalar_sisl, hlm_flag=.false.) call this%registry(r)%Register(key=hlm_fates_decomp_frac_temperature, & data=bc_in%t_scalar_sisl, hlm_flag=.false.) - + call this%registry(r)%Register(key=hlm_fates_effective_porosity, & + data=bc_in%eff_porosity_sl, hlm_flag=.false.) + call this%registry(r)%Register(key=hlm_fates_soil_water_saturation, & + data=bc_in%watsat_sl, hlm_flag=.false.) + call this%registry(r)%Register(key=hlm_fates_heterotrophic_respiration, & + data=bc_in%tot_het_resp, hlm_flag=.false.) + call this%registry(r)%Register(key=hlm_fates_snow_depth, & + data=bc_in%snow_depth, hlm_flag=.false.) + call this%registry(r)%Register(key=hlm_fates_snow_cover_frac, & + data=bc_in%frac_snow_eff, hlm_flag=.false.) + ! bc_out nlevdecomp = bc_in%nlevdecomp call this%registry(r)%Register(key=hlm_fates_litter_carbon_cellulose, & @@ -3154,6 +3156,33 @@ end subroutine UpdateInterfaceVariables ! ====================================================================================== +subroutine UpdateInterfaceVariablesTimeStep(this) + + ! This procedure handles updating the interface variables that need to be updated at + ! every model time step. + + ! Arguments + class(fates_interface_type), intent(inout) :: this + + ! Locals + integer :: n ! active registry index iterator + integer :: r ! registry index + + ! Set the registry active state + call this%SetRegistryActiveState() + + ! Loop through the active registries and update the litter fluxes + do n = 1, this%num_active_patches + r = this%filter_registry_active(n) + + call this%registry(r)%UpdateTimestep() + + end do + +end subroutine UpdateInterfaceVariablesTimeStep + +! ====================================================================================== + subroutine UpdateLitterFluxes(this, dtime) ! This procedure handles the updating of litter fluxes from FATES to the HLM. diff --git a/main/FatesInterfaceParametersMod.F90 b/main/FatesInterfaceParametersMod.F90 index 31e4d8c0ae..747e9bb072 100644 --- a/main/FatesInterfaceParametersMod.F90 +++ b/main/FatesInterfaceParametersMod.F90 @@ -11,6 +11,7 @@ module FatesInterfaceParametersMod character(len=*), parameter, public :: hlm_fates_soil_level = 'soil_level_number' character(len=*), parameter, public :: hlm_fates_decomp_frac_moisture = 'decomp_frac_moisture' character(len=*), parameter, public :: hlm_fates_decomp_frac_temperature = 'decomp_frac_temperature' + character(len=*), parameter, public :: hlm_fates_nlevground = 'nlevground' character(len=*), parameter, public :: hlm_fates_litter_carbon_cellulose = 'litter_carbon_cellulose' character(len=*), parameter, public :: hlm_fates_litter_carbon_lignin = 'litter_carbon_lignin' character(len=*), parameter, public :: hlm_fates_litter_carbon_labile = 'litter_carbon_labile' @@ -23,6 +24,11 @@ module FatesInterfaceParametersMod character(len=*), parameter, public :: hlm_fates_litter_nitrogen_lignin = 'litter_nitrogen_lignin' character(len=*), parameter, public :: hlm_fates_litter_nitrogen_labile = 'litter_nitrogen_labile' character(len=*), parameter, public :: hlm_fates_litter_nitrogen_total= 'litter_nitrogen_total' + character(len=*), parameter, public :: hlm_fates_effective_porosity = 'effective_porosity' + character(len=*), parameter, public :: hlm_fates_soil_water_saturation = 'soil_water_saturation' + character(len=*), parameter, public :: hlm_fates_heterotrophic_respiration = 'heterotrophic_respiration' + character(len=*), parameter, public :: hlm_fates_snow_depth = 'snow_depth' + character(len=*), parameter, public :: hlm_fates_snow_cover_frac = 'snow_cover_frac' ! Registry update frequency parameters integer, parameter, public :: registry_update_init_dims = 0 ! variable dimension that needs to be updated during initialization diff --git a/main/FatesInterfaceTypesMod.F90 b/main/FatesInterfaceTypesMod.F90 index 67567083e6..59c0dfe07d 100644 --- a/main/FatesInterfaceTypesMod.F90 +++ b/main/FatesInterfaceTypesMod.F90 @@ -392,6 +392,7 @@ module FatesInterfaceTypesMod ! Soil layer structure + integer :: nlevgrnd ! the number of ground layers in this column integer :: nlevsoil ! the number of soil layers in this column integer :: nlevdecomp ! the number of active soil layers in the column integer :: nlevdecomp_full ! the maximum possible soil layers for any column @@ -527,13 +528,11 @@ module FatesInterfaceTypesMod ! BGC Accounting real(r8) :: tot_het_resp ! total heterotrophic respiration (gC/m2/s) - real(r8) :: tot_somc ! total soil organic matter carbon (gc/m2) - real(r8) :: tot_litc ! total litter carbon tracked in the HLM (gc/m2) ! Canopy Structure - real(r8) :: snow_depth_si ! Depth of snow in snowy areas of site (m) - real(r8) :: frac_sno_eff_si ! Fraction of ground covered by snow (0-1) + real(r8) :: snow_depth ! Depth of snow in snowy areas of site (m) + real(r8) :: frac_snow_eff ! Fraction of ground covered by snow (0-1) ! Hydrology variables for BTRAN ! --------------------------------------------------------------------------------- @@ -930,6 +929,7 @@ module FatesInterfaceTypesMod procedure :: SetLastState procedure :: UpdateLitterFluxes procedure :: Update => UpdateInterfaceVariables + procedure :: UpdateTimeStep => UpdateInterfaceVariablesTimeStep generic :: Register => RegisterInterfaceVariables_0d, & RegisterInterfaceVariables_1d, & @@ -982,13 +982,19 @@ subroutine InitializeBCIn(this) allocate(this%dz_decomp_sisl(this%nlevdecomp_full)) allocate(this%w_scalar_sisl(this%nlevdecomp_full)) allocate(this%t_scalar_sisl(this%nlevdecomp_full)) + allocate(this%eff_porosity_sl(this%nlevgrnd)) + allocate(this%watsat_sl(this%nlevgrnd)) ! Unset variables this%decomp_id = fates_unset_int this%dz_decomp_sisl = nan this%w_scalar_sisl = nan this%t_scalar_sisl = nan + this%eff_porosity_sl = nan + this%watsat_sl = nan this%max_thaw_depth_index = fates_unset_int + this%snow_depth = nan + this%frac_snow_eff = nan end subroutine InitializeBCIn @@ -1166,6 +1172,8 @@ subroutine DefineInterfaceRegistry(this, initialize) ! Define the interface registry names and indices ! Variables that need to be updated during initialization and are necessary for other boundary conditions ! such as dimensions + call this%DefineInterfaceVariable(key=hlm_fates_nlevground, initialize=initialize, index=index, & + update_frequency=registry_update_init_dims) call this%DefineInterfaceVariable(key=hlm_fates_decomp_max, initialize=initialize, index=index, & update_frequency=registry_update_init_dims) call this%DefineInterfaceVariable(key=hlm_fates_decomp, initialize=initialize, index=index, & @@ -1191,6 +1199,16 @@ subroutine DefineInterfaceRegistry(this, initialize) update_frequency=registry_update_timestep, bc_dir=bc_out) call this%DefineInterfaceVariable(key=hlm_fates_litter_carbon_total, initialize=initialize, index=index, & update_frequency=registry_update_timestep, bc_dir=bc_out) + call this%DefineInterfaceVariable(key=hlm_fates_effective_porosity, initialize=initialize, index=index, & + update_frequency=registry_update_timestep, bc_dir=bc_in) + call this%DefineInterfaceVariable(key=hlm_fates_soil_water_saturation, initialize=initialize, index=index, & + update_frequency=registry_update_timestep, bc_dir=bc_in) + call this%DefineInterfaceVariable(key=hlm_fates_heterotrophic_respiration, initialize=initialize, index=index, & + update_frequency=registry_update_timestep, bc_dir=bc_in) + call this%DefineInterfaceVariable(key=hlm_fates_snow_depth, initialize=initialize, index=index, & + update_frequency=registry_update_timestep, bc_dir=bc_in) + call this%DefineInterfaceVariable(key=hlm_fates_snow_cover_frac, initialize=initialize, index=index, & + update_frequency=registry_update_timestep, bc_dir=bc_in) ! Define the N and P litter fluxes if in CNP mode ! We could define the interface variables always, even if not registered, but this helps reduce the memory needs @@ -1967,6 +1985,35 @@ subroutine UpdateLitterFluxes(this) end subroutine UpdateLitterFluxes + ! ====================================================================================== + + subroutine UpdateInterfaceVariablesTimeStep(this) + + ! This procedure updates the interface variables that are updated on the model time-step. + + class(fates_interface_registry_type), intent(inout) :: this ! registry being updated + + integer :: i ! update iterator + integer :: j ! variable index + + ! Update the boundary conditions necessary during time-step updates only + do i = 1, this%num_api_vars_update_timestep + + ! Get the variable index from the filter + j = this%filter_timestep(i) + + ! Skip updating litter flux variables as they are handled via a separate update call + if (.not. (any(this%filter_litter_flux(:) == j) .or. & + this%key(j) == hlm_fates_litter_carbon_total .or. & + this%key(j) == hlm_fates_litter_phosphorus_total .or. & + this%key(j) == hlm_fates_litter_nitrogen_total)) then + call this%fates_vars(j)%Update(this%hlm_vars(j)) + end if + + end do + + end subroutine UpdateInterfaceVariablesTimeStep + ! ====================================================================================== integer function GetRegistryVariableIndex(this, key) result(index) diff --git a/main/FatesRestartInterfaceMod.F90 b/main/FatesRestartInterfaceMod.F90 index 510defe04c..8a4380c5a7 100644 --- a/main/FatesRestartInterfaceMod.F90 +++ b/main/FatesRestartInterfaceMod.F90 @@ -111,7 +111,7 @@ module FatesRestartInterfaceMod integer :: ir_min_allowed_landuse_fraction_si integer :: ir_landuse_vector_gt_min_si integer :: ir_area_bareground_si - integer :: ir_snow_depth_si + integer :: ir_snow_depth_pa integer :: ir_landuse_config_si integer :: ir_gpp_acc_si integer :: ir_aresp_acc_si @@ -768,10 +768,6 @@ subroutine define_restart_vars(this, initialize_variables) long_name='minimum allowed land use fraction at each site', units='fraction', flushval = flushzero, & hlms='CLM:ALM', initialize=initialize_variables, ivar=ivar, index = ir_area_bareground_si ) - call this%set_restart_var(vname='fates_snow_depth_site', vtype=site_r8, & - long_name='average snow depth', units='m', flushval = flushzero, & - hlms='CLM:ALM', initialize=initialize_variables, ivar=ivar, index = ir_snow_depth_si ) - call this%set_restart_var(vname='fates_landuse_config_site', vtype=site_int, & long_name='hlm_use_potentialveg status of run that created this restart file', & units='kgC/m2', flushval = flushzero, & @@ -803,6 +799,10 @@ subroutine define_restart_vars(this, initialize_variables) long_name='Fraction of canopy covered in snow', units='unitless', flushval = flushinvalid, & hlms='CLM:ALM', initialize=initialize_variables, ivar=ivar, index = ir_fcansno_pa ) + call this%set_restart_var(vname='fates_snow_depth_pa', vtype=cohort_r8, & + long_name='average snow depth', units='m', flushval = flushzero, & + hlms='CLM:ALM', initialize=initialize_variables, ivar=ivar, index = ir_snow_depth_pa ) + call this%set_restart_var(vname='fates_solar_zenith_angle', vtype=site_r8, & long_name='the angle of the solar zenith', units='radians', flushval = flushinvalid, & hlms='CLM:ALM', initialize=initialize_variables, ivar=ivar, index = ir_solar_zenith_angle ) @@ -2294,7 +2294,7 @@ subroutine set_restart_vectors(this,nc,nsites,sites) rio_min_allowed_landuse_fraction_si => this%rvars(ir_min_allowed_landuse_fraction_si)%r81d, & rio_landuse_vector_gt_min_si => this%rvars(ir_landuse_vector_gt_min_si)%int1d, & rio_area_bareground_si => this%rvars(ir_area_bareground_si)%r81d, & - rio_snow_depth_si => this%rvars(ir_snow_depth_si)%r81d, & + rio_snow_depth_pa => this%rvars(ir_snow_depth_pa)%r81d, & rio_landuse_config_s => this%rvars(ir_landuse_config_si)%int1d, & rio_ncohort_pa => this%rvars(ir_ncohort_pa)%int1d, & rio_fcansno_pa => this%rvars(ir_fcansno_pa)%r81d, & @@ -2803,6 +2803,8 @@ subroutine set_restart_vectors(this,nc,nsites,sites) rio_fcansno_pa( io_idx_co_1st ) = cpatch%fcansno + rio_snow_depth_pa( io_idx_co_1st ) = cpatch%snow_depth + if ( debug ) then write(fates_log(),*) 'offsetNumCohorts III ' & ,io_idx_co,cohortsperpatch @@ -2983,7 +2985,6 @@ subroutine set_restart_vectors(this,nc,nsites,sites) rio_solar_zenith_angle(io_idx_si) = sites(s)%coszen rio_fireweather_index_si(io_idx_si) = sites(s)%fireWeather%fire_weather_index - rio_snow_depth_si(io_idx_si) = sites(s)%snow_depth ! land use flag rio_landuse_config_si(io_idx_si) = hlm_use_potentialveg @@ -3355,7 +3356,7 @@ subroutine get_restart_vectors(this, nc, nsites, sites) rio_min_allowed_landuse_fraction_si => this%rvars(ir_min_allowed_landuse_fraction_si)%r81d, & rio_landuse_vector_gt_min_si => this%rvars(ir_landuse_vector_gt_min_si)%int1d, & rio_area_bareground_si => this%rvars(ir_area_bareground_si)%r81d, & - rio_snow_depth_si => this%rvars(ir_snow_depth_si)%r81d, & + rio_snow_depth_pa => this%rvars(ir_snow_depth_pa)%r81d, & rio_landuse_config_si => this%rvars(ir_landuse_config_si)%int1d, & rio_ncohort_pa => this%rvars(ir_ncohort_pa)%int1d, & rio_fcansno_pa => this%rvars(ir_fcansno_pa)%r81d, & @@ -3820,6 +3821,7 @@ subroutine get_restart_vectors(this, nc, nsites, sites) cpatch%area = rio_area_pa(io_idx_co_1st) cpatch%age_class = get_age_class_index(cpatch%age) cpatch%fcansno = rio_fcansno_pa(io_idx_co_1st) + cpatch%snow_depth = rio_snow_depth_pa(io_idx_co_1st) ! Set zenith angle info @@ -4075,7 +4077,6 @@ subroutine get_restart_vectors(this, nc, nsites, sites) sites(s)%coszen = rio_solar_zenith_angle(io_idx_si) sites(s)%fireWeather%fire_weather_index = rio_fireweather_index_si(io_idx_si) - sites(s)%snow_depth = rio_snow_depth_si(io_idx_si) ! if needed, trigger the special procedure to initialize land use structure from a ! restart run that did not include land use. diff --git a/radiation/FatesRadiationDriveMod.F90 b/radiation/FatesRadiationDriveMod.F90 index a401388f7e..ecbd0219bd 100644 --- a/radiation/FatesRadiationDriveMod.F90 +++ b/radiation/FatesRadiationDriveMod.F90 @@ -185,7 +185,7 @@ subroutine FatesNormalizedCanopyRadiation(sites, bc_in, bc_out ) if(debug) then currentPatch%twostr%band(ib)%Rbeam_atm = 1._r8 currentPatch%twostr%band(ib)%Rdiff_atm = 1._r8 - call CheckPatchRadiationBalance(currentPatch, sites(s)%snow_depth, & + call CheckPatchRadiationBalance(currentPatch, currentPatch%snow_depth, & ib, bc_out(s)%fabd_parb(ifp,ib),bc_out(s)%fabi_parb(ifp,ib)) currentPatch%twostr%band(ib)%Rbeam_atm = fates_unset_r8 currentPatch%twostr%band(ib)%Rdiff_atm = fates_unset_r8 diff --git a/radiation/FatesTwoStreamUtilsMod.F90 b/radiation/FatesTwoStreamUtilsMod.F90 index f8d496c8b1..166c4e9083 100644 --- a/radiation/FatesTwoStreamUtilsMod.F90 +++ b/radiation/FatesTwoStreamUtilsMod.F90 @@ -225,7 +225,7 @@ subroutine FatesConstructRadElements(site) 0, & cohort%nv, & cohort%pft, & - site%snow_depth, & + patch%snow_depth, & vai_top, vai_bot, & elai_cohort,esai_cohort)