Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Active crown fire development #5

Closed
wants to merge 11 commits into from
Closed

Active crown fire development #5

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
11 commits
Select commit Hold shift + click to select a range
f908f8e
Introducing active crown fire equations
slevis-lmwg Sep 7, 2019
8ac2053
Merge remote-tracking branch 'jkshuman_repo/passive_crown_fire' into …
slevis-lmwg Sep 11, 2019
8eeb2c7
Merge remote-tracking branch 'jkshuman_repo/passive_crown_fire' into …
slevis-lmwg Sep 12, 2019
f00187f
One-line correction in active crown-fire code
slevis-lmwg Sep 16, 2019
12dc816
Merge remote-tracking branch 'jkshuman_repo/passive_crown_fire' into …
slevis-lmwg Oct 3, 2019
ce4933a
Use crown_fire param as real to scale effects of passive crown fires
slevis-lmwg Oct 4, 2019
2917215
Merge branch 'passive_crown_fire' into active_crown_fire
slevis-lmwg Oct 4, 2019
9ba50cb
Merge branch 'passive_crown_fire' into active_crown_fire
slevis-lmwg Oct 4, 2019
c3191bf
Merge remote-tracking branch 'jkshuman_repo/passive_crown_fire' into …
slevis-lmwg Oct 7, 2019
d377f27
Correction that fixes problem of no fire after year 1
slevis-lmwg Oct 8, 2019
402dfd0
Correct crown_fuel_bulkd formula to include n and c_area
slevis-lmwg Oct 9, 2019
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
14 changes: 8 additions & 6 deletions biogeochem/EDCohortDynamicsMod.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -534,11 +534,12 @@ subroutine nan_cohort(cc_p)
currentCohort%ddbhdt = nan ! time derivative of dbh

! FIRE
currentCohort%fraction_crown_burned = nan ! proportion of crown affected by fire
currentCohort%passive_crown_fire_flg = nan ! flag to indicate passive crown fire ignition
currentCohort%cambial_mort = nan ! probability that trees dies due to cambial char P&R (1986)
currentCohort%crownfire_mort = nan ! probability of tree post-fire mortality due to crown scorch
currentCohort%fire_mort = nan ! post-fire mortality from cambial and crown damage assuming two are independent
currentCohort%fraction_crown_burned = nan ! proportion of crown affected by fire
currentCohort%passive_crown_fire_flg = 9999 ! flag to indicate passive crown fire ignition
currentCohort%active_crown_fire_flg = 9999 ! flag to indicate active crown fire ignition
currentCohort%cambial_mort = nan ! probability that trees dies due to cambial char P&R (1986)
currentCohort%crownfire_mort = nan ! probability of tree post-fire mortality due to crown scorch
currentCohort%fire_mort = nan ! post-fire mortality from cambial and crown damage assuming two are independent

end subroutine nan_cohort

Expand Down Expand Up @@ -580,7 +581,8 @@ subroutine zero_cohort(cc_p)

currentcohort%year_net_uptake(:) = 999._r8 ! this needs to be 999, or trimming of new cohorts will break.
currentcohort%ts_net_uptake(:) = 0._r8
currentCohort%size_class = 1
currentcohort%fraction_crown_burned = 0._r8
currentCohort%size_class = 1
currentCohort%seed_prod = 0._r8
currentCohort%size_class_lasttimestep = 0
currentcohort%npp_acc_hold = 0._r8
Expand Down
98 changes: 69 additions & 29 deletions fire/SFMainMod.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -31,7 +31,6 @@ module SFMainMod
use EDtypesMod , only : TR_SF
use FatesLitterMod , only : litter_type

use PRTGenericMod, only : leaf_organ
use PRTGenericMod, only : carbon12_element
use PRTGenericMod, only : all_carbon_elements
use PRTGenericMod, only : leaf_organ
Expand Down Expand Up @@ -114,7 +113,7 @@ subroutine fire_danger_index ( currentSite, bc_in)
!*****************************************************************
! currentSite%acc_NI is the accumulated Nesterov fire danger index

