From 63359083fdfef00da3d38b5b7381f45a38f57783 Mon Sep 17 00:00:00 2001 From: ckoven Date: Fri, 13 Nov 2020 15:45:51 -0700 Subject: [PATCH 01/27] fixed error in fire intensity calculation --- fire/SFMainMod.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index dd417c9974..c48543b861 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -800,8 +800,8 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) W = currentPatch%TFC_ROS / 0.45_r8 !kgC/m2 to kgbiomass/m2 ! EQ 15 Thonicke et al 2010 - !units of fire intensity = (kJ/kg)*(kgBiomass/m2)*(m/min)*unitless_fraction - currentPatch%FI = SF_val_fuel_energy * W * ROS * currentPatch%frac_burnt !kj/m/s, or kW/m + !units of fire intensity = (kJ/kg)*(kgBiomass/m2)*(m/min) + currentPatch%FI = SF_val_fuel_energy * W * ROS !kj/m/s, or kW/m if(write_sf == itrue)then if( hlm_masterproc == itrue ) write(fates_log(),*) 'fire_intensity',currentPatch%fi,W,currentPatch%ROS_front From 41f5dea20812a971769e0cb238bb1da31ef9ed80 Mon Sep 17 00:00:00 2001 From: ckoven Date: Fri, 13 Nov 2020 16:09:55 -0700 Subject: [PATCH 02/27] added diagnostic to track total number of succesful iginitions, and added that to history --- fire/SFMainMod.F90 | 6 ++++-- main/EDInitMod.F90 | 2 ++ main/EDTypesMod.F90 | 1 + main/FatesHistoryInterfaceMod.F90 | 4 ++-- 4 files changed, 9 insertions(+), 4 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index c48543b861..8591b31b89 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -696,7 +696,7 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) ! ---initialize site parameters to zero--- currentSite%frac_burnt = 0.0_r8 - + currentSite%NF_successful = 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 @@ -810,7 +810,9 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) !'decide_fire' subroutine if (currentPatch%FI > SF_val_fire_threshold) then !track fires greater than kW/m energy threshold currentPatch%fire = 1 ! Fire... :D - + ! + currentSite%NF_successful = currentSite%NF_successful + & + currentSite%NF * currentSite%FDI * currentPatch%area / area else currentPatch%fire = 0 ! No fire... :-/ currentPatch%FD = 0.0_r8 diff --git a/main/EDInitMod.F90 b/main/EDInitMod.F90 index bec7a99537..caf93f8a7f 100644 --- a/main/EDInitMod.F90 +++ b/main/EDInitMod.F90 @@ -186,6 +186,7 @@ subroutine zero_site( site_in ) ! FIRE site_in%acc_ni = 0.0_r8 ! daily nesterov index accumulating over time. time unlimited theoretically. site_in%NF = 0.0_r8 ! daily lightning strikes per km2 + site_in%NF_successful = 0.0_r8 ! daily successful iginitions per km2 site_in%frac_burnt = 0.0_r8 ! burn area read in from external file do el=1,num_elements @@ -296,6 +297,7 @@ subroutine set_site_properties( nsites, sites,bc_in ) sites(s)%acc_NI = acc_NI sites(s)%NF = 0.0_r8 + sites(s)%NF_successful = 0.0_r8 sites(s)%frac_burnt = 0.0_r8 ! PLACEHOLDER FOR PFT AREA DATA MOVED ACROSS INTERFACE diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index 7288a4ef99..9a5fa28e6c 100644 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -730,6 +730,7 @@ module EDTypesMod real(r8) :: acc_ni ! daily nesterov index accumulating over time. real(r8) :: fdi ! daily probability an ignition event will start a fire real(r8) :: NF ! daily ignitions in km2 + real(r8) :: NF_successful ! daily ignitions in km2 that actually lead to fire real(r8) :: frac_burnt ! fraction of area burnt in this day. ! PLANT HYDRAULICS diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index ec02ecbfc1..f1789495af 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -2071,7 +2071,7 @@ subroutine update_history_dyn(this,nc,nsites,sites) ! site-level fire variables hio_nesterov_fire_danger_si(io_si) = sites(s)%acc_NI - hio_fire_nignitions_si(io_si) = sites(s)%NF + hio_fire_nignitions_si(io_si) = sites(s)%NF_successful hio_fire_fdi_si(io_si) = sites(s)%FDI ! If hydraulics are turned on, track the error terms @@ -4362,7 +4362,7 @@ subroutine define_history_vars(this, initialize_variables) ivar=ivar, initialize=initialize_variables, index = ih_nesterov_fire_danger_si) call this%set_history_var(vname='FIRE_IGNITIONS', units='number/km2/day', & - long='number of ignitions', use_default='active', & + long='number of successful ignitions', use_default='active', & avgflag='A', vtype=site_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & ivar=ivar, initialize=initialize_variables, index = ih_fire_nignitions_si) From 055a370358b8fb84ff18cff18227a4506c65d218 Mon Sep 17 00:00:00 2001 From: ckoven Date: Tue, 17 Nov 2020 10:24:09 -0700 Subject: [PATCH 03/27] reverted to ld frag scalar --- biogeochem/EDPhysiologyMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/biogeochem/EDPhysiologyMod.F90 b/biogeochem/EDPhysiologyMod.F90 index 26073ecd35..6d4aa312b8 100644 --- a/biogeochem/EDPhysiologyMod.F90 +++ b/biogeochem/EDPhysiologyMod.F90 @@ -2248,7 +2248,7 @@ subroutine fragmentation_scaler( currentPatch, bc_in) ! ! !LOCAL VARIABLES: logical :: use_century_tfunc = .false. - logical :: use_hlm_soil_scalar = .true. ! Use hlm input decomp fraction scalars + logical :: use_hlm_soil_scalar = .false. ! Use hlm input decomp fraction scalars integer :: j integer :: ifp ! Index of a FATES Patch "ifp" real(r8) :: t_scalar ! temperature scalar From 799dc650b7ef8e546e4610a064470793813f3512 Mon Sep 17 00:00:00 2001 From: ckoven Date: Sun, 22 Nov 2020 16:11:11 -0700 Subject: [PATCH 04/27] fixed frac_burnt unit error and added burned-area hist var to check order of operations --- fire/SFMainMod.F90 | 12 ++++++------ main/EDInitMod.F90 | 5 ++++- main/EDTypesMod.F90 | 4 ++-- main/FatesHistoryInterfaceMod.F90 | 11 +++++++++++ 4 files changed, 23 insertions(+), 9 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 8591b31b89..253559f7df 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -695,7 +695,7 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) real(r8),parameter :: km2_to_m2 = 1000000.0_r8 !area conversion for square km to square m ! ---initialize site parameters to zero--- - currentSite%frac_burnt = 0.0_r8 + currentSite%frac_burnt(:) = 0.0_r8 currentSite%NF_successful = 0._r8 ! Equation 7 from Venevsky et al GCB 2002 (modification of equation 8 in Thonicke et al. 2010) @@ -788,7 +788,7 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) AB = size_of_fire * currentSite%NF * currentSite%FDI !frac_burnt - currentPatch%frac_burnt = (min(0.99_r8, AB / km2_to_m2)) * currentPatch%area/area + currentPatch%frac_burnt = (min(0.99_r8, AB / km2_to_m2)) if(write_SF == itrue)then if ( hlm_masterproc == itrue ) write(fates_log(),*) 'frac_burnt',currentPatch%frac_burnt @@ -813,6 +813,10 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) ! currentSite%NF_successful = currentSite%NF_successful + & currentSite%NF * currentSite%FDI * currentPatch%area / area + ! + ! accumulate frac_burnt % at site level + currentSite%frac_burnt(cpatch%age_class) = currentSite%frac_burnt(cpatch%age_class) + currentPatch%frac_burnt + ! else currentPatch%fire = 0 ! No fire... :-/ currentPatch%FD = 0.0_r8 @@ -820,11 +824,7 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) endif endif! NF ignitions check - - ! accumulate frac_burnt % at site level - currentSite%frac_burnt = currentSite%frac_burnt + currentPatch%frac_burnt - currentPatch => currentPatch%younger enddo !end patch loop diff --git a/main/EDInitMod.F90 b/main/EDInitMod.F90 index caf93f8a7f..d59c5ffda0 100644 --- a/main/EDInitMod.F90 +++ b/main/EDInitMod.F90 @@ -44,6 +44,7 @@ module EDInitMod use FatesInterfaceTypesMod , only : nleafage use FatesInterfaceTypesMod , only : nlevsclass use FatesInterfaceTypesMod , only : nlevcoage + use FatesInterfaceTypesMod , only : nlevage use FatesAllometryMod , only : h2d_allom use FatesAllometryMod , only : bagw_allom use FatesAllometryMod , only : bbgw_allom @@ -119,6 +120,8 @@ subroutine init_site_vars( site_in, bc_in ) allocate(site_in%mass_balance(1:num_elements)) allocate(site_in%flux_diags(1:num_elements)) + allocate(site_in%frac_burnt(1:nlevage)) + site_in%nlevsoil = bc_in%nlevsoil allocate(site_in%rootfrac_scr(site_in%nlevsoil)) allocate(site_in%zi_soil(0:site_in%nlevsoil)) @@ -187,7 +190,7 @@ subroutine zero_site( site_in ) site_in%acc_ni = 0.0_r8 ! daily nesterov index accumulating over time. time unlimited theoretically. site_in%NF = 0.0_r8 ! daily lightning strikes per km2 site_in%NF_successful = 0.0_r8 ! daily successful iginitions per km2 - site_in%frac_burnt = 0.0_r8 ! burn area read in from external file + site_in%frac_burnt(:) = 0.0_r8 ! burn area do el=1,num_elements ! Zero the state variables used for checking mass conservation diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index 9a5fa28e6c..bca591f303 100644 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -555,7 +555,7 @@ module EDTypesMod ! FIRE EFFECTS real(r8) :: scorch_ht(maxpft) ! scorch height: m - real(r8) :: frac_burnt ! fraction burnt: frac gridcell/day + real(r8) :: frac_burnt ! fraction burnt: frac patch/day real(r8) :: tfc_ros ! total fuel consumed - no trunks. KgC/m2/day real(r8) :: burnt_frac_litter(nfsc) ! fraction of each litter pool burned:- @@ -731,7 +731,7 @@ module EDTypesMod real(r8) :: fdi ! daily probability an ignition event will start a fire real(r8) :: NF ! daily ignitions in km2 real(r8) :: NF_successful ! daily ignitions in km2 that actually lead to fire - real(r8) :: frac_burnt ! fraction of area burnt in this day. + real(r8), allocatable :: frac_burnt(:) ! fraction of gridcell area burnt in this day. indexed by patch age bins ! PLANT HYDRAULICS type(ed_site_hydr_type), pointer :: si_hydr diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index dd8d1287a2..f0cbaddef3 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -522,6 +522,7 @@ module FatesHistoryInterfaceMod integer :: ih_agesince_anthrodist_si_age integer :: ih_secondaryforest_area_si_age integer :: ih_area_burnt_si_age + integer :: ih_area_burnt2_si_age ! integer :: ih_fire_rate_of_spread_front_si_age integer :: ih_fire_intensity_si_age integer :: ih_fire_sum_fuel_si_age @@ -2001,6 +2002,7 @@ subroutine update_history_dyn(this,nc,nsites,sites) hio_agesince_anthrodist_si_age => this%hvars(ih_agesince_anthrodist_si_age)%r82d, & hio_secondaryforest_area_si_age => this%hvars(ih_secondaryforest_area_si_age)%r82d, & hio_area_burnt_si_age => this%hvars(ih_area_burnt_si_age)%r82d, & + hio_area_burnt2_si_age => this%hvars(ih_area_burnt2_si_age)%r82d, & ! hio_fire_rate_of_spread_front_si_age => this%hvars(ih_fire_rate_of_spread_front_si_age)%r82d, & hio_fire_intensity_si_age => this%hvars(ih_fire_intensity_si_age)%r82d, & hio_fire_sum_fuel_si_age => this%hvars(ih_fire_sum_fuel_si_age)%r82d, & @@ -2873,6 +2875,9 @@ subroutine update_history_dyn(this,nc,nsites,sites) hio_lai_si_age(io_si, ipa2) = 0._r8 hio_ncl_si_age(io_si, ipa2) = 0._r8 endif + + hio_area_burnt2_si_age(io_si,ipa2) = sites(s)%frac_burnt(ipa2) + end do ! pass the cohort termination mortality as a flux to the history, and then reset the termination mortality buffer @@ -4509,6 +4514,12 @@ subroutine define_history_vars(this, initialize_variables) avgflag='A', vtype=site_age_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & ivar=ivar, initialize=initialize_variables, index = ih_area_burnt_si_age ) + call this%set_history_var(vname='AREA_BURNT2_BY_PATCH_AGE', units='m2/m2', & + long='spitfire area burnt by patch age (divide by patch_area_by_age to get burnt fraction by age)', & + 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_burnt2_si_age ) + call this%set_history_var(vname='FIRE_INTENSITY_BY_PATCH_AGE', units='kJ/m/2', & long='product of fire intensity and burned area, resolved by patch age (so divide by AREA_BURNT_BY_PATCH_AGE to get burned-area-weighted-average intensity', & use_default='active', & From d793fe3a1b1a0e2f00543180bedc1f73969201e2 Mon Sep 17 00:00:00 2001 From: ckoven Date: Sun, 22 Nov 2020 16:39:23 -0700 Subject: [PATCH 05/27] bugfix on prior --- fire/SFMainMod.F90 | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 253559f7df..518aabb9a9 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -815,7 +815,8 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) currentSite%NF * currentSite%FDI * currentPatch%area / area ! ! accumulate frac_burnt % at site level - currentSite%frac_burnt(cpatch%age_class) = currentSite%frac_burnt(cpatch%age_class) + currentPatch%frac_burnt + currentSite%frac_burnt(currentPatch%age_class) = currentSite%frac_burnt(currentPatch%age_class) + & + currentPatch%frac_burnt ! else currentPatch%fire = 0 ! No fire... :-/ From 081bab968c2e9f05b9a2f84a690ab9a129143cb1 Mon Sep 17 00:00:00 2001 From: ckoven Date: Sun, 22 Nov 2020 16:45:45 -0700 Subject: [PATCH 06/27] unit fix on prior --- fire/SFMainMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 518aabb9a9..9cc7162464 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -816,7 +816,7 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) ! ! accumulate frac_burnt % at site level currentSite%frac_burnt(currentPatch%age_class) = currentSite%frac_burnt(currentPatch%age_class) + & - currentPatch%frac_burnt + currentPatch%frac_burnt * currentPatch%area / area ! else currentPatch%fire = 0 ! No fire... :-/ From cfc76321d5929a66a3410716d6441ba9ed691a78 Mon Sep 17 00:00:00 2001 From: ckoven Date: Mon, 23 Nov 2020 10:17:27 -0700 Subject: [PATCH 07/27] fix to vector-zeroing of site_in%frac_burnt in set_site_properties --- main/EDInitMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/main/EDInitMod.F90 b/main/EDInitMod.F90 index d59c5ffda0..bb09d13704 100644 --- a/main/EDInitMod.F90 +++ b/main/EDInitMod.F90 @@ -301,7 +301,7 @@ subroutine set_site_properties( nsites, sites,bc_in ) sites(s)%acc_NI = acc_NI sites(s)%NF = 0.0_r8 sites(s)%NF_successful = 0.0_r8 - sites(s)%frac_burnt = 0.0_r8 + sites(s)%frac_burnt(:) = 0.0_r8 ! PLACEHOLDER FOR PFT AREA DATA MOVED ACROSS INTERFACE if(hlm_use_fixed_biogeog.eq.itrue)then From 4f3bd4eb75e0f9fc39c0ca776d0c9c0b1ebb8490 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Tue, 24 Nov 2020 17:32:32 -0700 Subject: [PATCH 08/27] Update lb term for fire ellipse shape --- fire/SFMainMod.F90 | 36 ++++++++++++++++++++++-------------- 1 file changed, 22 insertions(+), 14 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index dd417c9974..4e3bbe0244 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -682,6 +682,7 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) real(r8) ROS !m/s real(r8) W !kgBiomass/m2 + real(r8) :: tree_fraction_patch ! patch level. no units real(r8) lb !length to breadth ratio of fire ellipse (unitless) real(r8) df !distance fire has travelled forward in m real(r8) db !distance fire has travelled backward in m @@ -692,7 +693,8 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) real(r8) anthro_ign_count ! anthropogenic ignition count/km2/day integer :: iofp ! index of oldest fates patch real(r8), parameter :: pot_hmn_ign_counts_alpha = 0.0035_r8 ! Potential human ignition counts (alpha in Li et al. 2012) (#/person/month) - real(r8),parameter :: km2_to_m2 = 1000000.0_r8 !area conversion for square km to square m + real(r8), parameter :: km2_to_m2 = 1000000.0_r8 !area conversion for square km to square m + real(r8), parameter :: wind_convert = 0.06_r8 ! convert wind speed from m/min to km/hr ! ---initialize site parameters to zero--- currentSite%frac_burnt = 0.0_r8 @@ -751,19 +753,25 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) !equation 15 in Arora and Boer CTEM model.Average fire is 1 day long. !currentPatch%FD = 60.0_r8 * 24.0_r8 !no minutes in a day - - ! The feedback between vegetation structure and ellipse size if turned off for now, - ! to reduce the positive feedback in the syste, - ! This will also be investigated by William Hoffmans proposal. - ! if (currentPatch%effect_wspeed < 16.67_r8) then !16.67m/min = 1km/hr - lb = 1.0_r8 - ! else - !FIX(RF,032414) FOR NO GRASS - ! lb = currentPatch%total_canopy_area/currentPatch%area*(1.0_r8)+(8.729_r8 * & - ! ((1.0_r8 -(exp(-0.03_r8 * 0.06_r8 * currentPatch%effect_wspeed)))**2.155_r8)) !& - !& +currentPatch%fpc_grass*(1.1_r8+((0.06_r8*currentPatch%effect_wspeed)**0.0464)) - - ! endif + tree_fraction_patch = 0.0_r8 + tree_fraction_patch = currentPatch%total_tree_area/currentPatch%area + + if(debug)then + write(fates_log(),*) 'SF currentPatch%area ',currentPatch%area + write(fates_log(),*) 'SF currentPatch%total_area ',currentPatch%total_tree_area + write(fates_log(),*) 'SF patch tree fraction ',tree_fraction_patch + write(fates_log(),*) 'SF AREA ',AREA + endif + + if (currentPatch%effect_wspeed < 16.67_r8) then !16.67m/min = 1km/hr + lb = 1.0_r8 + else + if (tree_fraction_patch > 0.55_r8) then !benchmark for forest cover per Staver 2010 + lb = (1.0_r8 + (8.729_r8 * & + ((1.0_r8 -(exp(-0.03_r8 * wind_convert * currentPatch%effect_wspeed)))**2.155_r8)) + else + lb = (1.1_r8+((wind_convert * currentPatch%effect_wspeed)**0.0464)) + endif ! if (lb > 8.0_r8)then ! lb = 8.0_r8 !Constraint Canadian Fire Behaviour System From 587e3aefb56eae41e1e3b25d416fb2c8cf04592d Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Tue, 24 Nov 2020 18:14:27 -0700 Subject: [PATCH 09/27] typo: close equation parenthesis --- fire/SFMainMod.F90 | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 4e3bbe0244..ab6fdf25b1 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -768,9 +768,10 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) else if (tree_fraction_patch > 0.55_r8) then !benchmark for forest cover per Staver 2010 lb = (1.0_r8 + (8.729_r8 * & - ((1.0_r8 -(exp(-0.03_r8 * wind_convert * currentPatch%effect_wspeed)))**2.155_r8)) + ((1.0_r8 -(exp(-0.03_r8 * wind_convert * currentPatch%effect_wspeed)))**2.155_r8))) else - lb = (1.1_r8+((wind_convert * currentPatch%effect_wspeed)**0.0464)) + lb = (1.1_r8+((wind_convert * currentPatch%effect_wspeed)**0.0464)) + endif endif ! if (lb > 8.0_r8)then From b6bdea99c3d1a3024098a9fba9c786c4c85d0ff3 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Wed, 25 Nov 2020 15:38:16 -0700 Subject: [PATCH 10/27] Add ref to original equations --- fire/SFMainMod.F90 | 11 ++++++----- 1 file changed, 6 insertions(+), 5 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index ab6fdf25b1..f043e7d31d 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -766,11 +766,12 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) if (currentPatch%effect_wspeed < 16.67_r8) then !16.67m/min = 1km/hr lb = 1.0_r8 else - if (tree_fraction_patch > 0.55_r8) then !benchmark for forest cover per Staver 2010 - lb = (1.0_r8 + (8.729_r8 * & - ((1.0_r8 -(exp(-0.03_r8 * wind_convert * currentPatch%effect_wspeed)))**2.155_r8))) - else - lb = (1.1_r8+((wind_convert * currentPatch%effect_wspeed)**0.0464)) + if (tree_fraction_patch > 0.55_r8) then !benchmark forest cover, Staver 2010 + ! EQ 79 forest fuels (Canadian Forest Fire Behavior Prediction System Ont.Inf.Rep. ST-X-3, 1992) + lb = (1.0_r8 + (8.729_r8 * & + ((1.0_r8 -(exp(-0.03_r8 * wind_convert * currentPatch%effect_wspeed)))**2.155_r8))) + else ! EQ 80 grass fuels (CFFBPS Ont.Inf.Rep. ST-X-3, 1992) + lb = (1.1_r8+((wind_convert * currentPatch%effect_wspeed)**0.0464)) endif endif From d1e75c925a6363b8ff93773db9eeed214688c80b Mon Sep 17 00:00:00 2001 From: ckoven Date: Tue, 1 Dec 2020 15:28:58 -0700 Subject: [PATCH 11/27] fixing typos in lb calc --- fire/SFMainMod.F90 | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index b81cb58bd1..ef48cf96f7 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -694,7 +694,7 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) integer :: iofp ! index of oldest fates patch real(r8), parameter :: pot_hmn_ign_counts_alpha = 0.0035_r8 ! Potential human ignition counts (alpha in Li et al. 2012) (#/person/month) real(r8), parameter :: km2_to_m2 = 1000000.0_r8 !area conversion for square km to square m - real(r8), parameter :: wind_convert = 0.06_r8 ! convert wind speed from m/min to km/hr + real(r8), parameter :: m_per_min__to__km_per_hour = 0.06_r8 ! convert wind speed from m/min to km/hr ! ---initialize site parameters to zero--- currentSite%frac_burnt(:) = 0.0_r8 @@ -763,15 +763,15 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) write(fates_log(),*) 'SF AREA ',AREA endif - if (currentPatch%effect_wspeed < 16.67_r8) then !16.67m/min = 1km/hr + if ((currentPatch%effect_wspeed*m_per_min__to__km_per_hour) < 1._r8) then !16.67m/min = 1km/hr lb = 1.0_r8 else if (tree_fraction_patch > 0.55_r8) then !benchmark forest cover, Staver 2010 ! EQ 79 forest fuels (Canadian Forest Fire Behavior Prediction System Ont.Inf.Rep. ST-X-3, 1992) lb = (1.0_r8 + (8.729_r8 * & - ((1.0_r8 -(exp(-0.03_r8 * wind_convert * currentPatch%effect_wspeed)))**2.155_r8))) - else ! EQ 80 grass fuels (CFFBPS Ont.Inf.Rep. ST-X-3, 1992) - lb = (1.1_r8+((wind_convert * currentPatch%effect_wspeed)**0.0464)) + ((1.0_r8 -(exp(-0.03_r8 * m_per_min__to__km_per_hour * currentPatch%effect_wspeed)))**2.155_r8))) + else ! EQ 80 grass fuels (CFFBPS Ont.Inf.Rep. ST-X-3, 1992, but with a correction in Information Report GLC-X-10 by Bottom et al., 2009) + lb = (1.1_r8*((m_per_min__to__km_per_hour * currentPatch%effect_wspeed)**0.464_r8)) endif endif From abcf30e7db4fe72be782e4fe2a1fab0c7875b851 Mon Sep 17 00:00:00 2001 From: ckoven Date: Tue, 1 Dec 2020 16:00:04 -0700 Subject: [PATCH 12/27] changing a pi-like number to pi --- fire/SFMainMod.F90 | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index ef48cf96f7..8cef4007f8 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -7,6 +7,7 @@ module SFMainMod use FatesConstantsMod , only : r8 => fates_r8 use FatesConstantsMod , only : itrue, ifalse + use FatesConstantsMod , only : pi_const use FatesInterfaceTypesMod , only : hlm_masterproc ! 1= master process, 0=not master process use EDTypesMod , only : numWaterMem use FatesGlobals , only : fates_log @@ -791,7 +792,7 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) ! AB = AB *3.0_r8 !size of fire = equation 14 Arora and Boer JGR 2005 - size_of_fire = ((3.1416_r8/(4.0_r8*lb))*((df+db)**2.0_r8)) + size_of_fire = ((pi_const/(4.0_r8*lb))*((df+db)**2.0_r8)) !AB = daily area burnt = size fires in m2 * num ignitions per day per km2 * prob ignition starts fire !AB = m2 per km2 per day From 671286f361426e2522a74708b023b3f51338a4d2 Mon Sep 17 00:00:00 2001 From: ckoven Date: Fri, 4 Dec 2020 17:29:15 -0700 Subject: [PATCH 13/27] some code cleanup and documentation --- fire/SFMainMod.F90 | 50 +++++++++++++++++++++++++++++++-------------- main/EDTypesMod.F90 | 4 ++-- 2 files changed, 37 insertions(+), 17 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 8cef4007f8..dccf1a0e62 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -451,10 +451,22 @@ subroutine rate_of_spread ( currentSite ) do while(associated(currentPatch)) - ! ---initialise parameters to zero.--- - beta_ratio = 0.0_r8; q_ig = 0.0_r8; eps = 0.0_r8; a = 0.0_r8; b = 0.0_r8; c = 0.0_r8; e = 0.0_r8 - phi_wind = 0.0_r8; xi = 0.0_r8; reaction_v_max = 0.0_r8; reaction_v_opt = 0.0_r8; mw_weight = 0.0_r8 - moist_damp = 0.0_r8; ir = 0.0_r8; a_beta = 0.0_r8; + ! ---initialise parameters to zero.--- + beta_ratio = 0.0_r8 + q_ig = 0.0_r8 + eps = 0.0_r8 + a = 0.0_r8 + b = 0.0_r8 + c = 0.0_r8 + e = 0.0_r8 + phi_wind = 0.0_r8 + xi = 0.0_r8 + reaction_v_max = 0.0_r8 + reaction_v_opt = 0.0_r8 + mw_weight = 0.0_r8 + moist_damp = 0.0_r8 + ir = 0.0_r8 + a_beta = 0.0_r8 currentPatch%ROS_front = 0.0_r8 ! remove mineral content from net fuel load per Thonicke 2010 for ir calculation @@ -493,11 +505,11 @@ subroutine rate_of_spread ( currentSite ) ! ---effective heating number--- ! Equation A3 in Thonicke et al. 2010. eps = exp(-4.528_r8 / currentPatch%fuel_sav) - ! Equation A7 in Thonicke et al. 2010 + ! Equation A7 in Thonicke et al. 2010 / eqn 49 from Rothermel 1972 b = 0.15988_r8 * (currentPatch%fuel_sav**0.54_r8) - ! Equation A8 in Thonicke et al. 2010 + ! Equation A8 in Thonicke et al. 2010 / eqn 48 from Rothermel 1972 c = 7.47_r8 * (exp(-0.8711_r8 * (currentPatch%fuel_sav**0.55_r8))) - ! Equation A9 in Thonicke et al. 2010. + ! Equation A9 in Thonicke et al. 2010. (which appears to have a typo, using the coefficient from Rothermel 1972 eqn. 50 instead) e = 0.715_r8 * (exp(-0.01094_r8 * currentPatch%fuel_sav)) if (debug) then @@ -587,7 +599,7 @@ subroutine ground_fuel_consumption ( currentSite ) real(r8) :: moist !effective fuel moisture real(r8) :: tau_b(nfsc) !lethal heating rates for each fuel class (min) - real(r8) :: fc_ground(nfsc) !proportion of fuel consumed + real(r8) :: fc_ground(nfsc) !total amount of fuel consumed per area of burned ground (kg C / m2 of burned area) integer :: c @@ -668,7 +680,7 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) !currentSite%FDI probability that an ignition will start a fire !currentSite%NF number of lighting strikes per day per km2 !currentPatch%ROS_front forward ROS (m/min) - !currentPatch%TFC_ROS total fuel consumed by flaming front (kgC/m2) + !currentPatch%TFC_ROS total fuel consumed by flaming front (kgC/m2 of burned area) use FatesInterfaceTypesMod, only : hlm_spitfire_mode use EDParamsMod, only : ED_val_nignitions @@ -791,24 +803,31 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) ! Equation 16 in arora and boer model JGR 2005 ! AB = AB *3.0_r8 - !size of fire = equation 14 Arora and Boer JGR 2005 + !size of fire = equation 14 Arora and Boer JGR 2005 (area of an ellipse) size_of_fire = ((pi_const/(4.0_r8*lb))*((df+db)**2.0_r8)) - !AB = daily area burnt = size fires in m2 * num ignitions per day per km2 * prob ignition starts fire - !AB = m2 per km2 per day + ! AB = daily area burnt = size fires in m2 * num ignitions per day per km2 * prob ignition starts fire + ! AB = m2 per km2 per day + ! the denominator in the units of currentSite%NF is total gridcell area, but since we assume that ignitions + ! are equally probable across patches, currentSite%NF is equivalently per area of a given patch + ! thus AB has units of m2 burned area per km2 patch area per day AB = size_of_fire * currentSite%NF * currentSite%FDI - !frac_burnt + ! frac_burnt + ! just a unit conversion from AB, to become area burned per area patch per day, + ! or just the fraction of the patch burned on that day currentPatch%frac_burnt = (min(0.99_r8, AB / km2_to_m2)) if(write_SF == itrue)then if ( hlm_masterproc == itrue ) write(fates_log(),*) 'frac_burnt',currentPatch%frac_burnt endif + else + currentPatch%frac_burnt = 0._r8 endif ! lb ROS = currentPatch%ROS_front / 60.0_r8 !m/min to m/sec - W = currentPatch%TFC_ROS / 0.45_r8 !kgC/m2 to kgbiomass/m2 + W = currentPatch%TFC_ROS / 0.45_r8 !kgC/m2 of burned area to kgbiomass/m2 of burned area ! EQ 15 Thonicke et al 2010 !units of fire intensity = (kJ/kg)*(kgBiomass/m2)*(m/min) @@ -825,7 +844,8 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) currentSite%NF_successful = currentSite%NF_successful + & currentSite%NF * currentSite%FDI * currentPatch%area / area ! - ! accumulate frac_burnt % at site level + ! accumulate frac_burnt % at site level. this is purely being tracked as a diagnostic + ! since patch%frac_burnt * patch%area changes due to the fire before history output currentSite%frac_burnt(currentPatch%age_class) = currentSite%frac_burnt(currentPatch%age_class) + & currentPatch%frac_burnt * currentPatch%area / area ! diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index bca591f303..9177fda013 100644 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -556,8 +556,8 @@ module EDTypesMod ! FIRE EFFECTS real(r8) :: scorch_ht(maxpft) ! scorch height: m real(r8) :: frac_burnt ! fraction burnt: frac patch/day - real(r8) :: tfc_ros ! total fuel consumed - no trunks. KgC/m2/day - real(r8) :: burnt_frac_litter(nfsc) ! fraction of each litter pool burned:- + real(r8) :: tfc_ros ! total intensity-relevant fuel consumed - no trunks. KgC/m2 of burned ground/day + real(r8) :: burnt_frac_litter(nfsc) ! fraction of each litter pool burned, conditional on it being burned ! PLANT HYDRAULICS (not currently used in hydraulics RGK 03-2018) From b9dd836c2fb3ec692f17046cc4f7321893a0131f Mon Sep 17 00:00:00 2001 From: ckoven Date: Fri, 4 Dec 2020 17:35:32 -0700 Subject: [PATCH 14/27] rerouting seed decay to roots instead of leaves --- biogeochem/EDPhysiologyMod.F90 | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/biogeochem/EDPhysiologyMod.F90 b/biogeochem/EDPhysiologyMod.F90 index 6d4aa312b8..e468365bf9 100644 --- a/biogeochem/EDPhysiologyMod.F90 +++ b/biogeochem/EDPhysiologyMod.F90 @@ -2201,19 +2201,20 @@ subroutine SeedDecayToFines(litt) ! ! !LOCAL VARIABLES: integer :: pft + integer, parameter :: seedlev = 2 - ! Add decaying seeds to the leaf litter + ! Add decaying seeds to a single level of root litter ! ----------------------------------------------------------------------------------- do pft = 1,numpft - litt%leaf_fines_in(ilabile) = litt%leaf_fines_in(ilabile) + & + litt%root_fines_in(ilabile,seedlev) = litt%root_fines_in(ilabile,seedlev) + & (litt%seed_decay(pft) + litt%seed_germ_decay(pft)) * EDPftvarcon_inst%lf_flab(pft) - litt%leaf_fines_in(icellulose) = litt%leaf_fines_in(icellulose) + & + litt%root_fines_in(icellulose,seedlev) = litt%root_fines_in(icellulose,seedlev) + & (litt%seed_decay(pft) + litt%seed_germ_decay(pft)) * EDPftvarcon_inst%lf_fcel(pft) - litt%leaf_fines_in(ilignin) = litt%leaf_fines_in(ilignin) + & + litt%root_fines_in(ilignin,seedlev) = litt%root_fines_in(ilignin,seedlev) + & (litt%seed_decay(pft) + litt%seed_germ_decay(pft)) * EDPftvarcon_inst%lf_flig(pft) enddo From 1766ba4f16a951fd3f8fa076159d33089cb8cd95 Mon Sep 17 00:00:00 2001 From: ckoven Date: Mon, 7 Dec 2020 20:53:06 -0700 Subject: [PATCH 15/27] turned on HLM decomp scalar again --- biogeochem/EDPhysiologyMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/biogeochem/EDPhysiologyMod.F90 b/biogeochem/EDPhysiologyMod.F90 index e468365bf9..0d08c2755e 100644 --- a/biogeochem/EDPhysiologyMod.F90 +++ b/biogeochem/EDPhysiologyMod.F90 @@ -2249,7 +2249,7 @@ subroutine fragmentation_scaler( currentPatch, bc_in) ! ! !LOCAL VARIABLES: logical :: use_century_tfunc = .false. - logical :: use_hlm_soil_scalar = .false. ! Use hlm input decomp fraction scalars + logical :: use_hlm_soil_scalar = .true. ! Use hlm input decomp fraction scalars integer :: j integer :: ifp ! Index of a FATES Patch "ifp" real(r8) :: t_scalar ! temperature scalar From 93da0f6c46ffbe7b8cf6a74a73678393ff152de9 Mon Sep 17 00:00:00 2001 From: ckoven Date: Tue, 8 Dec 2020 16:30:06 -0700 Subject: [PATCH 16/27] rerouting seed_decay flux straight to HLM surface litter pools --- biogeochem/EDPhysiologyMod.F90 | 40 ++---------------------------- biogeochem/FatesSoilBGCFluxMod.F90 | 22 +++++++++++++++- 2 files changed, 23 insertions(+), 39 deletions(-) diff --git a/biogeochem/EDPhysiologyMod.F90 b/biogeochem/EDPhysiologyMod.F90 index 0d08c2755e..32671e082e 100644 --- a/biogeochem/EDPhysiologyMod.F90 +++ b/biogeochem/EDPhysiologyMod.F90 @@ -227,10 +227,6 @@ subroutine PreDisturbanceLitterFluxes( currentSite, currentPatch, bc_in ) ! Calculate loss rate of viable seeds to litter call SeedDecay(litt) - ! Send those decaying seeds in the previous call - ! to the litter input flux - call SeedDecayToFines(litt) - ! Calculate seed germination rate, the status flags prevent ! germination from occuring when the site is in a drought ! (for drought deciduous) or too cold (for cold deciduous) @@ -254,7 +250,8 @@ subroutine PreDisturbanceLitterFluxes( currentSite, currentPatch, bc_in ) ! Fragmentation flux to soil decomposition model [kg/site/day] site_mass%frag_out = site_mass%frag_out + currentPatch%area * & ( sum(litt%ag_cwd_frag) + sum(litt%bg_cwd_frag) + & - sum(litt%leaf_fines_frag) + sum(litt%root_fines_frag)) + sum(litt%leaf_fines_frag) + sum(litt%root_fines_frag) + & + sum(litt%seed_decay) + sum(litt%seed_germ_decay)) end do @@ -2195,39 +2192,6 @@ end subroutine CWDInput ! ===================================================================================== - subroutine SeedDecayToFines(litt) - - type(litter_type) :: litt - ! - ! !LOCAL VARIABLES: - integer :: pft - integer, parameter :: seedlev = 2 - - ! Add decaying seeds to a single level of root litter - ! ----------------------------------------------------------------------------------- - - do pft = 1,numpft - - litt%root_fines_in(ilabile,seedlev) = litt%root_fines_in(ilabile,seedlev) + & - (litt%seed_decay(pft) + litt%seed_germ_decay(pft)) * EDPftvarcon_inst%lf_flab(pft) - - litt%root_fines_in(icellulose,seedlev) = litt%root_fines_in(icellulose,seedlev) + & - (litt%seed_decay(pft) + litt%seed_germ_decay(pft)) * EDPftvarcon_inst%lf_fcel(pft) - - litt%root_fines_in(ilignin,seedlev) = litt%root_fines_in(ilignin,seedlev) + & - (litt%seed_decay(pft) + litt%seed_germ_decay(pft)) * EDPftvarcon_inst%lf_flig(pft) - - enddo - - - return - end subroutine SeedDecayToFines - - - - - - ! ===================================================================================== subroutine fragmentation_scaler( currentPatch, bc_in) ! diff --git a/biogeochem/FatesSoilBGCFluxMod.F90 b/biogeochem/FatesSoilBGCFluxMod.F90 index 258c37e847..a09dc3725b 100644 --- a/biogeochem/FatesSoilBGCFluxMod.F90 +++ b/biogeochem/FatesSoilBGCFluxMod.F90 @@ -772,7 +772,6 @@ subroutine FluxIntoLitterPools(csite, bc_in, bc_out) ! into the soil/decomposition ! layers. It exponentially decays real(r8) :: surface_prof_tot ! normalizes the surface_prof array - integer :: ft ! PFT number integer :: nlev_eff_soil ! number of effective soil layers integer :: nlev_eff_decomp ! number of effective decomp layers real(r8) :: area_frac ! fraction of site's area of current patch @@ -782,6 +781,7 @@ subroutine FluxIntoLitterPools(csite, bc_in, bc_out) integer :: j ! Soil layer index integer :: id ! Decomposition layer index integer :: ic ! CWD type index + integer :: ipft ! PFT index ! NOTE(rgk, 201705) this parameter was brought over from SoilBiogeochemVerticalProfile ! how steep profile is for surface components (1/ e_folding depth) (1/m) @@ -930,6 +930,26 @@ subroutine FluxIntoLitterPools(csite, bc_in, bc_out) end do + + ! decaying seeds from the litter pool + do ipft = 1,numpft + do id = 1,nlev_eff_decomp + + flux_lab_si(id) = flux_lab_si(id) + & + (litt%seed_decay(ipft) + litt%seed_germ_decay(ipft)) * & + EDPftvarcon_inst%lf_flab(ipft) * area_frac* surface_prof(id) + + flux_cel_si(id) = flux_cel_si(id) + & + (litt%seed_decay(ipft) + litt%seed_germ_decay(ipft)) * & + EDPftvarcon_inst%lf_fcel(ipft) * area_frac* surface_prof(id) + + flux_lig_si(id) = flux_lig_si(id) + & + (litt%seed_decay(ipft) + litt%seed_germ_decay(ipft)) * & + EDPftvarcon_inst%lf_flig(ipft) * area_frac* surface_prof(id) + end do + end do + + do j = 1, nlev_eff_soil id = bc_in%decomp_id(j) From 6c02677facd2a77bd134d84336a7a6aaf93b5bcf Mon Sep 17 00:00:00 2001 From: ckoven Date: Tue, 8 Dec 2020 18:11:42 -0700 Subject: [PATCH 17/27] pulled grassland/forest tree coverage threshold as a named constant --- fire/SFMainMod.F90 | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index dccf1a0e62..4696442fd9 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -708,6 +708,7 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) real(r8), parameter :: pot_hmn_ign_counts_alpha = 0.0035_r8 ! Potential human ignition counts (alpha in Li et al. 2012) (#/person/month) real(r8), parameter :: km2_to_m2 = 1000000.0_r8 !area conversion for square km to square m real(r8), parameter :: m_per_min__to__km_per_hour = 0.06_r8 ! convert wind speed from m/min to km/hr + real(r8), parameter :: forest_grassland_lengthtobreadth_threshold = 0.55_r8 ! tree canopy cover below which to use grassland length-to-breadth eqn ! ---initialize site parameters to zero--- currentSite%frac_burnt(:) = 0.0_r8 @@ -779,11 +780,12 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) if ((currentPatch%effect_wspeed*m_per_min__to__km_per_hour) < 1._r8) then !16.67m/min = 1km/hr lb = 1.0_r8 else - if (tree_fraction_patch > 0.55_r8) then !benchmark forest cover, Staver 2010 + if (tree_fraction_patch > forest_grassland_lengthtobreadth_threshold) then !benchmark forest cover, Staver 2010 ! EQ 79 forest fuels (Canadian Forest Fire Behavior Prediction System Ont.Inf.Rep. ST-X-3, 1992) lb = (1.0_r8 + (8.729_r8 * & ((1.0_r8 -(exp(-0.03_r8 * m_per_min__to__km_per_hour * currentPatch%effect_wspeed)))**2.155_r8))) - else ! EQ 80 grass fuels (CFFBPS Ont.Inf.Rep. ST-X-3, 1992, but with a correction in Information Report GLC-X-10 by Bottom et al., 2009) + else ! EQ 80 grass fuels (CFFBPS Ont.Inf.Rep. ST-X-3, 1992, but with a correction from an errata published within + ! Information Report GLC-X-10 by Bottom et al., 2009 because there is a typo in CFFBPS Ont.Inf.Rep. ST-X-3, 1992) lb = (1.1_r8*((m_per_min__to__km_per_hour * currentPatch%effect_wspeed)**0.464_r8)) endif endif From f2ead9f622d88e701ed422be6cfd29d8b2a19908 Mon Sep 17 00:00:00 2001 From: ckoven Date: Wed, 9 Dec 2020 12:05:04 -0700 Subject: [PATCH 18/27] added pft-resolved gpp and npp hist vars, and fixed/improved units on some other hist vars --- main/FatesHistoryInterfaceMod.F90 | 32 +++++++++++++++++++++++++------ 1 file changed, 26 insertions(+), 6 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index f0cbaddef3..d9dab3a6c8 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -506,6 +506,8 @@ module FatesHistoryInterfaceMod integer :: ih_mortality_si_pft integer :: ih_crownarea_si_pft integer :: ih_canopycrownarea_si_pft + integer :: ih_gpp_si_pft + integer :: ih_npp_si_pft ! indices to (site x patch-age) variables integer :: ih_area_si_age @@ -1810,6 +1812,8 @@ subroutine update_history_dyn(this,nc,nsites,sites) hio_mortality_si_pft => this%hvars(ih_mortality_si_pft)%r82d, & hio_crownarea_si_pft => this%hvars(ih_crownarea_si_pft)%r82d, & hio_canopycrownarea_si_pft => this%hvars(ih_canopycrownarea_si_pft)%r82d, & + hio_gpp_si_pft => this%hvars(ih_gpp_si_pft)%r82d, & + hio_npp_si_pft => this%hvars(ih_npp_si_pft)%r82d, & hio_nesterov_fire_danger_si => this%hvars(ih_nesterov_fire_danger_si)%r81d, & hio_fire_nignitions_si => this%hvars(ih_fire_nignitions_si)%r81d, & hio_fire_fdi_si => this%hvars(ih_fire_fdi_si)%r81d, & @@ -2375,14 +2379,20 @@ subroutine update_history_dyn(this,nc,nsites,sites) ! Update PFT crown area hio_crownarea_si_pft(io_si, ft) = hio_crownarea_si_pft(io_si, ft) + & - ccohort%c_area + ccohort%c_area * AREA_INV if (ccohort%canopy_layer .eq. 1) then ! Update PFT canopy crown area hio_canopycrownarea_si_pft(io_si, ft) = hio_canopycrownarea_si_pft(io_si, ft) + & - ccohort%c_area + ccohort%c_area * AREA_INV end if + ! update pft-resolved NPP and GPP fluxes + hio_gpp_si_pft(io_si, ft) = hio_gpp_si_pft(io_si, ft) + & + ccohort%gpp_acc_hold * n_perm2 + + hio_npp_si_pft(io_si, ft) = hio_npp_si_pft(io_si, ft) + & + ccohort%npp_acc_hold * n_perm2 ! Site by Size-Class x PFT (SCPF) @@ -4202,12 +4212,12 @@ subroutine define_history_vars(this, initialize_variables) avgflag='A', vtype=site_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & ivar=ivar, initialize=initialize_variables, index = ih_trimming_si) - call this%set_history_var(vname='AREA_PLANT', units='m2', & + call this%set_history_var(vname='AREA_PLANT', units='m2/m2', & long='area occupied by all plants', use_default='active', & avgflag='A', vtype=site_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & ivar=ivar, initialize=initialize_variables, index = ih_area_plant_si) - call this%set_history_var(vname='AREA_TREES', units='m2', & + call this%set_history_var(vname='AREA_TREES', units='m2/m2', & long='area occupied by woody plants', use_default='active', & avgflag='A', vtype=site_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & ivar=ivar, initialize=initialize_variables, index = ih_area_trees_si) @@ -4293,16 +4303,26 @@ subroutine define_history_vars(this, initialize_variables) 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='PFTcrownarea', units='m2/ha', & + call this%set_history_var(vname='PFTcrownarea', units='m2/m2', & long='total PFT level crown area', use_default='inactive', & avgflag='A', vtype=site_pft_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & ivar=ivar, initialize=initialize_variables, index = ih_crownarea_si_pft ) - call this%set_history_var(vname='PFTcanopycrownarea', units='m2/ha', & + call this%set_history_var(vname='PFTcanopycrownarea', units='m2/m2', & long='total PFT-level canopy-layer crown area', use_default='inactive', & avgflag='A', vtype=site_pft_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & ivar=ivar, initialize=initialize_variables, index = ih_canopycrownarea_si_pft ) + call this%set_history_var(vname='PFTgpp', units='kg C m-2 y-1', & + long='total PFT-level GPP', use_default='active', & + avgflag='A', vtype=site_pft_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & + ivar=ivar, initialize=initialize_variables, index = ih_gpp_si_pft ) + + call this%set_history_var(vname='PFTnpp', units='kg C m-2 y-1', & + long='total PFT-level NPP', use_default='active', & + avgflag='A', vtype=site_pft_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & + ivar=ivar, initialize=initialize_variables, index = ih_npp_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=site_pft_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & From 13fb39824c3c9e3a0eef88953c1aa5c3182b84e6 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Fri, 11 Dec 2020 17:09:57 -0700 Subject: [PATCH 19/27] Comments and refs for MEF equation SPITFIRE --- fire/SFMainMod.F90 | 16 ++++++++++++---- 1 file changed, 12 insertions(+), 4 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 4696442fd9..3ab30e4094 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -244,7 +244,15 @@ subroutine charecteristics_of_fuel ( currentSite ) lg_sf,currentPatch%livegrass,currentPatch%sum_fuel endif - currentPatch%fuel_frac(lg_sf) = currentPatch%livegrass / currentPatch%sum_fuel + currentPatch%fuel_frac(lg_sf) = currentPatch%livegrass / currentPatch%sum_fuel + + ! MEF (moisure of extinction) depends on compactness of fuel, depth, particle size, wind, slope + ! MEF: pine needles=0.30 (text near EQ 28 Rothermal 1972) + ! Table II-1 NFFL mixed fuels models from Rothermal 1983 Gen. Tech. Rep. INT-143 + ! MEF: short grass=0.12,tall grass=0.25,chaparral=0.20,closed timber litter=0.30,hardwood litter=0.25 + ! Thonicke 2010 SAV give MEF:tw=0.45, sb=0.4874, lb=0.5245, tr=0.57, dg=0.404, lg=0.404 + ! no reference for MEF eqn. in Thonicke 2010 + ! Lasslop 2014 Table 1 MEF PFT level:grass=0.2,shrubs=0.3,TropEverGrnTree=0.2,TropDecid Tree=0.3, Extra-trop Tree=0.3 MEF(1:nfsc) = 0.524_r8 - 0.066_r8 * log10(SF_val_SAV(1:nfsc)) !--- weighted average of relative moisture content--- @@ -505,11 +513,11 @@ subroutine rate_of_spread ( currentSite ) ! ---effective heating number--- ! Equation A3 in Thonicke et al. 2010. eps = exp(-4.528_r8 / currentPatch%fuel_sav) - ! Equation A7 in Thonicke et al. 2010 / eqn 49 from Rothermel 1972 + ! Equation A7 in Thonicke et al. 2010 per eqn 49 from Rothermel 1972 b = 0.15988_r8 * (currentPatch%fuel_sav**0.54_r8) - ! Equation A8 in Thonicke et al. 2010 / eqn 48 from Rothermel 1972 + ! Equation A8 in Thonicke et al. 2010 per eqn 48 from Rothermel 1972 c = 7.47_r8 * (exp(-0.8711_r8 * (currentPatch%fuel_sav**0.55_r8))) - ! Equation A9 in Thonicke et al. 2010. (which appears to have a typo, using the coefficient from Rothermel 1972 eqn. 50 instead) + ! Equation A9 in Thonicke et al. 2010. (appears to have typo, using coefficient eqn.50 Rothermel 1972) e = 0.715_r8 * (exp(-0.01094_r8 * currentPatch%fuel_sav)) if (debug) then From 64f9e7b1f2ca86dcc63469831a9428179a330ca9 Mon Sep 17 00:00:00 2001 From: ckoven Date: Fri, 11 Dec 2020 18:10:45 -0700 Subject: [PATCH 20/27] removed site%frac_burnt variable as no tneeded --- fire/SFMainMod.F90 | 7 ------- main/EDInitMod.F90 | 4 ---- main/EDTypesMod.F90 | 1 - main/FatesHistoryInterfaceMod.F90 | 10 ---------- 4 files changed, 22 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 3ab30e4094..dec68185c0 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -719,7 +719,6 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) real(r8), parameter :: forest_grassland_lengthtobreadth_threshold = 0.55_r8 ! tree canopy cover below which to use grassland length-to-breadth eqn ! ---initialize site parameters to zero--- - currentSite%frac_burnt(:) = 0.0_r8 currentSite%NF_successful = 0._r8 ! Equation 7 from Venevsky et al GCB 2002 (modification of equation 8 in Thonicke et al. 2010) @@ -854,15 +853,9 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) currentSite%NF_successful = currentSite%NF_successful + & currentSite%NF * currentSite%FDI * currentPatch%area / area ! - ! accumulate frac_burnt % at site level. this is purely being tracked as a diagnostic - ! since patch%frac_burnt * patch%area changes due to the fire before history output - currentSite%frac_burnt(currentPatch%age_class) = currentSite%frac_burnt(currentPatch%age_class) + & - currentPatch%frac_burnt * currentPatch%area / area - ! else currentPatch%fire = 0 ! No fire... :-/ currentPatch%FD = 0.0_r8 - currentPatch%frac_burnt = 0.0_r8 endif endif! NF ignitions check diff --git a/main/EDInitMod.F90 b/main/EDInitMod.F90 index bb09d13704..27652d9ed9 100644 --- a/main/EDInitMod.F90 +++ b/main/EDInitMod.F90 @@ -120,8 +120,6 @@ subroutine init_site_vars( site_in, bc_in ) allocate(site_in%mass_balance(1:num_elements)) allocate(site_in%flux_diags(1:num_elements)) - allocate(site_in%frac_burnt(1:nlevage)) - site_in%nlevsoil = bc_in%nlevsoil allocate(site_in%rootfrac_scr(site_in%nlevsoil)) allocate(site_in%zi_soil(0:site_in%nlevsoil)) @@ -190,7 +188,6 @@ subroutine zero_site( site_in ) site_in%acc_ni = 0.0_r8 ! daily nesterov index accumulating over time. time unlimited theoretically. site_in%NF = 0.0_r8 ! daily lightning strikes per km2 site_in%NF_successful = 0.0_r8 ! daily successful iginitions per km2 - site_in%frac_burnt(:) = 0.0_r8 ! burn area do el=1,num_elements ! Zero the state variables used for checking mass conservation @@ -301,7 +298,6 @@ subroutine set_site_properties( nsites, sites,bc_in ) sites(s)%acc_NI = acc_NI sites(s)%NF = 0.0_r8 sites(s)%NF_successful = 0.0_r8 - sites(s)%frac_burnt(:) = 0.0_r8 ! PLACEHOLDER FOR PFT AREA DATA MOVED ACROSS INTERFACE if(hlm_use_fixed_biogeog.eq.itrue)then diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index 9177fda013..1a29410fa2 100644 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -731,7 +731,6 @@ module EDTypesMod real(r8) :: fdi ! daily probability an ignition event will start a fire real(r8) :: NF ! daily ignitions in km2 real(r8) :: NF_successful ! daily ignitions in km2 that actually lead to fire - real(r8), allocatable :: frac_burnt(:) ! fraction of gridcell area burnt in this day. indexed by patch age bins ! PLANT HYDRAULICS type(ed_site_hydr_type), pointer :: si_hydr diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index d9dab3a6c8..3dd391b9b1 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -524,7 +524,6 @@ module FatesHistoryInterfaceMod integer :: ih_agesince_anthrodist_si_age integer :: ih_secondaryforest_area_si_age integer :: ih_area_burnt_si_age - integer :: ih_area_burnt2_si_age ! integer :: ih_fire_rate_of_spread_front_si_age integer :: ih_fire_intensity_si_age integer :: ih_fire_sum_fuel_si_age @@ -2006,7 +2005,6 @@ subroutine update_history_dyn(this,nc,nsites,sites) hio_agesince_anthrodist_si_age => this%hvars(ih_agesince_anthrodist_si_age)%r82d, & hio_secondaryforest_area_si_age => this%hvars(ih_secondaryforest_area_si_age)%r82d, & hio_area_burnt_si_age => this%hvars(ih_area_burnt_si_age)%r82d, & - hio_area_burnt2_si_age => this%hvars(ih_area_burnt2_si_age)%r82d, & ! hio_fire_rate_of_spread_front_si_age => this%hvars(ih_fire_rate_of_spread_front_si_age)%r82d, & hio_fire_intensity_si_age => this%hvars(ih_fire_intensity_si_age)%r82d, & hio_fire_sum_fuel_si_age => this%hvars(ih_fire_sum_fuel_si_age)%r82d, & @@ -2886,8 +2884,6 @@ subroutine update_history_dyn(this,nc,nsites,sites) hio_ncl_si_age(io_si, ipa2) = 0._r8 endif - hio_area_burnt2_si_age(io_si,ipa2) = sites(s)%frac_burnt(ipa2) - end do ! pass the cohort termination mortality as a flux to the history, and then reset the termination mortality buffer @@ -4534,12 +4530,6 @@ subroutine define_history_vars(this, initialize_variables) avgflag='A', vtype=site_age_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & ivar=ivar, initialize=initialize_variables, index = ih_area_burnt_si_age ) - call this%set_history_var(vname='AREA_BURNT2_BY_PATCH_AGE', units='m2/m2', & - long='spitfire area burnt by patch age (divide by patch_area_by_age to get burnt fraction by age)', & - 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_burnt2_si_age ) - call this%set_history_var(vname='FIRE_INTENSITY_BY_PATCH_AGE', units='kJ/m/2', & long='product of fire intensity and burned area, resolved by patch age (so divide by AREA_BURNT_BY_PATCH_AGE to get burned-area-weighted-average intensity', & use_default='active', & From 372bfe4d55a113ca4fb75dfd955df28132c726c6 Mon Sep 17 00:00:00 2001 From: ckoven Date: Fri, 11 Dec 2020 18:37:54 -0700 Subject: [PATCH 21/27] fixed units on burned area --- main/FatesHistoryInterfaceMod.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index 3dd391b9b1..ed98301ee4 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -4472,7 +4472,7 @@ subroutine define_history_vars(this, initialize_variables) avgflag='A', vtype=site_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & ivar=ivar, initialize=initialize_variables, index = ih_fire_intensity_area_product_si ) - call this%set_history_var(vname='FIRE_AREA', units='fraction', & + call this%set_history_var(vname='FIRE_AREA', units='fraction/day', & long='spitfire fire area burn fraction', use_default='active', & avgflag='A', vtype=site_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & ivar=ivar, initialize=initialize_variables, index = ih_fire_area_si ) @@ -4524,7 +4524,7 @@ subroutine define_history_vars(this, initialize_variables) avgflag='A', vtype=site_agefuel_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & ivar=ivar, initialize=initialize_variables, index = ih_fuel_amount_age_fuel ) - call this%set_history_var(vname='AREA_BURNT_BY_PATCH_AGE', units='m2/m2', & + call this%set_history_var(vname='AREA_BURNT_BY_PATCH_AGE', units='m2/m2/day', & long='spitfire area burnt by patch age (divide by patch_area_by_age to get burnt fraction by age)', & use_default='active', & avgflag='A', vtype=site_age_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & From e590cd356d34990293334f046a0667809d2e1389 Mon Sep 17 00:00:00 2001 From: ckoven Date: Thu, 17 Dec 2020 11:18:06 -0700 Subject: [PATCH 22/27] updated MEF eqn and description --- fire/SFMainMod.F90 | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index dec68185c0..a508d27283 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -247,13 +247,14 @@ subroutine charecteristics_of_fuel ( currentSite ) currentPatch%fuel_frac(lg_sf) = currentPatch%livegrass / currentPatch%sum_fuel ! MEF (moisure of extinction) depends on compactness of fuel, depth, particle size, wind, slope + ! Eqn here is eqn 27 from Peterson and Ryan (1986) "Modeling Postfire Conifer Mortality for Long-Range Planning" + ! but lots of other approaches in use out there... ! MEF: pine needles=0.30 (text near EQ 28 Rothermal 1972) ! Table II-1 NFFL mixed fuels models from Rothermal 1983 Gen. Tech. Rep. INT-143 ! MEF: short grass=0.12,tall grass=0.25,chaparral=0.20,closed timber litter=0.30,hardwood litter=0.25 - ! Thonicke 2010 SAV give MEF:tw=0.45, sb=0.4874, lb=0.5245, tr=0.57, dg=0.404, lg=0.404 - ! no reference for MEF eqn. in Thonicke 2010 + ! Thonicke 2010 SAV values propagated thru P&R86 eqn below gives MEF:tw=0.355, sb=0.44, lb=0.525, tr=0.63, dg=0.248, lg=0.248 ! Lasslop 2014 Table 1 MEF PFT level:grass=0.2,shrubs=0.3,TropEverGrnTree=0.2,TropDecid Tree=0.3, Extra-trop Tree=0.3 - MEF(1:nfsc) = 0.524_r8 - 0.066_r8 * log10(SF_val_SAV(1:nfsc)) + MEF(1:nfsc) = 0.524_r8 - 0.066_r8 * log(SF_val_SAV(1:nfsc)) !--- weighted average of relative moisture content--- ! Equation 6 in Thonicke et al. 2010. across twig, small branch, large branch, and dead leaves From e01e8731cffc123bca48fd56bab18bdfefd2bd26 Mon Sep 17 00:00:00 2001 From: ckoven Date: Thu, 17 Dec 2020 13:45:58 -0700 Subject: [PATCH 23/27] added some units info to EDTypesMod.F90 and fates_params_default.cdl --- main/EDTypesMod.F90 | 4 ++-- parameter_files/fates_params_default.cdl | 10 +++++----- 2 files changed, 7 insertions(+), 7 deletions(-) diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index 1a29410fa2..110673f94a 100644 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -534,9 +534,9 @@ module EDTypesMod real(r8) :: sum_fuel ! total ground fuel related to ros (omits 1000hr fuels): KgC/m2 real(r8) :: fuel_frac(nfsc) ! fraction of each litter class in the ros_fuel:-. real(r8) :: livegrass ! total aboveground grass biomass in patch. KgC/m2 - real(r8) :: fuel_bulkd ! average fuel bulk density of the ground fuel + real(r8) :: fuel_bulkd ! average fuel bulk density of the ground fuel. kgBiomass/m3 ! (incl. live grasses. omits 1000hr fuels). KgC/m3 - real(r8) :: fuel_sav ! average surface area to volume ratio of the ground fuel + real(r8) :: fuel_sav ! average surface area to volume ratio of the ground fuel. cm-1 ! (incl. live grasses. omits 1000hr fuels). real(r8) :: fuel_mef ! average moisture of extinction factor ! of the ground fuel (incl. live grasses. omits 1000hr fuels). diff --git a/parameter_files/fates_params_default.cdl b/parameter_files/fates_params_default.cdl index 1f813c4b4e..600b97ac00 100644 --- a/parameter_files/fates_params_default.cdl +++ b/parameter_files/fates_params_default.cdl @@ -164,7 +164,7 @@ variables: fates_eca_vmax_ptase:units = "gP/m2/s" ; fates_eca_vmax_ptase:long_name = "maximum production rate for biochemical P (per m2) (ECA)" ; double fates_fire_alpha_SH(fates_pft) ; - fates_fire_alpha_SH:units = "NA" ; + fates_fire_alpha_SH:units = "m / (kw/m)**(2/3)" ; fates_fire_alpha_SH:long_name = "spitfire parameter, alpha scorch height, Equation 16 Thonicke et al 2010" ; double fates_fire_bark_scaler(fates_pft) ; fates_fire_bark_scaler:units = "fraction" ; @@ -501,8 +501,8 @@ variables: fates_z0mr:units = "unitless" ; fates_z0mr:long_name = "Ratio of momentum roughness length to canopy top height" ; double fates_fire_FBD(fates_litterclass) ; - fates_fire_FBD:units = "NA" ; - fates_fire_FBD:long_name = "spitfire parameter related to fuel bulk density, see SFMain.F90" ; + fates_fire_FBD:units = "kg Biomass/m3" ; + fates_fire_FBD:long_name = "fuel bulk density" ; double fates_fire_low_moisture_Coeff(fates_litterclass) ; fates_fire_low_moisture_Coeff:units = "NA" ; fates_fire_low_moisture_Coeff:long_name = "spitfire parameter, equation B1 Thonicke et al 2010" ; @@ -522,8 +522,8 @@ variables: fates_fire_min_moisture:units = "NA" ; fates_fire_min_moisture:long_name = "spitfire litter moisture threshold to be considered very dry" ; double fates_fire_SAV(fates_litterclass) ; - fates_fire_SAV:units = "NA" ; - fates_fire_SAV:long_name = "spitfire parameter related to surface area to volume ratio, see SFMain.F90" ; + fates_fire_SAV:units = "cm-1" ; + fates_fire_SAV:long_name = "fuel surface area to volume ratio" ; double fates_max_decomp(fates_litterclass) ; fates_max_decomp:units = "yr-1" ; fates_max_decomp:long_name = "maximum rate of litter & CWD transfer from non-decomposing class into decomposing class" ; From ea97d55c6de484863e8076ee1867c84d1571a28b Mon Sep 17 00:00:00 2001 From: ckoven Date: Thu, 17 Dec 2020 16:10:53 -0700 Subject: [PATCH 24/27] initing patch fire variables with nans instead of zeros, and removing some unnecessary zeroing --- biogeochem/EDPatchDynamicsMod.F90 | 36 +++++++++++++------------- fire/SFMainMod.F90 | 42 +++---------------------------- 2 files changed, 22 insertions(+), 56 deletions(-) diff --git a/biogeochem/EDPatchDynamicsMod.F90 b/biogeochem/EDPatchDynamicsMod.F90 index ea0b2918db..6f686f418d 100644 --- a/biogeochem/EDPatchDynamicsMod.F90 +++ b/biogeochem/EDPatchDynamicsMod.F90 @@ -2108,31 +2108,31 @@ subroutine zero_patch(cp_p) ! FIRE - currentPatch%litter_moisture(:) = 0.0_r8 ! litter moisture - currentPatch%fuel_eff_moist = 0.0_r8 ! average fuel moisture content of the ground fuel + currentPatch%litter_moisture(:) = nan ! litter moisture + currentPatch%fuel_eff_moist = nan ! average fuel moisture content of the ground fuel ! (incl. live grasses. omits 1000hr fuels) - currentPatch%livegrass = 0.0_r8 ! total ag grass biomass in patch. 1=c3 grass, 2=c4 grass. gc/m2 - currentPatch%sum_fuel = 0.0_r8 ! total ground fuel related to ros (omits 1000hr fuels). gc/m2 - currentPatch%fuel_bulkd = 0.0_r8 ! average fuel bulk density of the ground fuel + currentPatch%livegrass = nan ! total ag grass biomass in patch. 1=c3 grass, 2=c4 grass. gc/m2 + currentPatch%sum_fuel = nan ! total ground fuel related to ros (omits 1000hr fuels). gc/m2 + currentPatch%fuel_bulkd = nan ! average fuel bulk density of the ground fuel ! (incl. live grasses. omits 1000hr fuels). kgc/m3 - currentPatch%fuel_sav = 0.0_r8 ! average surface area to volume ratio of the ground fuel + currentPatch%fuel_sav = nan ! average surface area to volume ratio of the ground fuel ! (incl. live grasses. omits 1000hr fuels). - currentPatch%fuel_mef = 0.0_r8 ! average moisture of extinction factor of the ground fuel + currentPatch%fuel_mef = nan ! average moisture of extinction factor of the ground fuel ! (incl. live grasses. omits 1000hr fuels). - currentPatch%ros_front = 0.0_r8 ! average rate of forward spread of each fire in the patch. m/min. - currentPatch%effect_wspeed = 0.0_r8 ! dailywind modified by fraction of relative grass and tree cover. m/min. - currentPatch%tau_l = 0.0_r8 ! mins p&r(1986) - currentPatch%fuel_frac(:) = 0.0_r8 ! fraction of each litter class in the sum_fuel + currentPatch%ros_front = nan ! average rate of forward spread of each fire in the patch. m/min. + currentPatch%effect_wspeed = nan ! dailywind modified by fraction of relative grass and tree cover. m/min. + currentPatch%tau_l = nan ! mins p&r(1986) + currentPatch%fuel_frac(:) = nan ! fraction of each litter class in the sum_fuel !- for purposes of calculating weighted averages. - currentPatch%tfc_ros = 0.0_r8 ! used in fi calc - currentPatch%fi = 0._r8 ! average fire intensity of flaming front during day. + currentPatch%tfc_ros = nan ! used in fi calc + currentPatch%fi = nan ! average fire intensity of flaming front during day. ! backward ros plays no role. kj/m/s or kw/m. currentPatch%fire = 999 ! sr decide_fire.1=fire hot enough to proceed. 0=stop everything- no fires today - currentPatch%fd = 0.0_r8 ! fire duration (mins) - currentPatch%ros_back = 0.0_r8 ! backward ros (m/min) - currentPatch%scorch_ht(:) = 0.0_r8 ! scorch height of flames on a given PFT - currentPatch%frac_burnt = 0.0_r8 ! fraction burnt daily - currentPatch%burnt_frac_litter(:) = 0.0_r8 + currentPatch%fd = nan ! fire duration (mins) + currentPatch%ros_back = nan ! backward ros (m/min) + currentPatch%scorch_ht(:) = nan ! scorch height of flames on a given PFT + currentPatch%frac_burnt = nan ! fraction burnt daily + currentPatch%burnt_frac_litter(:) = nan currentPatch%btran_ft(:) = 0.0_r8 currentPatch%canopy_layer_tlai(:) = 0.0_r8 diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index a508d27283..01bbc84ac0 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -206,20 +206,11 @@ subroutine charecteristics_of_fuel ( currentSite ) ! NCWD =4 NFSC = 6 ! tw_sf = 1, lb_sf = 3, tr_sf = 4, dl_sf = 5, lg_sf = 6, - ! zero fire arrays. - currentPatch%fuel_eff_moist = 0.0_r8 - currentPatch%fuel_bulkd = 0.0_r8 !this is kgBiomass/m3 for use in rate of spread equations - currentPatch%fuel_sav = 0.0_r8 - currentPatch%fuel_frac(:) = 0.0_r8 - currentPatch%fuel_mef = 0.0_r8 - currentPatch%sum_fuel = 0.0_r8 - currentPatch%fuel_frac = 0.0_r8 if(write_sf == itrue)then if ( hlm_masterproc == itrue ) write(fates_log(),*) ' leaf_litter1 ',sum(litt_c%leaf_fines(:)) if ( hlm_masterproc == itrue ) write(fates_log(),*) ' leaf_litter2 ',sum(litt_c%ag_cwd(:)) if ( hlm_masterproc == itrue ) write(fates_log(),*) ' leaf_litter3 ',currentPatch%livegrass - if ( hlm_masterproc == itrue ) write(fates_log(),*) ' sum fuel', currentPatch%sum_fuel endif currentPatch%sum_fuel = sum(litt_c%leaf_fines(:)) + & @@ -238,8 +229,6 @@ subroutine charecteristics_of_fuel ( currentSite ) currentPatch%fuel_frac(tw_sf:tr_sf) = litt_c%ag_cwd(:) / currentPatch%sum_fuel if(write_sf == itrue)then - if ( hlm_masterproc == itrue ) write(fates_log(),*) 'ff1 ',currentPatch%fuel_frac - if ( hlm_masterproc == itrue ) write(fates_log(),*) 'ff2 ',currentPatch%fuel_frac if ( hlm_masterproc == itrue ) write(fates_log(),*) 'ff2a ', & lg_sf,currentPatch%livegrass,currentPatch%sum_fuel endif @@ -323,7 +312,6 @@ subroutine charecteristics_of_fuel ( currentSite ) currentPatch%fuel_frac(:) = 0.0000000001_r8 currentPatch%fuel_mef = 0.0000000001_r8 currentPatch%sum_fuel = 0.0000000001_r8 - currentPatch%fuel_frac = 0.0000000001_r8 endif ! check values. @@ -460,24 +448,6 @@ subroutine rate_of_spread ( currentSite ) do while(associated(currentPatch)) - ! ---initialise parameters to zero.--- - beta_ratio = 0.0_r8 - q_ig = 0.0_r8 - eps = 0.0_r8 - a = 0.0_r8 - b = 0.0_r8 - c = 0.0_r8 - e = 0.0_r8 - phi_wind = 0.0_r8 - xi = 0.0_r8 - reaction_v_max = 0.0_r8 - reaction_v_opt = 0.0_r8 - mw_weight = 0.0_r8 - moist_damp = 0.0_r8 - ir = 0.0_r8 - a_beta = 0.0_r8 - currentPatch%ROS_front = 0.0_r8 - ! remove mineral content from net fuel load per Thonicke 2010 for ir calculation currentPatch%sum_fuel = currentPatch%sum_fuel * (1.0_r8 - SF_val_miner_total) !net of minerals @@ -563,11 +533,6 @@ subroutine rate_of_spread ( currentSite ) moist_damp = max(0.0_r8,(1.0_r8 - (2.59_r8 * mw_weight) + (5.11_r8 * (mw_weight**2.0_r8)) - & (3.52_r8*(mw_weight**3.0_r8)))) - ! FIX(SPM, 040114) ask RF if this should be an endrun - ! if(write_SF == itrue)then - ! write(fates_log(),*) 'moist_damp' ,moist_damp,mw_weight,currentPatch%fuel_eff_moist,currentPatch%fuel_mef - ! endif - ! ir = reaction intenisty in kJ/m2/min ! currentPatch%sum_fuel converted from kgC/m2 to kgBiomass/m2 for ir calculation ir = reaction_v_opt*(currentPatch%sum_fuel/0.45_r8)*SF_val_fuel_energy*moist_damp*SF_val_miner_damp @@ -615,7 +580,7 @@ subroutine ground_fuel_consumption ( currentSite ) currentPatch => currentSite%oldest_patch; do while(associated(currentPatch)) - currentPatch%burnt_frac_litter = 1.0_r8 + currentPatch%burnt_frac_litter(:) = 1.0_r8 ! Calculate fraction of litter is burnt for all classes. ! Equation B1 in Thonicke et al. 2010--- do c = 1, nfsc !work out the burnt fraction for all pools, even if those pools dont exist. @@ -644,8 +609,9 @@ subroutine ground_fuel_consumption ( currentSite ) ! we can't ever kill -all- of the grass. currentPatch%burnt_frac_litter(lg_sf) = min(0.8_r8,currentPatch%burnt_frac_litter(lg_sf )) + ! reduce burnt amount for mineral content. - currentPatch%burnt_frac_litter = currentPatch%burnt_frac_litter * (1.0_r8-SF_val_miner_total) + currentPatch%burnt_frac_litter(:) = currentPatch%burnt_frac_litter(:) * (1.0_r8-SF_val_miner_total) !---Calculate amount of fuel burnt.--- @@ -758,11 +724,11 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) currentPatch => currentSite%oldest_patch; do while(associated(currentPatch)) ! ---initialize patch parameters to zero--- + currentPatch%FI = 0._r8 currentPatch%fire = 0 currentPatch%FD = 0.0_r8 currentPatch%frac_burnt = 0.0_r8 - if (currentSite%NF > 0.0_r8) then ! Equation 14 in Thonicke et al. 2010 From 32620d77d874357c7302253685d8e48206b4245a Mon Sep 17 00:00:00 2001 From: ckoven Date: Wed, 6 Jan 2021 21:47:43 -0700 Subject: [PATCH 25/27] zeroing patch fire vars for first timestep --- main/EDInitMod.F90 | 30 ++++++++++++++++++++++++++++++ 1 file changed, 30 insertions(+) diff --git a/main/EDInitMod.F90 b/main/EDInitMod.F90 index 27652d9ed9..19a9d2cb00 100644 --- a/main/EDInitMod.F90 +++ b/main/EDInitMod.F90 @@ -357,6 +357,7 @@ subroutine init_patches( nsites, sites, bc_in) type(ed_site_type), pointer :: sitep type(ed_patch_type), pointer :: newp + type(ed_patch_type), pointer :: currentPatch ! List out some nominal patch values that are used for Near Bear Ground initializations ! as well as initializing inventory @@ -437,6 +438,35 @@ subroutine init_patches( nsites, sites, bc_in) end if + ! zero all the patch fire variables for the first timestep + do s = 1, nsites + currentPatch => sites(s)%youngest_patch + do while(associated(currentPatch)) + + currentPatch%litter_moisture(:) = 0._r8 + currentPatch%fuel_eff_moist = 0._r8 + currentPatch%livegrass = 0._r8 + currentPatch%sum_fuel = 0._r8 + currentPatch%fuel_bulkd = 0._r8 + currentPatch%fuel_sav = 0._r8 + currentPatch%fuel_mef = 0._r8 + currentPatch%ros_front = 0._r8 + currentPatch%effect_wspeed = 0._r8 + currentPatch%tau_l = 0._r8 + currentPatch%fuel_frac(:) = 0._r8 + currentPatch%tfc_ros = 0._r8 + currentPatch%fi = 0._r8 + currentPatch%fire = 0 + currentPatch%fd = 0._r8 + currentPatch%ros_back = 0._r8 + currentPatch%scorch_ht(:) = 0._r8 + currentPatch%frac_burnt = 0._r8 + currentPatch%burnt_frac_litter(:) = 0._r8 + + currentPatch => currentPatch%older + enddo + enddo + ! This sets the rhizosphere shells based on the plant initialization ! The initialization of the plant-relevant hydraulics variables ! were set from a call inside of the init_cohorts()->create_cohort() subroutine From cc3b16faf5a05122b5d0dd4876da06630d14d59c Mon Sep 17 00:00:00 2001 From: ckoven Date: Wed, 6 Jan 2021 22:49:10 -0700 Subject: [PATCH 26/27] zeroing patch fire fluxes upon new patch creation --- biogeochem/EDPatchDynamicsMod.F90 | 17 +++++++++++++++++ 1 file changed, 17 insertions(+) diff --git a/biogeochem/EDPatchDynamicsMod.F90 b/biogeochem/EDPatchDynamicsMod.F90 index 6f686f418d..5a49695d15 100644 --- a/biogeochem/EDPatchDynamicsMod.F90 +++ b/biogeochem/EDPatchDynamicsMod.F90 @@ -2034,6 +2034,23 @@ subroutine create_patch(currentSite, new_patch, age, areap, label) new_patch%fabi_sha_z(:,:,:) = 0._r8 new_patch%scorch_ht(:) = 0._r8 new_patch%frac_burnt = 0._r8 + new_patch%litter_moisture(:) = 0._r8 + new_patch%fuel_eff_moist = 0._r8 + new_patch%livegrass = 0._r8 + new_patch%sum_fuel = 0._r8 + new_patch%fuel_bulkd = 0._r8 + new_patch%fuel_sav = 0._r8 + new_patch%fuel_mef = 0._r8 + new_patch%ros_front = 0._r8 + new_patch%effect_wspeed = 0._r8 + new_patch%tau_l = 0._r8 + new_patch%fuel_frac(:) = 0._r8 + new_patch%tfc_ros = 0._r8 + new_patch%fi = 0._r8 + new_patch%fd = 0._r8 + new_patch%ros_back = 0._r8 + new_patch%scorch_ht(:) = 0._r8 + new_patch%burnt_frac_litter(:) = 0._r8 new_patch%total_tree_area = 0.0_r8 new_patch%NCL_p = 1 From c0f3957e58288d27500f7aff909873a1cae2b74d Mon Sep 17 00:00:00 2001 From: ckoven Date: Tue, 12 Jan 2021 20:20:46 -0700 Subject: [PATCH 27/27] re-added patch%frac_burnt=0 when below threshold intensity --- fire/SFMainMod.F90 | 1 + 1 file changed, 1 insertion(+) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 01bbc84ac0..127dfa43f9 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -823,6 +823,7 @@ subroutine area_burnt_intensity ( currentSite, bc_in ) else currentPatch%fire = 0 ! No fire... :-/ currentPatch%FD = 0.0_r8 + currentPatch%frac_burnt = 0.0_r8 endif endif! NF ignitions check