use SFParamsMod, only : SF_val_fdi_a, SF_val_fdi_b
use SFParamsMod, only : SF_val_fdi_a, SF_val_fdi_b
use FatesConstantsMod , only : tfrz => t_water_freeze_k_1atm
use FatesConstantsMod , only : sec_per_day

Expand Down Expand Up @@ -167,7 +166,7 @@ subroutine charecteristics_of_fuel ( currentSite )

type(ed_patch_type), pointer :: currentPatch
type(ed_cohort_type), pointer :: currentCohort
type(litter_type), pointer :: litt_c
type(litter_type), pointer :: litt_c

real(r8) alpha_FMC(nfsc) ! Relative fuel moisture adjusted per drying ratio
real(r8) fuel_moisture(nfsc) ! Scaled moisture content of small litter fuels.
Expand Down Expand Up @@ -439,8 +438,8 @@ subroutine rate_of_spread ( currentSite )
real(r8) a_beta ! dummy variable for product of a* beta_ratio for react_v_opt equation
real(r8) a,b,c,e ! function of fuel sav

logical,parameter :: debug_windspeed = .false. !for debugging
real(r8),parameter :: q_dry = 581.0_r8 !heat of pre-ignition of dry fuels (kJ/kg)
logical,parameter :: debug_windspeed = .false. !for debugging
real(r8),parameter :: q_dry = 581.0_r8 !heat of pre-ignition of dry fuels (kJ/kg)

currentPatch=>currentSite%oldest_patch;

Expand Down Expand Up @@ -573,8 +572,8 @@ subroutine ground_fuel_consumption ( currentSite )
!returns the the hypothetic fuel consumed by the fire

use SFParamsMod, only : SF_val_miner_total, SF_val_min_moisture, &
SF_val_mid_moisture, SF_val_low_moisture_Coeff, SF_val_low_moisture_Slope, &
SF_val_mid_moisture_Coeff, SF_val_mid_moisture_Slope
SF_val_mid_moisture, SF_val_low_moisture_Coeff, SF_val_low_moisture_Slope, &
SF_val_mid_moisture_Coeff, SF_val_mid_moisture_Slope

type(ed_site_type) , intent(in), target :: currentSite
type(ed_patch_type), pointer :: currentPatch
Expand Down Expand Up @@ -663,8 +662,8 @@ subroutine fire_intensity ( currentSite )
!currentPatch%TFC_ROS total fuel consumed by flaming front (kgC/m2)

use FatesInterfaceMod, only : hlm_use_spitfire
use SFParamsMod, only : SF_val_fdi_alpha,SF_val_fuel_energy, &
SF_val_max_durat, SF_val_durat_slope
use SFParamsMod, only : SF_val_fdi_alpha,SF_val_fuel_energy, &
SF_val_max_durat, SF_val_durat_slope

type(ed_site_type), intent(inout), target :: currentSite

Expand All @@ -688,11 +687,11 @@ subroutine fire_intensity ( currentSite )
!'decide_fire' subroutine shortened and put in here...
if (currentPatch%FI >= fire_threshold) then ! 50kW/m is the threshold for a self-sustaining fire
currentPatch%fire = 1 ! Fire... :D

! 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
currentSite%FDI = 1.0_r8 - exp(-SF_val_fdi_alpha*currentSite%acc_NI)

! Equation 14 in Thonicke et al. 2010
! fire duration in minutes

Expand All @@ -706,7 +705,7 @@ subroutine fire_intensity ( currentSite )
!currentPatch%FD = 60.0_r8 * 24.0_r8 !no minutes in a day
else
currentPatch%fire = 0 ! No fire... :-/
currentPatch%FD = 0.0_r8
currentPatch%FD = 0.0_r8
endif
! FIX(SPM,032414) needs a refactor
! FIX(RF,032414) : should happen outside of SF loop - doing all spitfire code is inefficient otherwise.
Expand Down Expand Up @@ -737,7 +736,6 @@ subroutine area_burnt ( currentSite )
real(r8) AB !daily area burnt in m2 per km2
real(r8) size_of_fire !in m2
real(r8),parameter :: km2_to_m2 = 1000000.0_r8 !area conversion for square km to square m
integer g, p

! ---initialize site parameters to zero---
currentSite%frac_burnt = 0.0_r8
Expand Down Expand Up @@ -776,7 +774,7 @@ subroutine area_burnt ( currentSite )

! --- calculate area burnt---
if(lb > 0.0_r8) then

! Equation 1 in Thonicke et al. 2010
! To Do: Connect here with the Li & Levis GDP fire suppression algorithm.
! Equation 16 in arora and boer model JGR 2005
Expand Down Expand Up @@ -815,7 +813,7 @@ subroutine crown_scorching ( currentSite )

type(ed_site_type), intent(in), target :: currentSite

type(ed_patch_type), pointer :: currentPatch
type(ed_patch_type), pointer :: currentPatch
type(ed_cohort_type), pointer :: currentCohort

real(r8) :: f_ag_bmass ! fraction of tree cohort's above-ground biomass as a proportion of total patch ag tree biomass
Expand Down Expand Up @@ -889,14 +887,20 @@ subroutine crown_damage ( currentSite )

!returns the updated currentCohort%fraction_crown_burned for each tree cohort within each patch.
!currentCohort%fraction_crown_burned is the proportion of crown affected by fire



type(ed_site_type), intent(in), target :: currentSite

type(ed_patch_type) , pointer :: currentPatch
type(ed_cohort_type), pointer :: currentCohort

real(r8) :: leaf_c ! leaf carbon [kg]
real(r8) :: crown_fuel_bulkd ! leaf biomass per unit area divided by canopy depth [kg/m3] Van Wagner used this method according to Scott (2006) p. 12
real(r8), parameter :: low_heat_of_combustion = 12700.0_r8 ! [kJ/kg]
real(r8), parameter :: critical_mass_flow_rate = 0.05_r8 ! [kg/m2/s] value for conifer forests; if available for other vegetation, move to the params file?
real(r8) active_crown_FI ! critical fire intensity for active crown fire ignition (kW/m)
real(r8) ignite_active_crown ! ratio for ignition of active crown fire,EQ 14b Bessie & Johnson 1995
real(r8) :: crown_depth ! depth of crown (m)
real(r8) :: height_cbb ! clear branch bole height or crown base height (m)
real(r8) :: passive_crown_FI ! critical fire intensity for passive crown fire ignition (kW/m)
Expand All @@ -912,11 +916,13 @@ subroutine crown_damage ( currentSite )
do while(associated(currentCohort))
currentCohort%fraction_crown_burned = 0.0_r8


if (EDPftvarcon_inst%woody(currentCohort%pft) == 1) then !trees only

! height_cbb = clear branch bole height at base of crown (m)
! inst%crown = crown_depth_frac (PFT)
active_crown_FI = 0.0_r8 !critical fire intensity for active crown fire
ignite_active_crown = 0.0_r8 !ratio for ignition of active crown fire,EQ 14b Bessie & Johnson 1995
currentCohort%active_crown_fire_flg = 0 !flag for active crown fire ignition
crown_depth = currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft)
height_cbb = currentCohort%hite - crown_depth
passive_crown_FI = 0.0_r8
Expand All @@ -925,7 +931,7 @@ subroutine crown_damage ( currentSite )


! Evaluate for passive crown fire ignition
if (EDPftvarcon_inst%crown_fire(currentCohort%pft) > 0.0_r8) then
if (EDPftvarcon_inst%crown_fire(currentCohort%pft) > 0.001_r8) then

! Note: crown_ignition_energy to be calculated based on foliar moisture content from FATES-Hydro
! EQ 3 Van Wagner 1977
Expand All @@ -944,22 +950,56 @@ subroutine crown_damage ( currentSite )

if (ignite_crown >= 1.0_r8) then
currentCohort%passive_crown_fire_flg = 1 ! passive crown fire ignited
! "...in passive crown fires and high intensity surface fires trees can survive. Jack, red
! and white pine survive...scars used in dating fires"
write(fates_log(),*) 'SF currentCohort%passive_crown_fire_flg = ', currentCohort%passive_crown_fire_flg ! slevis diag
! "...in passive crown fires and high intensity surface
! fires trees can survive. Jack, red
! and white pine survive...scars used in dating fires"
! Johnson, E.A. 1992 Fire and Veg Dynamics: North American boreal forest. Cambridge Press
currentCohort%fraction_crown_burned = EDPftvarcon_inst%crown_fire(currentCohort%pft)
! else ! evaluate crown damage based on scorch height
currentCohort%fraction_crown_burned = EDPftvarcon_inst%crown_fire(currentCohort%pft)

! Critical intensity for active crowning (kW/m)
! EQ 12 Bessie and Johnson 1995
! Fuels / 0.45 to get biomass but note that the 0.45
! cancels out and could be removed. Also dividing
! critical_mass_flow_rate by 3.34, an empirical
! constant in Bessie & Johnson 1995
leaf_c = currentCohort%prt%GetState(leaf_organ, all_carbon_elements)
crown_fuel_bulkd = currentCohort%n * leaf_c / &
(0.45_r8 * currentCohort%c_area * crown_depth)
active_crown_FI = critical_mass_flow_rate * &
low_heat_of_combustion * currentPatch%sum_fuel / &
(0.45_r8 * 3.34_r8 * crown_fuel_bulkd)

! Initiate active crown fire?
! EQ 14b Bessie & Johnson 1995
ignite_active_crown = currentPatch%FI / active_crown_FI

if (ignite_active_crown >= 1.0_r8) then
currentCohort%active_crown_fire_flg = 1 ! active crown fire ignited
write(fates_log(),*) 'SF currentCohort%active_crown_fire_flg = ', currentCohort%active_crown_fire_flg ! slevis diag
currentCohort%fraction_crown_burned = 1.0_r8
! TODO Additional effects of active crown fire.
! Currently in code design phase; see
! https://github.com/NGEET/fates/issues/573
else
write(fates_log(),*) 'SF currentPatch%FI < active_crown_FI = ', currentPatch%FI, active_crown_FI ! slevis diag
endif ! ignite active crown fire
else ! crown damage based on scorch height done below
write(fates_log(),*) 'SF currentPatch%FI < passive_crown_FI = ', currentPatch%FI, passive_crown_FI ! slevis diag
endif ! ignite passive crown fire
! else no crown fire today
endif ! crown fire intensity
! else ! not crown fire plant
endif ! evaluate passive crown fire
else ! crown fire not possible
write(fates_log(),*) 'SF currentPatch%FI < crown_fire_threshold = ', currentPatch%FI, crown_fire_threshold ! slevis diag
endif ! fire intensity vs. crown fire threshold
else ! not a crown-fire pft
write(fates_log(),*) 'SF Not a crown-fire pft' ! slevis diag
endif ! crown-fire pft

! For surface fires, are flames in the canopy?
! height_cbb is clear branch bole height or height of bottom of canopy
! Equation 17 in Thonicke et al. 2010
if (currentCohort%hite > 0.0_r8 .and. currentPatch%SH > height_cbb &
.and. currentCohort%passive_crown_fire_flg == 0) then
.and. currentCohort%passive_crown_fire_flg == 0 &
.and. currentCohort%active_crown_fire_flg == 0) then

currentCohort%fraction_crown_burned = min(1.0_r8, ((currentPatch%SH - height_cbb)/crown_depth))
endif !SH frac crown burnt calculation
Expand Down Expand Up @@ -1002,7 +1042,7 @@ subroutine cambial_damage_kill ( currentSite )
if (currentPatch%fire == 1) then
currentCohort => currentPatch%tallest;
do while(associated(currentCohort))
!currentCohort%cambial_mort = 0.0_r8 !testing this removal
currentCohort%cambial_mort = 0.0_r8
if (EDPftvarcon_inst%woody(currentCohort%pft) == 1) then !trees only
! Equation 21 in Thonicke et al 2010
bt = EDPftvarcon_inst%bark_scaler(currentCohort%pft)*currentCohort%dbh ! bark thickness.
Expand Down Expand Up @@ -1059,7 +1099,7 @@ subroutine post_fire_mortality ( currentSite )
! Equation 22 in Thonicke et al. 2010.
currentCohort%crownfire_mort = EDPftvarcon_inst%crown_kill(currentCohort%pft)*currentCohort%fraction_crown_burned**3.0_r8
! Equation 18 in Thonicke et al. 2010.
currentCohort%fire_mort = max(0.0_r8,min(1.0_r8,currentCohort%crownfire_mort+currentCohort%cambial_mort- &
currentCohort%fire_mort = max(0._r8,min(1.0_r8,currentCohort%crownfire_mort+currentCohort%cambial_mort- &
(currentCohort%crownfire_mort*currentCohort%cambial_mort))) !joint prob.
else
currentCohort%fire_mort = 0.0_r8 !Set to zero. Grass mode of death is removal of leaves.
Expand Down
4 changes: 3 additions & 1 deletion main/EDTypesMod.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -344,7 +344,8 @@ module EDTypesMod

! FIRE
real(r8) :: fraction_crown_burned ! proportion of crown affected by fire:-
real(r8) :: passive_crown_fire_flg ! flag for passive crown fire ignition (1=ignition)
integer :: active_crown_fire_flg ! flag for active crown fire ignition
integer :: passive_crown_fire_flg ! flag for passive crown fire ignition (1=ignition)
real(r8) :: cambial_mort ! probability that trees dies due to cambial char
! (conditional on the tree being subjected to the fire)
real(r8) :: crownfire_mort ! probability of tree post-fire mortality
Expand Down Expand Up @@ -1011,6 +1012,7 @@ subroutine dump_cohort(ccohort)
write(fates_log(),*) 'co%ddbhdt = ', ccohort%ddbhdt
write(fates_log(),*) 'co%dbdeaddt = ', ccohort%dbdeaddt
write(fates_log(),*) 'co%fraction_crown_burned = ', ccohort%fraction_crown_burned
write(fates_log(),*) 'co%active_crown_fire_flg = ', ccohort%active_crown_fire_flg
write(fates_log(),*) 'co%passive_crown_fire_flg = ', ccohort%passive_crown_fire_flg
write(fates_log(),*) 'co%fire_mort = ', ccohort%fire_mort
write(fates_log(),*) 'co%crownfire_mort = ', ccohort%crownfire_mort
Expand Down
4 changes: 2 additions & 2 deletions main/FatesHistoryInterfaceMod.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4540,12 +4540,12 @@ subroutine define_history_vars(this, initialize_variables)
upfreq=1, ivar=ivar, initialize=initialize_variables, index = ih_m5_si_scpf )

call this%set_history_var(vname='CROWNFIREMORT_SCPF', units = 'N/ha/yr', &
long='crown fire mortality by pft/size',use_default='inactive', &
long='crown fire mortality by pft/size',use_default='active', &
avgflag='A', vtype=site_size_pft_r8, hlms='CLM:ALM', flushval=0.0_r8, &
upfreq=1, ivar=ivar, initialize=initialize_variables, index = ih_crownfiremort_si_scpf )

call this%set_history_var(vname='CAMBIALFIREMORT_SCPF', units = 'N/ha/yr', &
long='cambial fire mortality by pft/size',use_default='inactive', &
long='cambial fire mortality by pft/size',use_default='active', &
avgflag='A', vtype=site_size_pft_r8, hlms='CLM:ALM', flushval=0.0_r8, &
upfreq=1, ivar=ivar, initialize=initialize_variables, index = ih_cambialfiremort_si_scpf )

Expand Down