From 7193bc611992ab328923b1d953881d95f70ec045 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Wed, 4 Sep 2019 16:06:17 -0600 Subject: [PATCH 01/72] Passive Crown Fire development Add EQ to evaluate potential for crown fire Add EQ for passive crown fire igntion Set fraction of crown burnt based on this conditional Add crown fire threshold, crown fire flag, and crown fire PFT param --- fire/SFMainMod.F90 | 50 +++++++++++++++++++++--- fire/SFParamsMod.F90 | 9 +++++ main/EDPftvarcon.F90 | 14 ++++++- main/EDTypesMod.F90 | 1 + parameter_files/fates_params_default.cdl | 14 ++++++- 5 files changed, 78 insertions(+), 10 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index e7ce073c14..888cf52a8d 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -902,11 +902,17 @@ 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 + + !SF_val_crown_ignition_energy (kJ/kg) for crown foliage with 100% moisture (dry mass) + use SFParamsMod, only : SF_val_crown_ignition_energy !kJ/kg + type(ed_site_type), intent(in), target :: currentSite type(ed_patch_type) , pointer :: currentPatch type(ed_cohort_type), pointer :: currentCohort + real(r8) height_cbb ! lower crown base height (m) or clear branch bole height (m) + currentPatch => currentSite%oldest_patch do while(associated(currentPatch)) @@ -916,21 +922,53 @@ subroutine crown_damage ( currentSite ) do while(associated(currentCohort)) currentCohort%fraction_crown_burned = 0.0_r8 + currentCohort%crown_FI = 0.0_r8 !critical fire intensity for passive crown fire + currentCohort%crown_fire = 0 ! binary flag for passvie crown fire igntion + if (EDPftvarcon_inst%woody(currentCohort%pft) == 1) then !trees only + + ! Calculate clear branch bole height at base of crown + ! inst%crown = crown_depth_frac (PFT) + height_cbb = currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft) + + ! Evaluate for passive crown fire ignition + if (EDPftvarcon_inst%crown_fire(currentCohort%pft) == 1) then + + ! Crown fuel ignition potential, EQ 8 Bessie and Johnson 1995 + currentCohort%crown_FI = (0.01 * height_cbb * SF_val_crown_ignition_energy) + + if (currentCohort%crown_FI >= crown_fire_threshold) then ! 200 kW/m = threshold for crown fire potential + + ! Initiation of passive crown fire, EQ 9 Bessie and Johnson 1995 + currentCohort%passive_crown_fire = currentPatch%FI/currentCohort%crown_FI + if (currentCohort%passive_crown_fire > 1.0_r8) then + currentCohort%crown_fire = 1 ! passive crown fire ignited + currentCohort%fraction_crown_burned = 1.0_r8 + else ! evaluate crown damage based on scorch height + currentCohort%fraction_crown_burned = 0.0_r8 + endif ! passive crown fire + else + currentCohort%crown_fire = 0 ! no crown fire today + currentCohort%fraction_crown_burned = 0.0_r8 + endif ! crown fire intensity + else + currentCohort%crown_fire = 0 ! not crown fire plant + currentCohort%fraction_crown_burned = 0.0_r8 + endif + ! Flames lower than bottom of canopy. - ! c%hite is height of cohort - if (currentPatch%SH < (currentCohort%hite-currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft))) then + ! c%hite is height of cohort + if (currentPatch%SH < (currentCohort%hite- height_cbb)) .and. currentCohort%crown_fire = 0 then currentCohort%fraction_crown_burned = 0.0_r8 else ! Flames part of way up canopy. ! Equation 17 in Thonicke et al. 2010. ! flames over bottom of canopy but not over top. if ((currentCohort%hite > 0.0_r8).and.(currentPatch%SH >= & - (currentCohort%hite-currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft)))) then + (currentCohort%hite - height_cbb))) .and. currentCohort%crown_fire = 0 then - currentCohort%fraction_crown_burned = (currentPatch%SH-currentCohort%hite*(1.0_r8- & - EDPftvarcon_inst%crown(currentCohort%pft)))/(currentCohort%hite* & - EDPftvarcon_inst%crown(currentCohort%pft)) + currentCohort%fraction_crown_burned = (currentPatch%SH-currentCohort%hite*(1.0_r8- & + EDPftvarcon_inst%crown(currentCohort%pft)))/(height_cbb) else ! Flames over top of canopy. diff --git a/fire/SFParamsMod.F90 b/fire/SFParamsMod.F90 index 1c9f278c1e..dddef4d927 100644 --- a/fire/SFParamsMod.F90 +++ b/fire/SFParamsMod.F90 @@ -23,6 +23,7 @@ module SFParamsMod real(r8),protected, public :: SF_val_fdi_alpha real(r8),protected, public :: SF_val_miner_total real(r8),protected, public :: SF_val_fuel_energy + real(r8),protected, public :: SF_val_crown_ignition_energy real(r8),protected, public :: SF_val_part_dens real(r8),protected, public :: SF_val_miner_damp real(r8),protected, public :: SF_val_max_durat @@ -44,6 +45,7 @@ module SFParamsMod character(len=param_string_length),parameter :: SF_name_fdi_alpha = "fates_fire_fdi_alpha" character(len=param_string_length),parameter :: SF_name_miner_total = "fates_fire_miner_total" character(len=param_string_length),parameter :: SF_name_fuel_energy = "fates_fire_fuel_energy" + character(len=param_string_length),parameter :: SF_name_crown_ignition_energy = "fates_fire_crown_ignition_energy" character(len=param_string_length),parameter :: SF_name_part_dens = "fates_fire_part_dens" character(len=param_string_length),parameter :: SF_name_miner_damp = "fates_fire_miner_damp" character(len=param_string_length),parameter :: SF_name_max_durat = "fates_fire_max_durat" @@ -137,6 +139,7 @@ subroutine SpitFireParamsInit() SF_val_fdi_alpha = nan SF_val_miner_total = nan SF_val_fuel_energy = nan + SF_val_crown_ignition_energy = nan SF_val_part_dens = nan SF_val_miner_damp = nan SF_val_max_durat = nan @@ -215,6 +218,9 @@ subroutine SpitFireRegisterScalars(fates_params) call fates_params%RegisterParameter(name=SF_name_fuel_energy, dimension_shape=dimension_shape_scalar, & dimension_names=dim_names_scalar) + call fates_params%RegisterParameter(name=SF_name_crown_ignition_energy, dimension_shape=dimension_shape_scalar, & + dimension_names=dim_names_scalar) + call fates_params%RegisterParameter(name=SF_name_part_dens, dimension_shape=dimension_shape_scalar, & dimension_names=dim_names_scalar) @@ -256,6 +262,9 @@ subroutine SpitFireReceiveScalars(fates_params) call fates_params%RetreiveParameter(name=SF_name_fuel_energy, & data=SF_val_fuel_energy) + call fates_params%RetreiveParameter(name=SF_name_crown_ignition_energy, & + data=SF_val_crown_ignition_energy) + call fates_params%RetreiveParameter(name=SF_name_part_dens, & data=SF_val_part_dens) diff --git a/main/EDPftvarcon.F90 b/main/EDPftvarcon.F90 index 0e4a9f7720..fe06fc91a9 100644 --- a/main/EDPftvarcon.F90 +++ b/main/EDPftvarcon.F90 @@ -48,9 +48,10 @@ module EDPftvarcon real(r8), allocatable :: crown(:) ! fraction of the height of the plant ! that is occupied by crown. For fire model. real(r8), allocatable :: bark_scaler(:) ! scaler from dbh to bark thickness. For fire model. - real(r8), allocatable :: crown_kill(:) ! scaler on fire death. For fire model. + real(r8), allocatable :: crown_kill(:) ! resistance to fire. (1 = none) For fire model. + real(r8), allocatable :: crown_fire(:) ! Does plant have passive crown fire? (1=yes, 0=no) real(r8), allocatable :: initd(:) ! initial seedling density - real(r8), allocatable :: seed_suppl(:) ! seeds that come from outside the gridbox. + real(r8), allocatable :: seed_suppl(:) ! seeds that come from outside the gridbox. real(r8), allocatable :: BB_slope(:) ! ball berry slope parameter real(r8), allocatable :: seed_alloc_mature(:) ! fraction of carbon balance allocated to @@ -381,6 +382,10 @@ subroutine Register_PFT(this, fates_params) call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) + name = 'fates_fire_crown_fire' + call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & + dimension_names=dim_names, lower_bounds=dim_lower_bound) + name = 'fates_recruit_initd' call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) @@ -813,6 +818,10 @@ subroutine Receive_PFT(this, fates_params) call fates_params%RetreiveParameterAllocate(name=name, & data=this%crown_kill) + name = 'fates_fire_crown_fire' + call fates_params%RetreiveParameterAllocate(name=name, & + data=this%crown_fire) + name = 'fates_recruit_initd' call fates_params%RetreiveParameterAllocate(name=name, & data=this%initd) @@ -1722,6 +1731,7 @@ subroutine FatesReportPFTParams(is_master) write(fates_log(),fmt0) 'crown = ',EDPftvarcon_inst%crown write(fates_log(),fmt0) 'bark_scaler = ',EDPftvarcon_inst%bark_scaler write(fates_log(),fmt0) 'crown_kill = ',EDPftvarcon_inst%crown_kill + write(fates_log(),fmt0) 'crown_fire = ',EDPftvarcon_inst%crown_fire write(fates_log(),fmt0) 'initd = ',EDPftvarcon_inst%initd write(fates_log(),fmt0) 'seed_suppl = ',EDPftvarcon_inst%seed_suppl write(fates_log(),fmt0) 'BB_slope = ',EDPftvarcon_inst%BB_slope diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index 1850ea9ce9..4dfcede3d2 100644 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -146,6 +146,7 @@ module EDTypesMod integer, parameter, public :: lg_sf = 6 ! array index of live grass pool for spitfire real(r8), parameter, public :: fire_threshold = 50.0_r8 ! threshold for fires that spread or go out. KWm-2 (Pyne 1986) + real(r8), parameter, public :: crown_fire_threshold = 200.0_r8 ! threshold for passive crown fire ignition. KWm-2 (Bessie & Johnson 1995) ! PATCH FUSION real(r8), parameter, public :: force_patchfuse_min_biomass = 0.005_r8 ! min biomass (kg / m2 patch area) below which to force-fuse patches diff --git a/parameter_files/fates_params_default.cdl b/parameter_files/fates_params_default.cdl index 49b7662508..fa276bf42b 100644 --- a/parameter_files/fates_params_default.cdl +++ b/parameter_files/fates_params_default.cdl @@ -122,7 +122,7 @@ variables: fates_displar:long_name = "Ratio of displacement height to canopy top height" ; double fates_fire_alpha_SH(fates_pft) ; fates_fire_alpha_SH:units = "NA" ; - fates_fire_alpha_SH:long_name = "spitfire parameter, alpha scorch height, Equation 16 Thonicke et al 2010" ; + fates_fire_alpha_SH:long_name = "spitfire parameter, alpha scorch height, EQ 16 Thonicke et al 2010" ; double fates_fire_bark_scaler(fates_pft) ; fates_fire_bark_scaler:units = "fraction" ; fates_fire_bark_scaler:long_name = "the thickness of a cohorts bark as a fraction of its dbh" ; @@ -131,7 +131,10 @@ variables: fates_fire_crown_depth_frac:long_name = "the depth of a cohorts crown as a fraction of its height" ; double fates_fire_crown_kill(fates_pft) ; fates_fire_crown_kill:units = "NA" ; - fates_fire_crown_kill:long_name = "fire parameter, see equation 22 in Thonicke et al 2010" ; + fates_fire_crown_kill:long_name = "resistance to fire (1 = none), EQ 22 in Thonicke et al 2010" ; + double fates_fire_crown_fire(fates_pft) ; + fates_fire_crown_fire:units = "logical flag" ; + fates_fire_crown_fire:long_name = "Binary flag for passive crown fire"; double fates_fr_fcel(fates_pft) ; fates_fr_fcel:units = "fraction" ; fates_fr_fcel:long_name = "Fine root litter cellulose fraction" ; @@ -492,6 +495,9 @@ variables: double fates_fire_fuel_energy ; fates_fire_fuel_energy:units = "kJ/kg" ; fates_fire_fuel_energy:long_name = "spitfire parameter, heat content of fuel" ; + double fates_fire_crown_ignition_energy ; + fates_fire_crown_ignition_energy:units = "kJ/kg" ; + fates_fire_crown_ignition_energy:long_name = "spitfire parameter, heat of crown foliage ignition" ; double fates_fire_max_durat ; fates_fire_max_durat:units = "minutes" ; fates_fire_max_durat:long_name = "spitfire parameter, fire maximum duration, Equation 14 Thonicke et al 2010" ; @@ -718,6 +724,8 @@ data: fates_fire_crown_kill = 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775 ; +fates_fire_passive_crown_fire = 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 ; + fates_fr_fcel = 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5 ; fates_fr_flab = 0.25, 0.25, 0.25, 0.25, 0.25, 0.25, 0.25, 0.25, 0.25, 0.25, @@ -1116,6 +1124,8 @@ data: fates_fire_fuel_energy = 18000 ; + fates_fire_crown_ignition_energy = 3060 ; + fates_fire_max_durat = 240 ; fates_fire_miner_damp = 0.41739 ; From bad24e413d8bd119bbab669d407b72d4270025fc Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Wed, 4 Sep 2019 16:31:22 -0600 Subject: [PATCH 02/72] Correct EQ for crown_FI --- fire/SFMainMod.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 888cf52a8d..11dfe6246f 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -923,7 +923,7 @@ subroutine crown_damage ( currentSite ) do while(associated(currentCohort)) currentCohort%fraction_crown_burned = 0.0_r8 currentCohort%crown_FI = 0.0_r8 !critical fire intensity for passive crown fire - currentCohort%crown_fire = 0 ! binary flag for passvie crown fire igntion + currentCohort%crown_fire = 0 !flag for passvie crown fire ignition if (EDPftvarcon_inst%woody(currentCohort%pft) == 1) then !trees only @@ -935,7 +935,7 @@ subroutine crown_damage ( currentSite ) if (EDPftvarcon_inst%crown_fire(currentCohort%pft) == 1) then ! Crown fuel ignition potential, EQ 8 Bessie and Johnson 1995 - currentCohort%crown_FI = (0.01 * height_cbb * SF_val_crown_ignition_energy) + currentCohort%crown_FI = (0.01 * height_cbb * SF_val_crown_ignition_energy)**1.5 if (currentCohort%crown_FI >= crown_fire_threshold) then ! 200 kW/m = threshold for crown fire potential From 093db93ff2026722d9dfc7ef450602de1b4afeb4 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Thu, 5 Sep 2019 10:39:39 -0600 Subject: [PATCH 03/72] Update and correct variable names for passive crown fire --- fire/SFMainMod.F90 | 6 +++--- main/EDPftvarcon.F90 | 10 +++++----- parameter_files/fates_params_default.cdl | 10 +++++----- 3 files changed, 13 insertions(+), 13 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 11dfe6246f..bc85d80b91 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -940,8 +940,8 @@ subroutine crown_damage ( currentSite ) if (currentCohort%crown_FI >= crown_fire_threshold) then ! 200 kW/m = threshold for crown fire potential ! Initiation of passive crown fire, EQ 9 Bessie and Johnson 1995 - currentCohort%passive_crown_fire = currentPatch%FI/currentCohort%crown_FI - if (currentCohort%passive_crown_fire > 1.0_r8) then + currentCohort%ignite_crown = currentPatch%FI/currentCohort%crown_FI + if (currentCohort%ignite_crown > 1.0_r8) then currentCohort%crown_fire = 1 ! passive crown fire ignited currentCohort%fraction_crown_burned = 1.0_r8 else ! evaluate crown damage based on scorch height @@ -1069,7 +1069,7 @@ subroutine post_fire_mortality ( currentSite ) currentCohort%crownfire_mort = 0.0_r8 if (EDPftvarcon_inst%woody(currentCohort%pft) == 1) then ! Equation 22 in Thonicke et al. 2010. - currentCohort%crownfire_mort = EDPftvarcon_inst%crown_kill(currentCohort%pft)*currentCohort%fraction_crown_burned**3.0_r8 + currentCohort%crownfire_mort = EDPftvarcon_inst%crown_resist(currentCohort%pft)*currentCohort%fraction_crown_burned**3.0_r8 ! Equation 18 in Thonicke et al. 2010. currentCohort%fire_mort = max(0._r8,min(1.0_r8,currentCohort%crownfire_mort+currentCohort%cambial_mort- & (currentCohort%crownfire_mort*currentCohort%cambial_mort))) !joint prob. diff --git a/main/EDPftvarcon.F90 b/main/EDPftvarcon.F90 index fe06fc91a9..7fbda29121 100644 --- a/main/EDPftvarcon.F90 +++ b/main/EDPftvarcon.F90 @@ -48,7 +48,7 @@ module EDPftvarcon real(r8), allocatable :: crown(:) ! fraction of the height of the plant ! that is occupied by crown. For fire model. real(r8), allocatable :: bark_scaler(:) ! scaler from dbh to bark thickness. For fire model. - real(r8), allocatable :: crown_kill(:) ! resistance to fire. (1 = none) For fire model. + real(r8), allocatable :: crown_resist(:) ! resistance to fire. (1 = none) For fire model. real(r8), allocatable :: crown_fire(:) ! Does plant have passive crown fire? (1=yes, 0=no) real(r8), allocatable :: initd(:) ! initial seedling density real(r8), allocatable :: seed_suppl(:) ! seeds that come from outside the gridbox. @@ -378,7 +378,7 @@ subroutine Register_PFT(this, fates_params) call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) - name = 'fates_fire_crown_kill' + name = 'fates_fire_crown_resist' call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) @@ -814,9 +814,9 @@ subroutine Receive_PFT(this, fates_params) call fates_params%RetreiveParameterAllocate(name=name, & data=this%bark_scaler) - name = 'fates_fire_crown_kill' + name = 'fates_fire_crown_resist' call fates_params%RetreiveParameterAllocate(name=name, & - data=this%crown_kill) + data=this%crown_resist) name = 'fates_fire_crown_fire' call fates_params%RetreiveParameterAllocate(name=name, & @@ -1730,7 +1730,7 @@ subroutine FatesReportPFTParams(is_master) write(fates_log(),fmt0) 'leaf_stor_priority = ',EDPftvarcon_inst%leaf_stor_priority write(fates_log(),fmt0) 'crown = ',EDPftvarcon_inst%crown write(fates_log(),fmt0) 'bark_scaler = ',EDPftvarcon_inst%bark_scaler - write(fates_log(),fmt0) 'crown_kill = ',EDPftvarcon_inst%crown_kill + write(fates_log(),fmt0) 'crown_resist = ',EDPftvarcon_inst%crown_resist write(fates_log(),fmt0) 'crown_fire = ',EDPftvarcon_inst%crown_fire write(fates_log(),fmt0) 'initd = ',EDPftvarcon_inst%initd write(fates_log(),fmt0) 'seed_suppl = ',EDPftvarcon_inst%seed_suppl diff --git a/parameter_files/fates_params_default.cdl b/parameter_files/fates_params_default.cdl index fa276bf42b..d98fc2e1ea 100644 --- a/parameter_files/fates_params_default.cdl +++ b/parameter_files/fates_params_default.cdl @@ -129,9 +129,9 @@ variables: double fates_fire_crown_depth_frac(fates_pft) ; fates_fire_crown_depth_frac:units = "fraction" ; fates_fire_crown_depth_frac:long_name = "the depth of a cohorts crown as a fraction of its height" ; - double fates_fire_crown_kill(fates_pft) ; - fates_fire_crown_kill:units = "NA" ; - fates_fire_crown_kill:long_name = "resistance to fire (1 = none), EQ 22 in Thonicke et al 2010" ; + double fates_fire_crown_resist(fates_pft) ; + fates_fire_crown_resist:units = "NA" ; + fates_fire_crown_resist:long_name = "resistance to fire (1 = none), EQ 22 in Thonicke et al 2010" ; double fates_fire_crown_fire(fates_pft) ; fates_fire_crown_fire:units = "logical flag" ; fates_fire_crown_fire:long_name = "Binary flag for passive crown fire"; @@ -721,10 +721,10 @@ data: fates_fire_crown_depth_frac = 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.95, 0.95, 0.95, 1, 1, 1 ; - fates_fire_crown_kill = 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, + fates_fire_crown_resist = 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775 ; -fates_fire_passive_crown_fire = 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 ; +fates_fire_crown_fire = 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 ; fates_fr_fcel = 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5 ; From f3077dfda1ceba060cb8493bf868835d8bfe7cd9 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Thu, 5 Sep 2019 14:28:34 -0600 Subject: [PATCH 04/72] Update name for crown_fire flag --- fire/SFMainMod.F90 | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index bc85d80b91..f8f90402a9 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -923,7 +923,7 @@ subroutine crown_damage ( currentSite ) do while(associated(currentCohort)) currentCohort%fraction_crown_burned = 0.0_r8 currentCohort%crown_FI = 0.0_r8 !critical fire intensity for passive crown fire - currentCohort%crown_fire = 0 !flag for passvie crown fire ignition + currentCohort%crown_fire_flg = 0 !flag for passvie crown fire ignition if (EDPftvarcon_inst%woody(currentCohort%pft) == 1) then !trees only @@ -942,30 +942,30 @@ subroutine crown_damage ( currentSite ) ! Initiation of passive crown fire, EQ 9 Bessie and Johnson 1995 currentCohort%ignite_crown = currentPatch%FI/currentCohort%crown_FI if (currentCohort%ignite_crown > 1.0_r8) then - currentCohort%crown_fire = 1 ! passive crown fire ignited + currentCohort%crown_fire_flg = 1 ! passive crown fire ignited currentCohort%fraction_crown_burned = 1.0_r8 else ! evaluate crown damage based on scorch height currentCohort%fraction_crown_burned = 0.0_r8 endif ! passive crown fire else - currentCohort%crown_fire = 0 ! no crown fire today + currentCohort%crown_fire_flg = 0 ! no crown fire today currentCohort%fraction_crown_burned = 0.0_r8 endif ! crown fire intensity else - currentCohort%crown_fire = 0 ! not crown fire plant + currentCohort%crown_fire_flg = 0 ! not crown fire plant currentCohort%fraction_crown_burned = 0.0_r8 endif ! Flames lower than bottom of canopy. ! c%hite is height of cohort - if (currentPatch%SH < (currentCohort%hite- height_cbb)) .and. currentCohort%crown_fire = 0 then + if (currentPatch%SH < (currentCohort%hite- height_cbb)) .and. currentCohort%crown_fire_flg = 0 then currentCohort%fraction_crown_burned = 0.0_r8 else ! Flames part of way up canopy. ! Equation 17 in Thonicke et al. 2010. ! flames over bottom of canopy but not over top. if ((currentCohort%hite > 0.0_r8).and.(currentPatch%SH >= & - (currentCohort%hite - height_cbb))) .and. currentCohort%crown_fire = 0 then + (currentCohort%hite - height_cbb))) .and. currentCohort%crown_fire_flg = 0 then currentCohort%fraction_crown_burned = (currentPatch%SH-currentCohort%hite*(1.0_r8- & EDPftvarcon_inst%crown(currentCohort%pft)))/(height_cbb) From f28f546ca29580555042788841589a4a1d10c447 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Thu, 5 Sep 2019 14:47:53 -0600 Subject: [PATCH 05/72] Remove redundant details within evaluation for crown fire ignition --- fire/SFMainMod.F90 | 18 +++++++----------- 1 file changed, 7 insertions(+), 11 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index f8f90402a9..17c0aae191 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -939,22 +939,18 @@ subroutine crown_damage ( currentSite ) if (currentCohort%crown_FI >= crown_fire_threshold) then ! 200 kW/m = threshold for crown fire potential - ! Initiation of passive crown fire, EQ 9 Bessie and Johnson 1995 + ! Initiation of passive crown fire, EQ 9 Bessie and Johnson 1995 currentCohort%ignite_crown = currentPatch%FI/currentCohort%crown_FI + if (currentCohort%ignite_crown > 1.0_r8) then currentCohort%crown_fire_flg = 1 ! passive crown fire ignited currentCohort%fraction_crown_burned = 1.0_r8 - else ! evaluate crown damage based on scorch height - currentCohort%fraction_crown_burned = 0.0_r8 - endif ! passive crown fire - else - currentCohort%crown_fire_flg = 0 ! no crown fire today - currentCohort%fraction_crown_burned = 0.0_r8 + ! else ! evaluate crown damage based on scorch height + endif ! ignite passive crown fire + ! else no crown fire today endif ! crown fire intensity - else - currentCohort%crown_fire_flg = 0 ! not crown fire plant - currentCohort%fraction_crown_burned = 0.0_r8 - endif + ! else ! not crown fire plant + endif ! evaluate passive crown fire ! Flames lower than bottom of canopy. ! c%hite is height of cohort From 79a8f5eb9e948c9553652eee3fd01bed65499f1e Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Thu, 5 Sep 2019 16:11:39 -0600 Subject: [PATCH 06/72] Add q_dry as available to entire fire module --- fire/SFMainMod.F90 | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 17c0aae191..60ffde1a51 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -60,6 +60,8 @@ module SFMainMod integer :: write_SF = 0 ! for debugging logical :: debug = .false. ! for debugging + + real(r8),parameter :: q_dry = 581.0_r8 !heat of pre-ignition of dry fuels (kJ/kg) ! ============================================================================ ! ============================================================================ @@ -167,7 +169,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) timeav_swc real(r8) alpha_FMC(nfsc) ! Relative fuel moisture adjusted per drying ratio @@ -483,7 +485,7 @@ subroutine rate_of_spread ( currentSite ) ! Equation A4 in Thonicke et al. 2010 ! conversion of Rohtermal (1972) equation 12 in BTU/lb to current kJ/kg ! q_ig in kJ/kg - q_ig = 581.0_r8 +2594.0_r8 * currentPatch%fuel_eff_moist + q_ig = q_dry + 2594.0_r8 * currentPatch%fuel_eff_moist ! ---effective heating number--- ! Equation A3 in Thonicke et al. 2010. From ae0a5e7196fa9219b37712c140c55fa7a0c2d51b Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Thu, 5 Sep 2019 16:39:57 -0600 Subject: [PATCH 07/72] Change name for passive_crown fire intensity --- fire/SFMainMod.F90 | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 60ffde1a51..ef26fd81e1 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -924,7 +924,7 @@ subroutine crown_damage ( currentSite ) do while(associated(currentCohort)) currentCohort%fraction_crown_burned = 0.0_r8 - currentCohort%crown_FI = 0.0_r8 !critical fire intensity for passive crown fire + currentCohort%passive_crown_FI = 0.0_r8 !critical fire intensity for passive crown fire currentCohort%crown_fire_flg = 0 !flag for passvie crown fire ignition if (EDPftvarcon_inst%woody(currentCohort%pft) == 1) then !trees only @@ -937,12 +937,12 @@ subroutine crown_damage ( currentSite ) if (EDPftvarcon_inst%crown_fire(currentCohort%pft) == 1) then ! Crown fuel ignition potential, EQ 8 Bessie and Johnson 1995 - currentCohort%crown_FI = (0.01 * height_cbb * SF_val_crown_ignition_energy)**1.5 + currentCohort%passive_crown_FI = (0.01 * height_cbb * SF_val_crown_ignition_energy)**1.5 - if (currentCohort%crown_FI >= crown_fire_threshold) then ! 200 kW/m = threshold for crown fire potential + if (currentCohort%passive_crown_FI >= crown_fire_threshold) then ! 200 kW/m = threshold for crown fire potential - ! Initiation of passive crown fire, EQ 9 Bessie and Johnson 1995 - currentCohort%ignite_crown = currentPatch%FI/currentCohort%crown_FI + ! Initiation of passive crown fire, EQ 14 Bessie and Johnson 1995 + currentCohort%ignite_crown = currentPatch%FI/currentCohort%passive_crown_FI if (currentCohort%ignite_crown > 1.0_r8) then currentCohort%crown_fire_flg = 1 ! passive crown fire ignited From cafdac5295dff07a545fa0326971450b73cade60 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Thu, 5 Sep 2019 18:32:01 -0600 Subject: [PATCH 08/72] Revert q_dry to be local within fire ROS subroutine --- fire/SFMainMod.F90 | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index ef26fd81e1..81d285472e 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -60,8 +60,6 @@ module SFMainMod integer :: write_SF = 0 ! for debugging logical :: debug = .false. ! for debugging - - real(r8),parameter :: q_dry = 581.0_r8 !heat of pre-ignition of dry fuels (kJ/kg) ! ============================================================================ ! ============================================================================ @@ -443,7 +441,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 + 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; From 555838103d385cd5884ce27937a9c38d86682b70 Mon Sep 17 00:00:00 2001 From: jkshuman Date: Mon, 9 Sep 2019 10:12:24 -0600 Subject: [PATCH 09/72] Update fire/SFMainMod.F90 Co-Authored-By: Samuel Levis --- fire/SFMainMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 81d285472e..d27dfcb6ad 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -943,7 +943,7 @@ subroutine crown_damage ( currentSite ) ! Initiation of passive crown fire, EQ 14 Bessie and Johnson 1995 currentCohort%ignite_crown = currentPatch%FI/currentCohort%passive_crown_FI - if (currentCohort%ignite_crown > 1.0_r8) then + if (currentCohort%ignite_crown >= 1.0_r8) then currentCohort%crown_fire_flg = 1 ! passive crown fire ignited currentCohort%fraction_crown_burned = 1.0_r8 ! else ! evaluate crown damage based on scorch height From 1d2367137d742de35f415f2c8fae89ed96041d34 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Mon, 9 Sep 2019 14:54:37 -0600 Subject: [PATCH 10/72] Clean up calculation of fraction crown burnt --- fire/SFMainMod.F90 | 17 +++++++++-------- 1 file changed, 9 insertions(+), 8 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index d27dfcb6ad..5c2a7ba11f 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -928,9 +928,10 @@ subroutine crown_damage ( currentSite ) if (EDPftvarcon_inst%woody(currentCohort%pft) == 1) then !trees only - ! Calculate clear branch bole height at base of crown + ! height_cbb = clear branch bole height at base of crown (m) ! inst%crown = crown_depth_frac (PFT) - height_cbb = currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft) + height_cbb = currentCohort%hite * (1.0_r8 - EDPftvarcon_inst%crown(currentCohort%pft) + crown_depth = currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft) ! Evaluate for passive crown fire ignition if (EDPftvarcon_inst%crown_fire(currentCohort%pft) == 1) then @@ -954,18 +955,18 @@ subroutine crown_damage ( currentSite ) endif ! evaluate passive crown fire ! Flames lower than bottom of canopy. - ! c%hite is height of cohort - if (currentPatch%SH < (currentCohort%hite- height_cbb)) .and. currentCohort%crown_fire_flg = 0 then + ! cohort%hite is height of cohort + if (currentPatch%SH <= height_cbb .and. currentCohort%crown_fire_flg = 0) then currentCohort%fraction_crown_burned = 0.0_r8 else ! Flames part of way up canopy. ! Equation 17 in Thonicke et al. 2010. ! flames over bottom of canopy but not over top. - if ((currentCohort%hite > 0.0_r8).and.(currentPatch%SH >= & - (currentCohort%hite - height_cbb))) .and. currentCohort%crown_fire_flg = 0 then + if ((currentCohort%hite > 0.0_r8).and.(currentPatch%SH > height_cbb) & + .and. (currentPatch%SH <= currentCohort%hite) & + .and. currentCohort%crown_fire_flg = 0 then - currentCohort%fraction_crown_burned = (currentPatch%SH-currentCohort%hite*(1.0_r8- & - EDPftvarcon_inst%crown(currentCohort%pft)))/(height_cbb) + currentCohort%fraction_crown_burned = (currentPatch%SH - height_cbb)/crown_depth else ! Flames over top of canopy. From 764d6c46d7cae23883a2c0c1aca0f57cd5476bb9 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Mon, 9 Sep 2019 17:04:32 -0600 Subject: [PATCH 11/72] Re-order for clearer calculation of clear branch bole height --- fire/SFMainMod.F90 | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 5c2a7ba11f..e542e72d7b 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -930,8 +930,9 @@ subroutine crown_damage ( currentSite ) ! height_cbb = clear branch bole height at base of crown (m) ! inst%crown = crown_depth_frac (PFT) - height_cbb = currentCohort%hite * (1.0_r8 - EDPftvarcon_inst%crown(currentCohort%pft) crown_depth = currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft) + height_cbb = currentCohort%hite - crown_depth + ! Evaluate for passive crown fire ignition if (EDPftvarcon_inst%crown_fire(currentCohort%pft) == 1) then From a1ebec924c92c749cc145d698786deccf635d932 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Tue, 10 Sep 2019 09:46:51 -0600 Subject: [PATCH 12/72] Clean up else statements for scorch height fraction crown burnt --- fire/SFMainMod.F90 | 46 ++++++++++++++++++++-------------------------- 1 file changed, 20 insertions(+), 26 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index e542e72d7b..51dccdd464 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -114,7 +114,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 @@ -574,8 +574,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 @@ -664,8 +664,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 @@ -828,7 +828,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 @@ -902,7 +902,6 @@ 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 - !SF_val_crown_ignition_energy (kJ/kg) for crown foliage with 100% moisture (dry mass) use SFParamsMod, only : SF_val_crown_ignition_energy !kJ/kg @@ -912,7 +911,9 @@ subroutine crown_damage ( currentSite ) type(ed_patch_type) , pointer :: currentPatch type(ed_cohort_type), pointer :: currentCohort - real(r8) height_cbb ! lower crown base height (m) or clear branch bole height (m) + real(r8) crown_depth ! depth of crown (m) + real(r8) height_cbb ! clear branch bole height or crown base height (m) + currentPatch => currentSite%oldest_patch @@ -955,26 +956,19 @@ subroutine crown_damage ( currentSite ) ! else ! not crown fire plant endif ! evaluate passive crown fire - ! Flames lower than bottom of canopy. - ! cohort%hite is height of cohort + ! Flames lower than bottom of canopy. + ! height_cbb is clear branch bole height or height of bottom of canopy if (currentPatch%SH <= height_cbb .and. currentCohort%crown_fire_flg = 0) then currentCohort%fraction_crown_burned = 0.0_r8 - else - ! Flames part of way up canopy. - ! Equation 17 in Thonicke et al. 2010. - ! flames over bottom of canopy but not over top. - if ((currentCohort%hite > 0.0_r8).and.(currentPatch%SH > height_cbb) & - .and. (currentPatch%SH <= currentCohort%hite) & - .and. currentCohort%crown_fire_flg = 0 then - - currentCohort%fraction_crown_burned = (currentPatch%SH - height_cbb)/crown_depth - - else - ! Flames over top of canopy. - currentCohort%fraction_crown_burned = 1.0_r8 - endif - endif + ! Flames part way into canopy + ! Equation 17 in Thonicke et al. 2010 + elseif ((currentCohort%hite > 0.0_r8).and.(currentPatch%SH > height_cbb) & + .and. currentCohort%crown_fire_flg = 0 then + + currentCohort%fraction_crown_burned = max(0.0_r8, & + min(1.0_r8, (currentPatch%SH - height_cbb)/crown_depth)) + endif !SH frac crown burnt calculation ! Check for strange values. currentCohort%fraction_crown_burned = min(1.0_r8, max(0.0_r8,currentCohort%fraction_crown_burned)) endif !trees only @@ -1070,7 +1064,7 @@ subroutine post_fire_mortality ( currentSite ) ! Equation 22 in Thonicke et al. 2010. currentCohort%crownfire_mort = EDPftvarcon_inst%crown_resist(currentCohort%pft)*currentCohort%fraction_crown_burned**3.0_r8 ! Equation 18 in Thonicke et al. 2010. - currentCohort%fire_mort = max(0._r8,min(1.0_r8,currentCohort%crownfire_mort+currentCohort%cambial_mort- & + currentCohort%fire_mort = max(0.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. From a0e59a297f746d3733eceedd029d085147d9372c Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Tue, 10 Sep 2019 10:46:45 -0600 Subject: [PATCH 13/72] Add Cohort variables for crown damage to CohortDynamics --- biogeochem/EDCohortDynamicsMod.F90 | 2 ++ 1 file changed, 2 insertions(+) diff --git a/biogeochem/EDCohortDynamicsMod.F90 b/biogeochem/EDCohortDynamicsMod.F90 index 943c5abf37..79d2171fcc 100644 --- a/biogeochem/EDCohortDynamicsMod.F90 +++ b/biogeochem/EDCohortDynamicsMod.F90 @@ -535,6 +535,8 @@ subroutine nan_cohort(cc_p) ! FIRE currentCohort%fraction_crown_burned = nan ! proportion of crown affected by fire + currentCohort%passive_crown_FI = nan ! critical fire intensity for passive crown fire + currentCohort%crown_fire_flg = nan ! flag for 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 From 1365bf4ec754272a241ac1a3cb2194d31c074663 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Tue, 10 Sep 2019 10:59:29 -0600 Subject: [PATCH 14/72] Define variable name for crown fire damage --- biogeochem/EDCohortDynamicsMod.F90 | 1 + fire/SFMainMod.F90 | 11 +++++++---- 2 files changed, 8 insertions(+), 4 deletions(-) diff --git a/biogeochem/EDCohortDynamicsMod.F90 b/biogeochem/EDCohortDynamicsMod.F90 index 79d2171fcc..455633592e 100644 --- a/biogeochem/EDCohortDynamicsMod.F90 +++ b/biogeochem/EDCohortDynamicsMod.F90 @@ -537,6 +537,7 @@ subroutine nan_cohort(cc_p) currentCohort%fraction_crown_burned = nan ! proportion of crown affected by fire currentCohort%passive_crown_FI = nan ! critical fire intensity for passive crown fire currentCohort%crown_fire_flg = nan ! flag for passive crown fire ignition + currentCohort%ignite_crown = nan ! ratio to determine passive crown fire ignition, EQ14 Bessie&Johnson 1995 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 diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 51dccdd464..2ded57639a 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -912,7 +912,9 @@ subroutine crown_damage ( currentSite ) type(ed_cohort_type), pointer :: currentCohort real(r8) crown_depth ! depth of crown (m) - real(r8) height_cbb ! clear branch bole height or crown base height (m) + real(r8) height_cbb ! clear branch bole height or crown base height (m) + + real(r8), parameter :: crown_fire_threshold = 200.0_r8 ! threshold for crown fire potential (kW/m) currentPatch => currentSite%oldest_patch @@ -926,6 +928,7 @@ subroutine crown_damage ( currentSite ) currentCohort%fraction_crown_burned = 0.0_r8 currentCohort%passive_crown_FI = 0.0_r8 !critical fire intensity for passive crown fire currentCohort%crown_fire_flg = 0 !flag for passvie crown fire ignition + currentCohort%ignite_crown = 0.0_r8 !ratio of ignition of passive crown fire,EQ 14 Bessie & Johnson 1995 if (EDPftvarcon_inst%woody(currentCohort%pft) == 1) then !trees only @@ -939,7 +942,7 @@ subroutine crown_damage ( currentSite ) if (EDPftvarcon_inst%crown_fire(currentCohort%pft) == 1) then ! Crown fuel ignition potential, EQ 8 Bessie and Johnson 1995 - currentCohort%passive_crown_FI = (0.01 * height_cbb * SF_val_crown_ignition_energy)**1.5 + currentCohort%passive_crown_FI = (0.01_r8 * height_cbb * SF_val_crown_ignition_energy)**1.5_r8 if (currentCohort%passive_crown_FI >= crown_fire_threshold) then ! 200 kW/m = threshold for crown fire potential @@ -958,13 +961,13 @@ subroutine crown_damage ( currentSite ) ! Flames lower than bottom of canopy. ! height_cbb is clear branch bole height or height of bottom of canopy - if (currentPatch%SH <= height_cbb .and. currentCohort%crown_fire_flg = 0) then + if (currentPatch%SH <= height_cbb .and. currentCohort%crown_fire_flg = 0) then & currentCohort%fraction_crown_burned = 0.0_r8 ! Flames part way into canopy ! Equation 17 in Thonicke et al. 2010 elseif ((currentCohort%hite > 0.0_r8).and.(currentPatch%SH > height_cbb) & - .and. currentCohort%crown_fire_flg = 0 then + .and. currentCohort%crown_fire_flg = 0 then & currentCohort%fraction_crown_burned = max(0.0_r8, & min(1.0_r8, (currentPatch%SH - height_cbb)/crown_depth)) From 45377c91dfd681f0ab6ec91024120e0f7e79f9bd Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Tue, 10 Sep 2019 11:21:30 -0600 Subject: [PATCH 15/72] Zero new cohort values for crown damage --- biogeochem/EDCohortDynamicsMod.F90 | 5 ++++- fire/SFMainMod.F90 | 4 ++-- 2 files changed, 6 insertions(+), 3 deletions(-) diff --git a/biogeochem/EDCohortDynamicsMod.F90 b/biogeochem/EDCohortDynamicsMod.F90 index 455633592e..c43c8f3871 100644 --- a/biogeochem/EDCohortDynamicsMod.F90 +++ b/biogeochem/EDCohortDynamicsMod.F90 @@ -582,7 +582,10 @@ 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%fraction_crown_burned = 0._r8 + currentcohort%fraction_crown_burned = 0._r8 + currentCohort%passive_crown_FI = 0._r8 + currentCohort%crown_fire_flg = 0 + currentCohort%ignite_crown = 0.0_r8 currentCohort%size_class = 1 currentCohort%seed_prod = 0._r8 currentCohort%size_class_lasttimestep = 0 diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 2ded57639a..c23d7373d1 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -961,13 +961,13 @@ subroutine crown_damage ( currentSite ) ! Flames lower than bottom of canopy. ! height_cbb is clear branch bole height or height of bottom of canopy - if (currentPatch%SH <= height_cbb .and. currentCohort%crown_fire_flg = 0) then & + if (currentPatch%SH <= height_cbb .and. currentCohort%crown_fire_flg == 0) then & currentCohort%fraction_crown_burned = 0.0_r8 ! Flames part way into canopy ! Equation 17 in Thonicke et al. 2010 elseif ((currentCohort%hite > 0.0_r8).and.(currentPatch%SH > height_cbb) & - .and. currentCohort%crown_fire_flg = 0 then & + .and. currentCohort%crown_fire_flg == 0 then & currentCohort%fraction_crown_burned = max(0.0_r8, & min(1.0_r8, (currentPatch%SH - height_cbb)/crown_depth)) From 73c876294e1e72f03a79f2ba006bc75ea2a8dd7c Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Tue, 10 Sep 2019 11:32:01 -0600 Subject: [PATCH 16/72] Add missing parentheses crown damage conditional --- fire/SFMainMod.F90 | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index c23d7373d1..57896f2fd4 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -961,16 +961,16 @@ subroutine crown_damage ( currentSite ) ! Flames lower than bottom of canopy. ! height_cbb is clear branch bole height or height of bottom of canopy - if (currentPatch%SH <= height_cbb .and. currentCohort%crown_fire_flg == 0) then & + if (currentPatch%SH <= height_cbb) .and. (currentCohort%crown_fire_flg == 0) then & currentCohort%fraction_crown_burned = 0.0_r8 ! Flames part way into canopy ! Equation 17 in Thonicke et al. 2010 elseif ((currentCohort%hite > 0.0_r8).and.(currentPatch%SH > height_cbb) & - .and. currentCohort%crown_fire_flg == 0 then & + .and. (currentCohort%crown_fire_flg == 0) then & currentCohort%fraction_crown_burned = max(0.0_r8, & - min(1.0_r8, (currentPatch%SH - height_cbb)/crown_depth)) + min(1.0_r8, ((currentPatch%SH - height_cbb)/crown_depth))) endif !SH frac crown burnt calculation ! Check for strange values. currentCohort%fraction_crown_burned = min(1.0_r8, max(0.0_r8,currentCohort%fraction_crown_burned)) From c288d76d7724235bd87c95a47d63692b82580c0a Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Tue, 10 Sep 2019 12:31:06 -0600 Subject: [PATCH 17/72] Clean up crown damage subroutine --- fire/SFMainMod.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 57896f2fd4..f86486de05 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -961,13 +961,13 @@ subroutine crown_damage ( currentSite ) ! Flames lower than bottom of canopy. ! height_cbb is clear branch bole height or height of bottom of canopy - if (currentPatch%SH <= height_cbb) .and. (currentCohort%crown_fire_flg == 0) then & + if ((currentPatch%SH <= height_cbb) .and. (currentCohort%crown_fire_flg == 0)) then & currentCohort%fraction_crown_burned = 0.0_r8 ! Flames part way into canopy ! Equation 17 in Thonicke et al. 2010 elseif ((currentCohort%hite > 0.0_r8).and.(currentPatch%SH > height_cbb) & - .and. (currentCohort%crown_fire_flg == 0) then & + .and. (currentCohort%crown_fire_flg == 0)) then & currentCohort%fraction_crown_burned = max(0.0_r8, & min(1.0_r8, ((currentPatch%SH - height_cbb)/crown_depth))) From 0fc0c147080b0e958d4e9113a86bdb9327844bed Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Tue, 10 Sep 2019 12:42:31 -0600 Subject: [PATCH 18/72] Remove unnecssary line character crown damage routine --- fire/SFMainMod.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index f86486de05..f559213ec4 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -961,13 +961,13 @@ subroutine crown_damage ( currentSite ) ! Flames lower than bottom of canopy. ! height_cbb is clear branch bole height or height of bottom of canopy - if ((currentPatch%SH <= height_cbb) .and. (currentCohort%crown_fire_flg == 0)) then & + if ((currentPatch%SH <= height_cbb) .and. (currentCohort%crown_fire_flg == 0)) then currentCohort%fraction_crown_burned = 0.0_r8 ! Flames part way into canopy ! Equation 17 in Thonicke et al. 2010 elseif ((currentCohort%hite > 0.0_r8).and.(currentPatch%SH > height_cbb) & - .and. (currentCohort%crown_fire_flg == 0)) then & + .and. (currentCohort%crown_fire_flg == 0)) then currentCohort%fraction_crown_burned = max(0.0_r8, & min(1.0_r8, ((currentPatch%SH - height_cbb)/crown_depth))) From fa1973fafb9ec0d57e1151cd5e23f10c0f185d1b Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Tue, 10 Sep 2019 13:38:24 -0600 Subject: [PATCH 19/72] Define new crown damage variables in EDTypes --- main/EDTypesMod.F90 | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index 4dfcede3d2..ba5455c9e8 100644 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -344,6 +344,9 @@ module EDTypesMod ! FIRE real(r8) :: fraction_crown_burned ! proportion of crown affected by fire:- + real(r8) :: passive_crown_FI ! critical fire intensity for passive crown fire + real(r8) :: crown_fire_flg ! flag for passive crown fire ignition + real(r8) :: ignite_crown ! ratio to determine passive crown fire 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 @@ -1012,6 +1015,9 @@ 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%passive_crown_FI = ', ccohort%passive_crown_FI + write(fates_log(),*) 'co%crown_fire_flg = ', ccohort%crown_fire_flg + write(fates_log(),*) 'co%ignite_crown = ', ccohort%ignite_crown write(fates_log(),*) 'co%fire_mort = ', ccohort%fire_mort write(fates_log(),*) 'co%crownfire_mort = ', ccohort%crownfire_mort write(fates_log(),*) 'co%cambial_mort = ', ccohort%cambial_mort From ee2e576daa0a7408f45e7459f9ca24af4396ef9b Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Tue, 10 Sep 2019 13:41:42 -0600 Subject: [PATCH 20/72] Update to only "nan" crown damage variables in EDCohortDynamics --- biogeochem/EDCohortDynamicsMod.F90 | 6 +----- 1 file changed, 1 insertion(+), 5 deletions(-) diff --git a/biogeochem/EDCohortDynamicsMod.F90 b/biogeochem/EDCohortDynamicsMod.F90 index c43c8f3871..f780b46207 100644 --- a/biogeochem/EDCohortDynamicsMod.F90 +++ b/biogeochem/EDCohortDynamicsMod.F90 @@ -582,11 +582,7 @@ 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%fraction_crown_burned = 0._r8 - currentCohort%passive_crown_FI = 0._r8 - currentCohort%crown_fire_flg = 0 - currentCohort%ignite_crown = 0.0_r8 - currentCohort%size_class = 1 + currentCohort%size_class = 1 currentCohort%seed_prod = 0._r8 currentCohort%size_class_lasttimestep = 0 currentcohort%npp_acc_hold = 0._r8 From 39a39f3e68db272c19f781ca094f52d5f0cff89c Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Tue, 10 Sep 2019 14:37:26 -0600 Subject: [PATCH 21/72] Remove extra parentheses crown damage --- fire/SFMainMod.F90 | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index f559213ec4..abfe894400 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -961,13 +961,13 @@ subroutine crown_damage ( currentSite ) ! Flames lower than bottom of canopy. ! height_cbb is clear branch bole height or height of bottom of canopy - if ((currentPatch%SH <= height_cbb) .and. (currentCohort%crown_fire_flg == 0)) then + if (currentPatch%SH <= height_cbb .and. currentCohort%crown_fire_flg == 0) then currentCohort%fraction_crown_burned = 0.0_r8 ! Flames part way into canopy ! Equation 17 in Thonicke et al. 2010 - elseif ((currentCohort%hite > 0.0_r8).and.(currentPatch%SH > height_cbb) & - .and. (currentCohort%crown_fire_flg == 0)) then + elseif (currentCohort%hite > 0.0_r8 .and. currentPatch%SH > height_cbb & + .and. currentCohort%crown_fire_flg == 0) then currentCohort%fraction_crown_burned = max(0.0_r8, & min(1.0_r8, ((currentPatch%SH - height_cbb)/crown_depth))) From 703d87f4fbf903bc916466040f037fcc3cbc785b Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Tue, 10 Sep 2019 15:36:38 -0600 Subject: [PATCH 22/72] Add crown_fire_threshold read from EDTypes --- fire/SFMainMod.F90 | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index abfe894400..fe53e2e99b 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -22,6 +22,7 @@ module SFMainMod use EDtypesMod , only : AREA use EDtypesMod , only : DL_SF use EDtypesMod , only : FIRE_THRESHOLD + use EDtypesMod , only : crown_fire_threshold use EDTypesMod , only : TW_SF use EDtypesMod , only : LB_SF use EDtypesMod , only : LG_SF @@ -914,9 +915,6 @@ subroutine crown_damage ( currentSite ) real(r8) crown_depth ! depth of crown (m) real(r8) height_cbb ! clear branch bole height or crown base height (m) - real(r8), parameter :: crown_fire_threshold = 200.0_r8 ! threshold for crown fire potential (kW/m) - - currentPatch => currentSite%oldest_patch do while(associated(currentPatch)) From 7e691e01ef05fb99b19c2c3947584fd68fd8c108 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Tue, 10 Sep 2019 16:04:47 -0600 Subject: [PATCH 23/72] Update crown damage for clarity with surface fires --- fire/SFMainMod.F90 | 10 +++------- 1 file changed, 3 insertions(+), 7 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index fe53e2e99b..edede03bd6 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -957,18 +957,14 @@ subroutine crown_damage ( currentSite ) ! else ! not crown fire plant endif ! evaluate passive crown fire - ! Flames lower than bottom of canopy. + ! For surface fires, are flames in the canopy? ! height_cbb is clear branch bole height or height of bottom of canopy - if (currentPatch%SH <= height_cbb .and. currentCohort%crown_fire_flg == 0) then - currentCohort%fraction_crown_burned = 0.0_r8 - - ! Flames part way into canopy ! Equation 17 in Thonicke et al. 2010 - elseif (currentCohort%hite > 0.0_r8 .and. currentPatch%SH > height_cbb & + if (currentCohort%hite > 0.0_r8 .and. currentPatch%SH > height_cbb & .and. currentCohort%crown_fire_flg == 0) then currentCohort%fraction_crown_burned = max(0.0_r8, & - min(1.0_r8, ((currentPatch%SH - height_cbb)/crown_depth))) + min(1.0_r8, ((currentPatch%SH - height_cbb)/crown_depth))) endif !SH frac crown burnt calculation ! Check for strange values. currentCohort%fraction_crown_burned = min(1.0_r8, max(0.0_r8,currentCohort%fraction_crown_burned)) From 52b073460b5168e540f583e97755dc814233f1f6 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Wed, 11 Sep 2019 11:41:56 -0600 Subject: [PATCH 24/72] Update crown damage variables to be local to subroutine --- biogeochem/EDCohortDynamicsMod.F90 | 2 -- fire/SFMainMod.F90 | 25 ++++++++++++++----------- 2 files changed, 14 insertions(+), 13 deletions(-) diff --git a/biogeochem/EDCohortDynamicsMod.F90 b/biogeochem/EDCohortDynamicsMod.F90 index f780b46207..3b70366793 100644 --- a/biogeochem/EDCohortDynamicsMod.F90 +++ b/biogeochem/EDCohortDynamicsMod.F90 @@ -535,9 +535,7 @@ subroutine nan_cohort(cc_p) ! FIRE currentCohort%fraction_crown_burned = nan ! proportion of crown affected by fire - currentCohort%passive_crown_FI = nan ! critical fire intensity for passive crown fire currentCohort%crown_fire_flg = nan ! flag for passive crown fire ignition - currentCohort%ignite_crown = nan ! ratio to determine passive crown fire ignition, EQ14 Bessie&Johnson 1995 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 diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index edede03bd6..e78b496e9c 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -912,8 +912,10 @@ subroutine crown_damage ( currentSite ) type(ed_patch_type) , pointer :: currentPatch type(ed_cohort_type), pointer :: currentCohort - real(r8) crown_depth ! depth of crown (m) - real(r8) height_cbb ! clear branch bole height or crown base height (m) + 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) + real(r8) ignite_crown ! ratio for ignition of passive crown fire,EQ 14 Bessie & Johnson 1995 currentPatch => currentSite%oldest_patch @@ -924,30 +926,31 @@ subroutine crown_damage ( currentSite ) do while(associated(currentCohort)) currentCohort%fraction_crown_burned = 0.0_r8 - currentCohort%passive_crown_FI = 0.0_r8 !critical fire intensity for passive crown fire - currentCohort%crown_fire_flg = 0 !flag for passvie crown fire ignition - currentCohort%ignite_crown = 0.0_r8 !ratio of ignition of passive crown fire,EQ 14 Bessie & Johnson 1995 + 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) - crown_depth = currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft) - height_cbb = currentCohort%hite - crown_depth + crown_depth = currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft) + height_cbb = currentCohort%hite - crown_depth + passive_crown_FI = 0.0_r8 !critical fire intensity for passive crown fire + ignite_crown = 0.0_r8 !ratio for ignition of passive crown fire,EQ 14 Bessie & Johnson 1995 + currentCohort%crown_fire_flg = 0 !flag for passvie crown fire ignition ! Evaluate for passive crown fire ignition if (EDPftvarcon_inst%crown_fire(currentCohort%pft) == 1) then ! Crown fuel ignition potential, EQ 8 Bessie and Johnson 1995 - currentCohort%passive_crown_FI = (0.01_r8 * height_cbb * SF_val_crown_ignition_energy)**1.5_r8 + passive_crown_FI = (0.01_r8 * height_cbb * SF_val_crown_ignition_energy)**1.5_r8 - if (currentCohort%passive_crown_FI >= crown_fire_threshold) then ! 200 kW/m = threshold for crown fire potential + if (passive_crown_FI >= crown_fire_threshold) then ! 200 kW/m = threshold for crown fire potential ! Initiation of passive crown fire, EQ 14 Bessie and Johnson 1995 - currentCohort%ignite_crown = currentPatch%FI/currentCohort%passive_crown_FI + ignite_crown = currentPatch%FI/passive_crown_FI - if (currentCohort%ignite_crown >= 1.0_r8) then + if (ignite_crown >= 1.0_r8) then currentCohort%crown_fire_flg = 1 ! passive crown fire ignited currentCohort%fraction_crown_burned = 1.0_r8 ! else ! evaluate crown damage based on scorch height From cce9213a7f0795d555adae3895d0f0a683f6c723 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Wed, 11 Sep 2019 12:11:45 -0600 Subject: [PATCH 25/72] Remove local crown damage variables from EDTypes --- main/EDTypesMod.F90 | 2 -- 1 file changed, 2 deletions(-) diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index ba5455c9e8..ab54367fc4 100644 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -344,9 +344,7 @@ module EDTypesMod ! FIRE real(r8) :: fraction_crown_burned ! proportion of crown affected by fire:- - real(r8) :: passive_crown_FI ! critical fire intensity for passive crown fire real(r8) :: crown_fire_flg ! flag for passive crown fire ignition - real(r8) :: ignite_crown ! ratio to determine passive crown fire 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 From eb6e6d153eabd068e056c710937ab97d61790fc1 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Wed, 11 Sep 2019 12:20:12 -0600 Subject: [PATCH 26/72] Remove local crown damage variables from EDTypes, another location --- main/EDTypesMod.F90 | 2 -- 1 file changed, 2 deletions(-) diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index ab54367fc4..2c25d6eeda 100644 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -1013,9 +1013,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%passive_crown_FI = ', ccohort%passive_crown_FI write(fates_log(),*) 'co%crown_fire_flg = ', ccohort%crown_fire_flg - write(fates_log(),*) 'co%ignite_crown = ', ccohort%ignite_crown write(fates_log(),*) 'co%fire_mort = ', ccohort%fire_mort write(fates_log(),*) 'co%crownfire_mort = ', ccohort%crownfire_mort write(fates_log(),*) 'co%cambial_mort = ', ccohort%cambial_mort From dc1ba209135c54f18b99921f8ca73a77962aa9c9 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Thu, 12 Sep 2019 09:04:55 -0600 Subject: [PATCH 27/72] Zero cambial_mort inside cambial_damage_kill not zero_cohort --- biogeochem/EDCohortDynamicsMod.F90 | 2 -- fire/SFMainMod.F90 | 5 ++--- 2 files changed, 2 insertions(+), 5 deletions(-) diff --git a/biogeochem/EDCohortDynamicsMod.F90 b/biogeochem/EDCohortDynamicsMod.F90 index 3b70366793..4d7bba06bb 100644 --- a/biogeochem/EDCohortDynamicsMod.F90 +++ b/biogeochem/EDCohortDynamicsMod.F90 @@ -595,8 +595,6 @@ subroutine zero_cohort(cc_p) currentCohort%leaf_cost = 0._r8 currentcohort%excl_weight = 0._r8 currentcohort%prom_weight = 0._r8 - currentcohort%crownfire_mort = 0._r8 - currentcohort%cambial_mort = 0._r8 currentCohort%c13disc_clm = 0._r8 currentCohort%c13disc_acc = 0._r8 diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index e78b496e9c..6c2a29c8a1 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -1007,7 +1007,8 @@ subroutine cambial_damage_kill ( currentSite ) if (currentPatch%fire == 1) then currentCohort => currentPatch%tallest; - do while(associated(currentCohort)) + do while(associated(currentCohort)) + 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. @@ -1019,8 +1020,6 @@ subroutine cambial_damage_kill ( currentSite ) else if ((currentPatch%tau_l/tau_c) > 0.22_r8) then currentCohort%cambial_mort = (0.563_r8*(currentPatch%tau_l/tau_c)) - 0.125_r8 - else - currentCohort%cambial_mort = 0.0_r8 endif endif endif !trees From e11602f9d5ef7c184851f9dd1ffcd8962752b7ed Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Thu, 12 Sep 2019 15:48:22 -0600 Subject: [PATCH 28/72] Update drying_ratio value in params file --- parameter_files/fates_params_default.cdl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/parameter_files/fates_params_default.cdl b/parameter_files/fates_params_default.cdl index d98fc2e1ea..8abab650ae 100644 --- a/parameter_files/fates_params_default.cdl +++ b/parameter_files/fates_params_default.cdl @@ -1112,7 +1112,7 @@ fates_fire_crown_fire = 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 ; fates_cwd_flig = 0.24 ; - fates_fire_drying_ratio = 13000 ; + fates_fire_drying_ratio = 66000 ; fates_fire_durat_slope = -11.06 ; From 5587782b138f48b350413fad8946bae47c9ba333 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Wed, 25 Sep 2019 16:46:58 -0600 Subject: [PATCH 29/72] Update crown_ignite_energy to PFT level --- fire/SFMainMod.F90 | 6 +++--- fire/SFParamsMod.F90 | 9 --------- main/EDPftvarcon.F90 | 10 ++++++++++ parameter_files/fates_params_default.cdl | 12 ++++++------ 4 files changed, 19 insertions(+), 18 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 90dffb71a4..6d367c1937 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -890,8 +890,7 @@ 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 - !SF_val_crown_ignition_energy (kJ/kg) for crown foliage with 100% moisture (dry mass) - use SFParamsMod, only : SF_val_crown_ignition_energy !kJ/kg + type(ed_site_type), intent(in), target :: currentSite @@ -929,7 +928,8 @@ subroutine crown_damage ( currentSite ) if (EDPftvarcon_inst%crown_fire(currentCohort%pft) == 1) then ! Crown fuel ignition potential, EQ 8 Bessie and Johnson 1995 - passive_crown_FI = (0.01_r8 * height_cbb * SF_val_crown_ignition_energy)**1.5_r8 + ! Note: future crown_ignition_energy to be calculated based on foliar moisture content from FATES-Hydro + passive_crown_FI = (0.01_r8 * height_cbb *EDPftvarcon_inst%crown_ignite_energy(currentCohort%pft)**1.5_r8 if (passive_crown_FI >= crown_fire_threshold) then ! 200 kW/m = threshold for crown fire potential diff --git a/fire/SFParamsMod.F90 b/fire/SFParamsMod.F90 index dddef4d927..1c9f278c1e 100644 --- a/fire/SFParamsMod.F90 +++ b/fire/SFParamsMod.F90 @@ -23,7 +23,6 @@ module SFParamsMod real(r8),protected, public :: SF_val_fdi_alpha real(r8),protected, public :: SF_val_miner_total real(r8),protected, public :: SF_val_fuel_energy - real(r8),protected, public :: SF_val_crown_ignition_energy real(r8),protected, public :: SF_val_part_dens real(r8),protected, public :: SF_val_miner_damp real(r8),protected, public :: SF_val_max_durat @@ -45,7 +44,6 @@ module SFParamsMod character(len=param_string_length),parameter :: SF_name_fdi_alpha = "fates_fire_fdi_alpha" character(len=param_string_length),parameter :: SF_name_miner_total = "fates_fire_miner_total" character(len=param_string_length),parameter :: SF_name_fuel_energy = "fates_fire_fuel_energy" - character(len=param_string_length),parameter :: SF_name_crown_ignition_energy = "fates_fire_crown_ignition_energy" character(len=param_string_length),parameter :: SF_name_part_dens = "fates_fire_part_dens" character(len=param_string_length),parameter :: SF_name_miner_damp = "fates_fire_miner_damp" character(len=param_string_length),parameter :: SF_name_max_durat = "fates_fire_max_durat" @@ -139,7 +137,6 @@ subroutine SpitFireParamsInit() SF_val_fdi_alpha = nan SF_val_miner_total = nan SF_val_fuel_energy = nan - SF_val_crown_ignition_energy = nan SF_val_part_dens = nan SF_val_miner_damp = nan SF_val_max_durat = nan @@ -218,9 +215,6 @@ subroutine SpitFireRegisterScalars(fates_params) call fates_params%RegisterParameter(name=SF_name_fuel_energy, dimension_shape=dimension_shape_scalar, & dimension_names=dim_names_scalar) - call fates_params%RegisterParameter(name=SF_name_crown_ignition_energy, dimension_shape=dimension_shape_scalar, & - dimension_names=dim_names_scalar) - call fates_params%RegisterParameter(name=SF_name_part_dens, dimension_shape=dimension_shape_scalar, & dimension_names=dim_names_scalar) @@ -262,9 +256,6 @@ subroutine SpitFireReceiveScalars(fates_params) call fates_params%RetreiveParameter(name=SF_name_fuel_energy, & data=SF_val_fuel_energy) - call fates_params%RetreiveParameter(name=SF_name_crown_ignition_energy, & - data=SF_val_crown_ignition_energy) - call fates_params%RetreiveParameter(name=SF_name_part_dens, & data=SF_val_part_dens) diff --git a/main/EDPftvarcon.F90 b/main/EDPftvarcon.F90 index 7fbda29121..0fe6c4a97a 100644 --- a/main/EDPftvarcon.F90 +++ b/main/EDPftvarcon.F90 @@ -50,6 +50,7 @@ module EDPftvarcon real(r8), allocatable :: bark_scaler(:) ! scaler from dbh to bark thickness. For fire model. real(r8), allocatable :: crown_resist(:) ! resistance to fire. (1 = none) For fire model. real(r8), allocatable :: crown_fire(:) ! Does plant have passive crown fire? (1=yes, 0=no) + real(r8), allocatable :: crown_ignite_energy(:) ! heat of crown foliage ignition [kJ/kg] real(r8), allocatable :: initd(:) ! initial seedling density real(r8), allocatable :: seed_suppl(:) ! seeds that come from outside the gridbox. real(r8), allocatable :: BB_slope(:) ! ball berry slope parameter @@ -386,6 +387,10 @@ subroutine Register_PFT(this, fates_params) call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) + name = 'fates_fire_crown_ignition_energy' + call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & + dimension_names=dim_names, lower_bounds=dim_lower_bound) + name = 'fates_recruit_initd' call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) @@ -822,6 +827,10 @@ subroutine Receive_PFT(this, fates_params) call fates_params%RetreiveParameterAllocate(name=name, & data=this%crown_fire) + name = 'fates_fire_crown_ignition_energy' + call fates_params%RetreiveParameterAllocate(name=name, & + data=this%crown_ignition_energy) + name = 'fates_recruit_initd' call fates_params%RetreiveParameterAllocate(name=name, & data=this%initd) @@ -1732,6 +1741,7 @@ subroutine FatesReportPFTParams(is_master) write(fates_log(),fmt0) 'bark_scaler = ',EDPftvarcon_inst%bark_scaler write(fates_log(),fmt0) 'crown_resist = ',EDPftvarcon_inst%crown_resist write(fates_log(),fmt0) 'crown_fire = ',EDPftvarcon_inst%crown_fire + write(fates_log(),fmt0) 'crown_ignition_energy = ',EDPftvarcon_inst%crown_ignition_energy write(fates_log(),fmt0) 'initd = ',EDPftvarcon_inst%initd write(fates_log(),fmt0) 'seed_suppl = ',EDPftvarcon_inst%seed_suppl write(fates_log(),fmt0) 'BB_slope = ',EDPftvarcon_inst%BB_slope diff --git a/parameter_files/fates_params_default.cdl b/parameter_files/fates_params_default.cdl index 179dc01b35..b32045b6e3 100644 --- a/parameter_files/fates_params_default.cdl +++ b/parameter_files/fates_params_default.cdl @@ -135,6 +135,9 @@ variables: double fates_fire_crown_fire(fates_pft) ; fates_fire_crown_fire:units = "logical flag" ; fates_fire_crown_fire:long_name = "Binary flag for passive crown fire"; + double fates_fire_crown_ignite_energy(fates_pft) ; + fates_fire_crown_ignite_energy:units = "kJ/kg" ; + fates_fire_crown_ignite_energy:long_name = "spitfire parameter, heat of crown foliage ignition" ; double fates_fr_fcel(fates_pft) ; fates_fr_fcel:units = "fraction" ; fates_fr_fcel:long_name = "Fine root litter cellulose fraction" ; @@ -495,9 +498,6 @@ variables: double fates_fire_fuel_energy ; fates_fire_fuel_energy:units = "kJ/kg" ; fates_fire_fuel_energy:long_name = "spitfire parameter, heat content of fuel" ; - double fates_fire_crown_ignition_energy ; - fates_fire_crown_ignition_energy:units = "kJ/kg" ; - fates_fire_crown_ignition_energy:long_name = "spitfire parameter, heat of crown foliage ignition" ; double fates_fire_max_durat ; fates_fire_max_durat:units = "minutes" ; fates_fire_max_durat:long_name = "spitfire parameter, fire maximum duration, Equation 14 Thonicke et al 2010" ; @@ -724,7 +724,9 @@ data: fates_fire_crown_resist = 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775 ; -fates_fire_crown_fire = 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 ; + fates_fire_crown_fire = 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 ; + + fates_fire_crown_ignite_energy = 3060, 3060, 3060, 3060, 3060, 3060, 3060, 3060, 3060, 3060, 3060, 3060 ; fates_fr_fcel = 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5 ; @@ -1124,8 +1126,6 @@ fates_fire_crown_fire = 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 ; fates_fire_fuel_energy = 18000 ; - fates_fire_crown_ignition_energy = 3060 ; - fates_fire_max_durat = 240 ; fates_fire_miner_damp = 0.41739 ; From 76ed1fade6e8d3b15bc62a9295218084fad2c95b Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Wed, 25 Sep 2019 17:00:14 -0600 Subject: [PATCH 30/72] Update crown_ignite name in EDPftvarcon --- main/EDPftvarcon.F90 | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/main/EDPftvarcon.F90 b/main/EDPftvarcon.F90 index 0fe6c4a97a..38b35f996c 100644 --- a/main/EDPftvarcon.F90 +++ b/main/EDPftvarcon.F90 @@ -387,7 +387,7 @@ subroutine Register_PFT(this, fates_params) call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) - name = 'fates_fire_crown_ignition_energy' + name = 'fates_fire_crown_ignite_energy' call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) @@ -827,9 +827,9 @@ subroutine Receive_PFT(this, fates_params) call fates_params%RetreiveParameterAllocate(name=name, & data=this%crown_fire) - name = 'fates_fire_crown_ignition_energy' + name = 'fates_fire_crown_ignite_energy' call fates_params%RetreiveParameterAllocate(name=name, & - data=this%crown_ignition_energy) + data=this%crown_ignite_energy) name = 'fates_recruit_initd' call fates_params%RetreiveParameterAllocate(name=name, & @@ -1741,7 +1741,7 @@ subroutine FatesReportPFTParams(is_master) write(fates_log(),fmt0) 'bark_scaler = ',EDPftvarcon_inst%bark_scaler write(fates_log(),fmt0) 'crown_resist = ',EDPftvarcon_inst%crown_resist write(fates_log(),fmt0) 'crown_fire = ',EDPftvarcon_inst%crown_fire - write(fates_log(),fmt0) 'crown_ignition_energy = ',EDPftvarcon_inst%crown_ignition_energy + write(fates_log(),fmt0) 'crown_ignite_energy = ',EDPftvarcon_inst%crown_ignite_energy write(fates_log(),fmt0) 'initd = ',EDPftvarcon_inst%initd write(fates_log(),fmt0) 'seed_suppl = ',EDPftvarcon_inst%seed_suppl write(fates_log(),fmt0) 'BB_slope = ',EDPftvarcon_inst%BB_slope From 2feddcb86e402519c4849108584e3220879fb9f4 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Wed, 25 Sep 2019 17:05:34 -0600 Subject: [PATCH 31/72] Add paran in SFMain --- fire/SFMainMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 6d367c1937..6050f0d062 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -929,7 +929,7 @@ subroutine crown_damage ( currentSite ) ! Crown fuel ignition potential, EQ 8 Bessie and Johnson 1995 ! Note: future crown_ignition_energy to be calculated based on foliar moisture content from FATES-Hydro - passive_crown_FI = (0.01_r8 * height_cbb *EDPftvarcon_inst%crown_ignite_energy(currentCohort%pft)**1.5_r8 + passive_crown_FI = (0.01_r8 * height_cbb *EDPftvarcon_inst%crown_ignite_energy(currentCohort%pft))**1.5_r8 if (passive_crown_FI >= crown_fire_threshold) then ! 200 kW/m = threshold for crown fire potential From 0aef4a978b022f3483e088081ba69e06ce610ee4 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Mon, 30 Sep 2019 10:44:29 -0600 Subject: [PATCH 32/72] Add comment with equation for crown_ignite_energy per VanWagner 1977 --- fire/SFMainMod.F90 | 6 +++++- 1 file changed, 5 insertions(+), 1 deletion(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 6050f0d062..bb5c264e9a 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -927,8 +927,12 @@ subroutine crown_damage ( currentSite ) ! Evaluate for passive crown fire ignition if (EDPftvarcon_inst%crown_fire(currentCohort%pft) == 1) then + ! Note: crown_ignition_energy to be calculated based on foliar moisture content from FATES-Hydro + ! EQ 3 Van Wagner 1977 + ! crown ignite_energy (kJ/kg), m = foliar moisture content based on dry fuel (%) + ! crown_ignite_energy = 460 + 26 * m + ! Crown fuel ignition potential, EQ 8 Bessie and Johnson 1995 - ! Note: future crown_ignition_energy to be calculated based on foliar moisture content from FATES-Hydro passive_crown_FI = (0.01_r8 * height_cbb *EDPftvarcon_inst%crown_ignite_energy(currentCohort%pft))**1.5_r8 if (passive_crown_FI >= crown_fire_threshold) then ! 200 kW/m = threshold for crown fire potential From 8d582b74b760ff573e76d5fb07e94b99d48d95a8 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Mon, 30 Sep 2019 10:58:52 -0600 Subject: [PATCH 33/72] Add comment with details for empirical value C in passive_crown_FI --- fire/SFMainMod.F90 | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index bb5c264e9a..e49c68cacb 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -929,10 +929,12 @@ subroutine crown_damage ( currentSite ) ! Note: crown_ignition_energy to be calculated based on foliar moisture content from FATES-Hydro ! EQ 3 Van Wagner 1977 - ! crown ignite_energy (kJ/kg), m = foliar moisture content based on dry fuel (%) + ! h = crown_ignite_energy (kJ/kg), m = foliar moisture content based on dry fuel (%) ! crown_ignite_energy = 460 + 26 * m ! Crown fuel ignition potential, EQ 8 Bessie and Johnson 1995 + ! Van Wagner EQ4 FI = (Czh)**3/2 where z=crown base height,h=heat crown ignite energy, FI=fire intensity + ! 0.01 = C from Van Wagner 1977 EQ4 for crown of base height 6m, 100% FMC, and FI 2500kW/m passive_crown_FI = (0.01_r8 * height_cbb *EDPftvarcon_inst%crown_ignite_energy(currentCohort%pft))**1.5_r8 if (passive_crown_FI >= crown_fire_threshold) then ! 200 kW/m = threshold for crown fire potential From 29fd4ff373d3b3cf3a4bce53ba6b780be99b7096 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Mon, 30 Sep 2019 15:09:00 -0600 Subject: [PATCH 34/72] Update flag to logical, clean up conditional --- biogeochem/EDCohortDynamicsMod.F90 | 10 ++++----- fire/SFMainMod.F90 | 26 ++++++++++++------------ main/EDPftvarcon.F90 | 20 +++++++++--------- parameter_files/fates_params_default.cdl | 10 ++++----- 4 files changed, 33 insertions(+), 33 deletions(-) diff --git a/biogeochem/EDCohortDynamicsMod.F90 b/biogeochem/EDCohortDynamicsMod.F90 index 4d7bba06bb..c129aace5d 100644 --- a/biogeochem/EDCohortDynamicsMod.F90 +++ b/biogeochem/EDCohortDynamicsMod.F90 @@ -534,11 +534,11 @@ 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%crown_fire_flg = nan ! flag for 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 = nan ! flag for 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 end subroutine nan_cohort diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index e49c68cacb..a80d72542c 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -897,10 +897,11 @@ subroutine crown_damage ( currentSite ) type(ed_patch_type) , pointer :: currentPatch type(ed_cohort_type), pointer :: currentCohort - 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) - real(r8) ignite_crown ! ratio for ignition of passive crown fire,EQ 14 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) + real(r8) :: ignite_crown ! ratio for ignition of passive crown fire,EQ 14 Bessie & Johnson 1995 + logical :: passive_crown_flg ! flag to check for passive crown fire ignition currentPatch => currentSite%oldest_patch @@ -917,11 +918,11 @@ subroutine crown_damage ( currentSite ) ! height_cbb = clear branch bole height at base of crown (m) ! inst%crown = crown_depth_frac (PFT) - crown_depth = currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft) - height_cbb = currentCohort%hite - crown_depth - passive_crown_FI = 0.0_r8 !critical fire intensity for passive crown fire - ignite_crown = 0.0_r8 !ratio for ignition of passive crown fire,EQ 14 Bessie & Johnson 1995 - currentCohort%crown_fire_flg = 0 !flag for passvie crown fire ignition + crown_depth = currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft) + height_cbb = currentCohort%hite - crown_depth + passive_crown_FI = 0.0_r8 + ignite_crown = 0.0_r8 + currentCohort%passive_crown_fire_flg = .false. ! Evaluate for passive crown fire ignition @@ -943,7 +944,7 @@ subroutine crown_damage ( currentSite ) ignite_crown = currentPatch%FI/passive_crown_FI if (ignite_crown >= 1.0_r8) then - currentCohort%crown_fire_flg = 1 ! passive crown fire ignited + currentCohort%passive_crown_fire_flg = .true. ! passive crown fire ignited currentCohort%fraction_crown_burned = 1.0_r8 ! else ! evaluate crown damage based on scorch height endif ! ignite passive crown fire @@ -956,10 +957,9 @@ subroutine crown_damage ( currentSite ) ! 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%crown_fire_flg == 0) then + .and. currentCohort%passive_crown_fire_flg = .true.) then - currentCohort%fraction_crown_burned = max(0.0_r8, & - min(1.0_r8, ((currentPatch%SH - height_cbb)/crown_depth))) + currentCohort%fraction_crown_burned = min(1.0_r8, ((currentPatch%SH - height_cbb)/crown_depth)) endif !SH frac crown burnt calculation ! Check for strange values. currentCohort%fraction_crown_burned = min(1.0_r8, max(0.0_r8,currentCohort%fraction_crown_burned)) diff --git a/main/EDPftvarcon.F90 b/main/EDPftvarcon.F90 index 38b35f996c..7d9aa3d900 100644 --- a/main/EDPftvarcon.F90 +++ b/main/EDPftvarcon.F90 @@ -48,8 +48,8 @@ module EDPftvarcon real(r8), allocatable :: crown(:) ! fraction of the height of the plant ! that is occupied by crown. For fire model. real(r8), allocatable :: bark_scaler(:) ! scaler from dbh to bark thickness. For fire model. - real(r8), allocatable :: crown_resist(:) ! resistance to fire. (1 = none) For fire model. - real(r8), allocatable :: crown_fire(:) ! Does plant have passive crown fire? (1=yes, 0=no) + real(r8), allocatable :: crown_kill(:) ! crown resistance to fire. (1 = none) For fire model. + real(r8), allocatable :: passive_crown_fire(:) ! Is plant susceptible to passive crown fire?(1=yes,0=no) real(r8), allocatable :: crown_ignite_energy(:) ! heat of crown foliage ignition [kJ/kg] real(r8), allocatable :: initd(:) ! initial seedling density real(r8), allocatable :: seed_suppl(:) ! seeds that come from outside the gridbox. @@ -379,11 +379,11 @@ subroutine Register_PFT(this, fates_params) call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) - name = 'fates_fire_crown_resist' + name = 'fates_fire_crown_kill' call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) - name = 'fates_fire_crown_fire' + name = 'fates_fire_passive_crown_fire' call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) @@ -819,13 +819,13 @@ subroutine Receive_PFT(this, fates_params) call fates_params%RetreiveParameterAllocate(name=name, & data=this%bark_scaler) - name = 'fates_fire_crown_resist' + name = 'fates_fire_crown_kill' call fates_params%RetreiveParameterAllocate(name=name, & - data=this%crown_resist) + data=this%crown_kill) - name = 'fates_fire_crown_fire' + name = 'fates_fire_passive_crown_fire' call fates_params%RetreiveParameterAllocate(name=name, & - data=this%crown_fire) + data=this%passive_crown_fire) name = 'fates_fire_crown_ignite_energy' call fates_params%RetreiveParameterAllocate(name=name, & @@ -1739,8 +1739,8 @@ subroutine FatesReportPFTParams(is_master) write(fates_log(),fmt0) 'leaf_stor_priority = ',EDPftvarcon_inst%leaf_stor_priority write(fates_log(),fmt0) 'crown = ',EDPftvarcon_inst%crown write(fates_log(),fmt0) 'bark_scaler = ',EDPftvarcon_inst%bark_scaler - write(fates_log(),fmt0) 'crown_resist = ',EDPftvarcon_inst%crown_resist - write(fates_log(),fmt0) 'crown_fire = ',EDPftvarcon_inst%crown_fire + write(fates_log(),fmt0) 'crown_kill = ',EDPftvarcon_inst%crown_kill + write(fates_log(),fmt0) 'passive_crown_fire = ',EDPftvarcon_inst%passive_crown_fire write(fates_log(),fmt0) 'crown_ignite_energy = ',EDPftvarcon_inst%crown_ignite_energy write(fates_log(),fmt0) 'initd = ',EDPftvarcon_inst%initd write(fates_log(),fmt0) 'seed_suppl = ',EDPftvarcon_inst%seed_suppl diff --git a/parameter_files/fates_params_default.cdl b/parameter_files/fates_params_default.cdl index b32045b6e3..e6d56f9f07 100644 --- a/parameter_files/fates_params_default.cdl +++ b/parameter_files/fates_params_default.cdl @@ -129,9 +129,9 @@ variables: double fates_fire_crown_depth_frac(fates_pft) ; fates_fire_crown_depth_frac:units = "fraction" ; fates_fire_crown_depth_frac:long_name = "the depth of a cohorts crown as a fraction of its height" ; - double fates_fire_crown_resist(fates_pft) ; - fates_fire_crown_resist:units = "NA" ; - fates_fire_crown_resist:long_name = "resistance to fire (1 = none), EQ 22 in Thonicke et al 2010" ; + double fates_fire_crown_kill(fates_pft) ; + fates_fire_crown_kill:units = "NA" ; + fates_fire_crown_kill:long_name = "crown resistance to fire (1 = none), EQ 22 in Thonicke et al 2010" ; double fates_fire_crown_fire(fates_pft) ; fates_fire_crown_fire:units = "logical flag" ; fates_fire_crown_fire:long_name = "Binary flag for passive crown fire"; @@ -721,10 +721,10 @@ data: fates_fire_crown_depth_frac = 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.95, 0.95, 0.95, 1, 1, 1 ; - fates_fire_crown_resist = 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, + fates_fire_crown_kill = 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775 ; - fates_fire_crown_fire = 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 ; + fates_fire_passive_crown_fire = 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 ; fates_fire_crown_ignite_energy = 3060, 3060, 3060, 3060, 3060, 3060, 3060, 3060, 3060, 3060, 3060, 3060 ; From cb67c913c637f004259f6e31eaf6368b6b3a9d60 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Mon, 30 Sep 2019 18:25:05 -0600 Subject: [PATCH 35/72] Revert passive crown fire flag to binary --- biogeochem/EDCohortDynamicsMod.F90 | 2 +- fire/SFMainMod.F90 | 15 +++++++-------- 2 files changed, 8 insertions(+), 9 deletions(-) diff --git a/biogeochem/EDCohortDynamicsMod.F90 b/biogeochem/EDCohortDynamicsMod.F90 index c129aace5d..88a0f02eb7 100644 --- a/biogeochem/EDCohortDynamicsMod.F90 +++ b/biogeochem/EDCohortDynamicsMod.F90 @@ -535,7 +535,7 @@ subroutine nan_cohort(cc_p) ! FIRE currentCohort%fraction_crown_burned = nan ! proportion of crown affected by fire - currentCohort%passive_crown_fire_flg = nan ! flag for passive crown fire ignition + 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 diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index a80d72542c..1067aafb2d 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -901,7 +901,6 @@ subroutine crown_damage ( currentSite ) 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) real(r8) :: ignite_crown ! ratio for ignition of passive crown fire,EQ 14 Bessie & Johnson 1995 - logical :: passive_crown_flg ! flag to check for passive crown fire ignition currentPatch => currentSite%oldest_patch @@ -918,11 +917,11 @@ subroutine crown_damage ( currentSite ) ! height_cbb = clear branch bole height at base of crown (m) ! inst%crown = crown_depth_frac (PFT) - crown_depth = currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft) - height_cbb = currentCohort%hite - crown_depth - passive_crown_FI = 0.0_r8 - ignite_crown = 0.0_r8 - currentCohort%passive_crown_fire_flg = .false. + crown_depth = currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft) + height_cbb = currentCohort%hite - crown_depth + passive_crown_FI = 0.0_r8 + ignite_crown = 0.0_r8 + currentCohort%passive_crown_fire_flg = 0 ! Evaluate for passive crown fire ignition @@ -944,7 +943,7 @@ subroutine crown_damage ( currentSite ) ignite_crown = currentPatch%FI/passive_crown_FI if (ignite_crown >= 1.0_r8) then - currentCohort%passive_crown_fire_flg = .true. ! passive crown fire ignited + currentCohort%passive_crown_fire_flg = 1 ! passive crown fire ignited currentCohort%fraction_crown_burned = 1.0_r8 ! else ! evaluate crown damage based on scorch height endif ! ignite passive crown fire @@ -957,7 +956,7 @@ subroutine crown_damage ( currentSite ) ! 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 = .true.) then + .and. currentCohort%passive_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 From 5b843566b287893697e1fd99949a2f2e331dfa84 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Tue, 1 Oct 2019 14:43:36 -0600 Subject: [PATCH 36/72] Update variable names in EDTypes --- fire/SFMainMod.F90 | 2 +- main/EDTypesMod.F90 | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 1067aafb2d..6e87575051 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -1052,7 +1052,7 @@ subroutine post_fire_mortality ( currentSite ) currentCohort%crownfire_mort = 0.0_r8 if (EDPftvarcon_inst%woody(currentCohort%pft) == 1) then ! Equation 22 in Thonicke et al. 2010. - currentCohort%crownfire_mort = EDPftvarcon_inst%crown_resist(currentCohort%pft)*currentCohort%fraction_crown_burned**3.0_r8 + 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%crownfire_mort*currentCohort%cambial_mort))) !joint prob. diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index 0594906bfe..70110951c2 100644 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -1011,7 +1011,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%crown_fire_flg = ', ccohort%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 write(fates_log(),*) 'co%cambial_mort = ', ccohort%cambial_mort From 188e6838bb3fb295fa31b4b596eb738c491c0053 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Tue, 1 Oct 2019 14:51:52 -0600 Subject: [PATCH 37/72] Update passive_crown_fire in EDTypes --- main/EDTypesMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index 70110951c2..131931b5c2 100644 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -344,7 +344,7 @@ module EDTypesMod ! FIRE real(r8) :: fraction_crown_burned ! proportion of crown affected by fire:- - real(r8) :: crown_fire_flg ! flag for passive crown fire ignition + real(r8) :: 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 From 4f017a2785beea12fd233a478d2febb6ffc6518e Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Tue, 1 Oct 2019 15:05:31 -0600 Subject: [PATCH 38/72] Update variable name and comments crown_fire flag --- main/EDPftvarcon.F90 | 10 +++++----- parameter_files/fates_params_default.cdl | 4 ++-- 2 files changed, 7 insertions(+), 7 deletions(-) diff --git a/main/EDPftvarcon.F90 b/main/EDPftvarcon.F90 index 7d9aa3d900..359c66f5ce 100644 --- a/main/EDPftvarcon.F90 +++ b/main/EDPftvarcon.F90 @@ -49,7 +49,7 @@ module EDPftvarcon ! that is occupied by crown. For fire model. real(r8), allocatable :: bark_scaler(:) ! scaler from dbh to bark thickness. For fire model. real(r8), allocatable :: crown_kill(:) ! crown resistance to fire. (1 = none) For fire model. - real(r8), allocatable :: passive_crown_fire(:) ! Is plant susceptible to passive crown fire?(1=yes,0=no) + real(r8), allocatable :: crown_fire(:) ! Is plant susceptible to passive crown fire?(1=yes,0=no) real(r8), allocatable :: crown_ignite_energy(:) ! heat of crown foliage ignition [kJ/kg] real(r8), allocatable :: initd(:) ! initial seedling density real(r8), allocatable :: seed_suppl(:) ! seeds that come from outside the gridbox. @@ -383,7 +383,7 @@ subroutine Register_PFT(this, fates_params) call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) - name = 'fates_fire_passive_crown_fire' + name = 'fates_fire_crown_fire' call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) @@ -823,9 +823,9 @@ subroutine Receive_PFT(this, fates_params) call fates_params%RetreiveParameterAllocate(name=name, & data=this%crown_kill) - name = 'fates_fire_passive_crown_fire' + name = 'fates_fire_crown_fire' call fates_params%RetreiveParameterAllocate(name=name, & - data=this%passive_crown_fire) + data=this%crown_fire) name = 'fates_fire_crown_ignite_energy' call fates_params%RetreiveParameterAllocate(name=name, & @@ -1740,7 +1740,7 @@ subroutine FatesReportPFTParams(is_master) write(fates_log(),fmt0) 'crown = ',EDPftvarcon_inst%crown write(fates_log(),fmt0) 'bark_scaler = ',EDPftvarcon_inst%bark_scaler write(fates_log(),fmt0) 'crown_kill = ',EDPftvarcon_inst%crown_kill - write(fates_log(),fmt0) 'passive_crown_fire = ',EDPftvarcon_inst%passive_crown_fire + write(fates_log(),fmt0) 'crown_fire = ',EDPftvarcon_inst%crown_fire write(fates_log(),fmt0) 'crown_ignite_energy = ',EDPftvarcon_inst%crown_ignite_energy write(fates_log(),fmt0) 'initd = ',EDPftvarcon_inst%initd write(fates_log(),fmt0) 'seed_suppl = ',EDPftvarcon_inst%seed_suppl diff --git a/parameter_files/fates_params_default.cdl b/parameter_files/fates_params_default.cdl index e6d56f9f07..7789cd2986 100644 --- a/parameter_files/fates_params_default.cdl +++ b/parameter_files/fates_params_default.cdl @@ -133,8 +133,8 @@ variables: fates_fire_crown_kill:units = "NA" ; fates_fire_crown_kill:long_name = "crown resistance to fire (1 = none), EQ 22 in Thonicke et al 2010" ; double fates_fire_crown_fire(fates_pft) ; - fates_fire_crown_fire:units = "logical flag" ; - fates_fire_crown_fire:long_name = "Binary flag for passive crown fire"; + fates_fire_crown_fire:units = "binary flag" ; + fates_fire_crown_fire:long_name = "binary flag for plants susceptible to passive crown fire (1=yes)"; double fates_fire_crown_ignite_energy(fates_pft) ; fates_fire_crown_ignite_energy:units = "kJ/kg" ; fates_fire_crown_ignite_energy:long_name = "spitfire parameter, heat of crown foliage ignition" ; From 3e0d8f7e5061621f43e79acae0fdbf3c4232af44 Mon Sep 17 00:00:00 2001 From: jkshuman Date: Thu, 3 Oct 2019 16:58:05 -0600 Subject: [PATCH 39/72] Update fire/SFMainMod.F90 Co-Authored-By: Samuel Levis --- fire/SFMainMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 6e87575051..16addd4c69 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -937,7 +937,7 @@ subroutine crown_damage ( currentSite ) ! 0.01 = C from Van Wagner 1977 EQ4 for crown of base height 6m, 100% FMC, and FI 2500kW/m passive_crown_FI = (0.01_r8 * height_cbb *EDPftvarcon_inst%crown_ignite_energy(currentCohort%pft))**1.5_r8 - if (passive_crown_FI >= crown_fire_threshold) then ! 200 kW/m = threshold for crown fire potential + if (currentPatch%FI >= crown_fire_threshold) then ! 200 kW/m = threshold for crown fire potential ! Initiation of passive crown fire, EQ 14 Bessie and Johnson 1995 ignite_crown = currentPatch%FI/passive_crown_FI From aeefb59cbbab848cf5ccde95fd1f992b7a8dcacf Mon Sep 17 00:00:00 2001 From: jkshuman Date: Fri, 4 Oct 2019 12:10:02 -0600 Subject: [PATCH 40/72] Update fire/SFMainMod.F90 Co-Authored-By: Samuel Levis --- fire/SFMainMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 16addd4c69..18e293a840 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -925,7 +925,7 @@ subroutine crown_damage ( currentSite ) ! Evaluate for passive crown fire ignition - if (EDPftvarcon_inst%crown_fire(currentCohort%pft) == 1) then + if (EDPftvarcon_inst%crown_fire(currentCohort%pft) > 0.0_r8) then ! Note: crown_ignition_energy to be calculated based on foliar moisture content from FATES-Hydro ! EQ 3 Van Wagner 1977 From 5a6fcc00239f65d39b3d40d5ff292696dab4fd9f Mon Sep 17 00:00:00 2001 From: jkshuman Date: Fri, 4 Oct 2019 12:13:49 -0600 Subject: [PATCH 41/72] Update fire/SFMainMod.F90 Co-Authored-By: Samuel Levis --- fire/SFMainMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 18e293a840..3e9001dc9e 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -944,7 +944,7 @@ subroutine crown_damage ( currentSite ) if (ignite_crown >= 1.0_r8) then currentCohort%passive_crown_fire_flg = 1 ! passive crown fire ignited - currentCohort%fraction_crown_burned = 1.0_r8 + currentCohort%fraction_crown_burned = EDPftvarcon_inst%crown_fire(currentCohort%pft) ! else ! evaluate crown damage based on scorch height endif ! ignite passive crown fire ! else no crown fire today From 9997422e49770535425a756fca9003ae1858265b Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Fri, 4 Oct 2019 13:03:18 -0600 Subject: [PATCH 42/72] Add comment about passive crown fire survival --- fire/SFMainMod.F90 | 3 +++ 1 file changed, 3 insertions(+) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 3e9001dc9e..2120b4d8ab 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -944,6 +944,9 @@ 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" + ! 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 endif ! ignite passive crown fire From 012e363d35184b4ef723d0cb7c9d07a54296da01 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Mon, 7 Oct 2019 13:41:17 -0600 Subject: [PATCH 43/72] Revert zero variables to again occur in EDCohortDynamics --- biogeochem/EDCohortDynamicsMod.F90 | 2 ++ fire/SFMainMod.F90 | 4 +++- 2 files changed, 5 insertions(+), 1 deletion(-) diff --git a/biogeochem/EDCohortDynamicsMod.F90 b/biogeochem/EDCohortDynamicsMod.F90 index 88a0f02eb7..9efacbf416 100644 --- a/biogeochem/EDCohortDynamicsMod.F90 +++ b/biogeochem/EDCohortDynamicsMod.F90 @@ -595,6 +595,8 @@ subroutine zero_cohort(cc_p) currentCohort%leaf_cost = 0._r8 currentcohort%excl_weight = 0._r8 currentcohort%prom_weight = 0._r8 + currentcohort%crownfire_mort = 0._r8 + currentcohort%cambial_mort = 0._r8 currentCohort%c13disc_clm = 0._r8 currentCohort%c13disc_acc = 0._r8 diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 2120b4d8ab..a47e93f903 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -1002,7 +1002,7 @@ subroutine cambial_damage_kill ( currentSite ) if (currentPatch%fire == 1) then currentCohort => currentPatch%tallest; do while(associated(currentCohort)) - currentCohort%cambial_mort = 0.0_r8 + !currentCohort%cambial_mort = 0.0_r8 !testing this removal 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. @@ -1014,6 +1014,8 @@ subroutine cambial_damage_kill ( currentSite ) else if ((currentPatch%tau_l/tau_c) > 0.22_r8) then currentCohort%cambial_mort = (0.563_r8*(currentPatch%tau_l/tau_c)) - 0.125_r8 + else + currentCohort%cambial_mort = 0.0_r8 endif endif endif !trees From 277352c529c62674ca5bb36e78493a2b1d1fc68d Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Wed, 6 Nov 2019 13:10:50 -0700 Subject: [PATCH 44/72] CG scaler on lightning, NF>0 for fire event --- fire/SFMainMod.F90 | 32 +++++++++++++++++--------------- 1 file changed, 17 insertions(+), 15 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index ca1077e5bd..519bc06191 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -438,7 +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 + 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; @@ -481,7 +482,7 @@ subroutine rate_of_spread ( currentSite ) ! Rothermal EQ12= 250 Btu/lb + 1116 Btu/lb * fuel_eff_moist ! conversion of Rothermal (1972) EQ12 in BTU/lb to current kJ/kg ! q_ig in kJ/kg - q_ig = 581.0_r8 +2594.0_r8 * currentPatch%fuel_eff_moist + q_ig = q_dry +2594.0_r8 * currentPatch%fuel_eff_moist ! ---effective heating number--- ! Equation A3 in Thonicke et al. 2010. @@ -661,6 +662,8 @@ subroutine fire_intensity ( currentSite ) !currentPatch%TFC_ROS total fuel consumed by flaming front (kgC/m2) use FatesInterfaceMod, only : hlm_use_spitfire + use EDParamsMod, only : ED_val_nignitions + use FatesConstantsMod, only : years_per_day use SFParamsMod, only : SF_val_fdi_alpha,SF_val_fuel_energy, & SF_val_max_durat, SF_val_durat_slope @@ -670,6 +673,15 @@ subroutine fire_intensity ( currentSite ) real(r8) ROS !m/s real(r8) W !kgBiomass/m2 + real(r8) NF !number of lighting strikes per day per km2 + real(r8),parameter :: CG_strikes = .20_r8 !cloud to ground lightning strikes + !Latham and Williams (2001) + + !NF = number of lighting strikes per day per km2 + NF = ED_val_nignitions * years_per_day * CG_strikes + + ! If there are 15 lightning strikes per year, per km2. (approx from NASA product for S.A.) + ! then there are 15 * 1/365 strikes/km2 each day. currentPatch => currentSite%oldest_patch; @@ -684,7 +696,7 @@ subroutine fire_intensity ( currentSite ) if( hlm_masterproc == itrue ) write(fates_log(),*) 'fire_intensity',currentPatch%fi,W,currentPatch%ROS_front endif !'decide_fire' subroutine shortened and put in here... - if (currentPatch%FI >= fire_threshold) then ! 50kW/m is the threshold for a self-sustaining fire + if (currentPatch%FI >= fire_threshold) .and. (NF > 0._r8) then !50kW/m is 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) @@ -722,30 +734,20 @@ end subroutine fire_intensity subroutine area_burnt ( currentSite ) !***************************************************************** - use EDParamsMod, only : ED_val_nignitions - use FatesConstantsMod, only : years_per_day - type(ed_site_type), intent(inout), target :: currentSite type(ed_patch_type), pointer :: currentPatch 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 - real(r8) NF !number of lightning strikes per day per km2 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 + 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 - !NF = number of lighting strikes per day per km2 - NF = ED_val_nignitions * years_per_day - - ! If there are 15 lightning strikes per year, per km2. (approx from NASA product for S.A.) - ! then there are 15 * 1/365 strikes/km2 each day. - currentPatch => currentSite%oldest_patch; do while(associated(currentPatch)) ! ---initialize patch parameters to zero--- From 51ca42c83892d716bd881e3f23d36aa33ca35fee Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Wed, 6 Nov 2019 14:00:43 -0700 Subject: [PATCH 45/72] Make NF site level variable --- fire/SFMainMod.F90 | 15 ++++++++------- main/EDInitMod.F90 | 3 +++ 2 files changed, 11 insertions(+), 7 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 519bc06191..f4a5bdba80 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -658,30 +658,30 @@ subroutine fire_intensity ( currentSite ) !currentPatch%FI avg fire intensity of flaming front during day. Backward ROS plays no role here. kJ/m/s or kW/m. !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) use FatesInterfaceMod, only : hlm_use_spitfire use EDParamsMod, only : ED_val_nignitions use FatesConstantsMod, only : years_per_day - use SFParamsMod, only : SF_val_fdi_alpha,SF_val_fuel_energy, & + 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 - type(ed_patch_type), pointer :: currentPatch real(r8) ROS !m/s real(r8) W !kgBiomass/m2 - real(r8) NF !number of lighting strikes per day per km2 real(r8),parameter :: CG_strikes = .20_r8 !cloud to ground lightning strikes !Latham and Williams (2001) + currentSite%NF = 0.0_r8 !NF = number of lighting strikes per day per km2 - NF = ED_val_nignitions * years_per_day * CG_strikes + currentSite%NF = ED_val_nignitions * years_per_day * CG_strikes ! If there are 15 lightning strikes per year, per km2. (approx from NASA product for S.A.) - ! then there are 15 * 1/365 strikes/km2 each day. + ! then there are 15 * 1/365 strikes/km2 each day currentPatch => currentSite%oldest_patch; @@ -696,7 +696,7 @@ subroutine fire_intensity ( currentSite ) if( hlm_masterproc == itrue ) write(fates_log(),*) 'fire_intensity',currentPatch%fi,W,currentPatch%ROS_front endif !'decide_fire' subroutine shortened and put in here... - if (currentPatch%FI >= fire_threshold) .and. (NF > 0._r8) then !50kW/m is threshold for a self-sustaining fire + if (currentPatch%FI >= fire_threshold .and. currentSite%NF > 0._r8) then !50kW/m threshold for 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) @@ -741,6 +741,7 @@ subroutine area_burnt ( currentSite ) real(r8) df !distance fire has travelled forward in m real(r8) db !distance fire has travelled backward in m 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 @@ -787,7 +788,7 @@ subroutine area_burnt ( currentSite ) !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 = size_of_fire * NF * currentSite%FDI + AB = size_of_fire * currentSite%NF * currentSite%FDI !frac_burnt in units of m2 here. currentPatch%frac_burnt = min(0.99_r8, AB / km2_to_m2) diff --git a/main/EDInitMod.F90 b/main/EDInitMod.F90 index 9613d79957..105d31581b 100644 --- a/main/EDInitMod.F90 +++ b/main/EDInitMod.F90 @@ -174,6 +174,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%frac_burnt = 0.0_r8 ! burn area read in from external file do el=1,num_elements @@ -255,6 +256,7 @@ subroutine set_site_properties( nsites, sites ) cleafoff = 300 cstat = phen_cstat_notcold ! Leaves are on acc_NI = 0.0_r8 + NF = 0.0_r8 dstat = phen_dstat_moiston ! Leaves are on dleafoff = 300 dleafon = 100 @@ -280,6 +282,7 @@ subroutine set_site_properties( nsites, sites ) sites(s)%dstatus = dstat sites(s)%acc_NI = acc_NI + sites(s)%NF = NF sites(s)%frac_burnt = 0.0_r8 From ace251b1d34c9342ba9824a55905c64634d5213f Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Wed, 6 Nov 2019 15:14:05 -0700 Subject: [PATCH 46/72] Update to NF and fire_intensity subroutine --- fire/SFMainMod.F90 | 5 ----- main/EDTypesMod.F90 | 3 ++- 2 files changed, 2 insertions(+), 6 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index f4a5bdba80..203da064e3 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -718,11 +718,6 @@ subroutine fire_intensity ( currentSite ) currentPatch%fire = 0 ! No fire... :-/ 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. - if( hlm_use_spitfire == ifalse )then - currentPatch%fire = 0 !fudge to turn fire off - endif currentPatch => currentPatch%younger; enddo !end patch loop diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index 6ed7598d63..17a7a3f16a 100644 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -687,7 +687,8 @@ module EDTypesMod real(r8) :: wind ! daily wind in m/min for Spitfire units real(r8) :: acc_ni ! daily nesterov index accumulating over time. real(r8) :: fdi ! daily probability an ignition event will start a fire - real(r8) :: frac_burnt ! fraction of soil burnt in this day. + real(r8) :: NF ! daily ignitions in km2 + real(r8) :: frac_burnt ! fraction of area burnt in this day. ! PLANT HYDRAULICS type(ed_site_hydr_type), pointer :: si_hydr From 20944a616455a032baeccea1bb0da9c00aec9bc5 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Wed, 6 Nov 2019 15:21:36 -0700 Subject: [PATCH 47/72] Update to NF in EDInit --- main/EDInitMod.F90 | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/main/EDInitMod.F90 b/main/EDInitMod.F90 index 105d31581b..7e89c71ade 100644 --- a/main/EDInitMod.F90 +++ b/main/EDInitMod.F90 @@ -256,7 +256,6 @@ subroutine set_site_properties( nsites, sites ) cleafoff = 300 cstat = phen_cstat_notcold ! Leaves are on acc_NI = 0.0_r8 - NF = 0.0_r8 dstat = phen_dstat_moiston ! Leaves are on dleafoff = 300 dleafon = 100 @@ -282,7 +281,7 @@ subroutine set_site_properties( nsites, sites ) sites(s)%dstatus = dstat sites(s)%acc_NI = acc_NI - sites(s)%NF = NF + sites(s)%NF = 0.0_r8 sites(s)%frac_burnt = 0.0_r8 From c79a97712be7aebff8a969d167962b5bf96f4553 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Thu, 7 Nov 2019 15:52:36 -0700 Subject: [PATCH 48/72] Remove zero of NF in SFMain --- fire/SFMainMod.F90 | 1 - 1 file changed, 1 deletion(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 203da064e3..3a67db965b 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -676,7 +676,6 @@ subroutine fire_intensity ( currentSite ) real(r8),parameter :: CG_strikes = .20_r8 !cloud to ground lightning strikes !Latham and Williams (2001) - currentSite%NF = 0.0_r8 !NF = number of lighting strikes per day per km2 currentSite%NF = ED_val_nignitions * years_per_day * CG_strikes From 281bce11ebe825db765cfbf033109ea6986bedca Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Fri, 8 Nov 2019 13:39:26 -0700 Subject: [PATCH 49/72] Update scorch height to cohort level SFMain,EDTypes --- fire/SFMainMod.F90 | 19 ++++++++++--------- main/EDTypesMod.F90 | 2 +- 2 files changed, 11 insertions(+), 10 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 3a67db965b..190112895c 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -805,8 +805,9 @@ end subroutine area_burnt !***************************************************************** subroutine crown_scorching ( currentSite ) !***************************************************************** - !currentPatch%SH !average scorch height for the patch(m) - !currentPatch%FI average fire intensity of flaming front during day. kW/m. + + !currentPatch%FI average fire intensity of flaming front during day. kW/m. + !currentCohort%SH scorch height for the cohort(m) type(ed_site_type), intent(in), target :: currentSite @@ -846,9 +847,10 @@ subroutine crown_scorching ( currentSite ) !This loop weights the scorch height for the contribution of each cohort to the overall biomass. ! does this do anything? I think it might be redundant? RF. - currentPatch%SH = 0.0_r8 + currentCohort => currentPatch%tallest; do while(associated(currentCohort)) + currentCohort%SH = 0.0_r8 if (EDPftvarcon_inst%woody(currentCohort%pft) == 1 & .and. (tree_ag_biomass > 0.0_r8)) then !trees only @@ -862,11 +864,10 @@ subroutine crown_scorching ( currentSite ) !equation 16 in Thonicke et al. 2010 if(write_SF == itrue)then - if ( hlm_masterproc == itrue ) write(fates_log(),*) 'currentPatch%SH',currentPatch%SH,f_ag_bmass + if ( hlm_masterproc == itrue ) write(fates_log(),*) 'currentCohort%SH',currentCohort%SH,f_ag_bmass endif !2/3 Byram (1959) - currentPatch%SH = currentPatch%SH + f_ag_bmass * & - EDPftvarcon_inst%fire_alpha_SH(currentCohort%pft) * (currentPatch%FI**0.667_r8) + currentCohort%SH = EDPftvarcon_inst%fire_alpha_SH(currentCohort%pft) * (currentPatch%FI**0.667_r8) endif !trees only currentCohort=>currentCohort%shorter; @@ -902,16 +903,16 @@ subroutine crown_damage ( currentSite ) if (EDPftvarcon_inst%woody(currentCohort%pft) == 1) then !trees only ! Flames lower than bottom of canopy. ! c%hite is height of cohort - if (currentPatch%SH < (currentCohort%hite-currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft))) then + if (currentCohort%SH < (currentCohort%hite-currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft))) then currentCohort%fraction_crown_burned = 0.0_r8 else ! Flames part of way up canopy. ! Equation 17 in Thonicke et al. 2010. ! flames over bottom of canopy but not over top. - if ((currentCohort%hite > 0.0_r8).and.(currentPatch%SH >= & + if ((currentCohort%hite > 0.0_r8).and.(currentCohort%SH >= & (currentCohort%hite-currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft)))) then - currentCohort%fraction_crown_burned = (currentPatch%SH-currentCohort%hite*(1.0_r8- & + currentCohort%fraction_crown_burned = (currentCohort%SH-currentCohort%hite*(1.0_r8- & EDPftvarcon_inst%crown(currentCohort%pft)))/(currentCohort%hite* & EDPftvarcon_inst%crown(currentCohort%pft)) diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index 17a7a3f16a..4a74971c45 100644 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -342,6 +342,7 @@ module EDTypesMod real(r8) :: dbdeaddt ! time derivative of dead biomass : KgC/year ! FIRE + real(r8) :: sh ! scorch height: m real(r8) :: fraction_crown_burned ! proportion of crown affected by fire:- real(r8) :: cambial_mort ! probability that trees dies due to cambial char ! (conditional on the tree being subjected to the fire) @@ -519,7 +520,6 @@ module EDTypesMod real(r8) :: fi ! average fire intensity of flaming front: kj/m/s or kw/m integer :: fire ! Is there a fire? 1=yes 0=no real(r8) :: fd ! fire duration: mins - real(r8) :: sh ! average scorch height: m ! FIRE EFFECTS real(r8) :: frac_burnt ! fraction burnt: frac gridcell/day From c3f7c230210560a4eb53ab1b5aa1786f6e30432c Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Fri, 8 Nov 2019 14:00:30 -0700 Subject: [PATCH 50/72] Turn off patch level scorch height history var --- main/FatesHistoryInterfaceMod.F90 | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index b1870beee0..fc7b8704d7 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -141,7 +141,7 @@ module FatesHistoryInterfaceMod integer :: ih_TFC_ROS_pa integer :: ih_fire_intensity_pa integer :: ih_fire_area_pa - integer :: ih_scorch_height_pa + ! integer :: ih_scorch_height_pa integer :: ih_fire_fuel_bulkd_pa integer :: ih_fire_fuel_eff_moist_pa integer :: ih_fire_fuel_sav_pa @@ -2438,7 +2438,7 @@ subroutine update_history_dyn(this,nc,nsites,sites) hio_tfc_ros_pa(io_pa) = cpatch%TFC_ROS hio_fire_intensity_pa(io_pa) = cpatch%FI hio_fire_area_pa(io_pa) = cpatch%frac_burnt - hio_scorch_height_pa(io_pa) = cpatch%SH + !hio_scorch_height_pa(io_pa) = cpatch%SH hio_fire_fuel_bulkd_pa(io_pa) = cpatch%fuel_bulkd hio_fire_fuel_eff_moist_pa(io_pa) = cpatch%fuel_eff_moist hio_fire_fuel_sav_pa(io_pa) = cpatch%fuel_sav @@ -3906,10 +3906,10 @@ subroutine define_history_vars(this, initialize_variables) avgflag='A', vtype=patch_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & ivar=ivar, initialize=initialize_variables, index = ih_fire_area_pa ) - call this%set_history_var(vname='SCORCH_HEIGHT', units='m', & - long='spitfire flame height:m', use_default='active', & - avgflag='A', vtype=patch_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & - ivar=ivar, initialize=initialize_variables, index = ih_scorch_height_pa ) + !call this%set_history_var(vname='SCORCH_HEIGHT', units='m', & + ! long='spitfire flame height:m', use_default='active', & + ! avgflag='A', vtype=patch_r8, hlms='CLM:ALM', flushval=0.0_r8, upfreq=1, & + ! ivar=ivar, initialize=initialize_variables, index = ih_scorch_height_pa ) call this%set_history_var(vname='fire_fuel_mef', units='m', & long='spitfire fuel moisture', use_default='active', & From 999fa6bcbae0ead23c7e57284de5f0da6672d32f Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Fri, 8 Nov 2019 14:05:05 -0700 Subject: [PATCH 51/72] Turn off SH history variable --- main/FatesHistoryInterfaceMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/main/FatesHistoryInterfaceMod.F90 b/main/FatesHistoryInterfaceMod.F90 index fc7b8704d7..57d30aa4c3 100644 --- a/main/FatesHistoryInterfaceMod.F90 +++ b/main/FatesHistoryInterfaceMod.F90 @@ -1631,7 +1631,7 @@ subroutine update_history_dyn(this,nc,nsites,sites) hio_effect_wspeed_pa => this%hvars(ih_effect_wspeed_pa)%r81d, & hio_fire_intensity_pa => this%hvars(ih_fire_intensity_pa)%r81d, & hio_fire_area_pa => this%hvars(ih_fire_area_pa)%r81d, & - hio_scorch_height_pa => this%hvars(ih_scorch_height_pa)%r81d, & + ! hio_scorch_height_pa => this%hvars(ih_scorch_height_pa)%r81d, & hio_fire_fuel_bulkd_pa => this%hvars(ih_fire_fuel_bulkd_pa)%r81d, & hio_fire_fuel_eff_moist_pa => this%hvars(ih_fire_fuel_eff_moist_pa)%r81d, & hio_fire_fuel_sav_pa => this%hvars(ih_fire_fuel_sav_pa)%r81d, & From a98f7d65cce2e149d179177ebf30c27d123107f4 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Fri, 8 Nov 2019 14:19:39 -0700 Subject: [PATCH 52/72] Add SH to EDCohortDynamics --- biogeochem/EDCohortDynamicsMod.F90 | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/biogeochem/EDCohortDynamicsMod.F90 b/biogeochem/EDCohortDynamicsMod.F90 index 943c5abf37..b9fb9f003d 100644 --- a/biogeochem/EDCohortDynamicsMod.F90 +++ b/biogeochem/EDCohortDynamicsMod.F90 @@ -534,6 +534,7 @@ subroutine nan_cohort(cc_p) currentCohort%ddbhdt = nan ! time derivative of dbh ! FIRE + currentCohort%sh = nan ! scorch height affecting crown (m) currentCohort%fraction_crown_burned = nan ! proportion of crown affected by fire 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 @@ -579,6 +580,7 @@ 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%sh = 0._r8 currentcohort%fraction_crown_burned = 0._r8 currentCohort%size_class = 1 currentCohort%seed_prod = 0._r8 @@ -1652,7 +1654,8 @@ subroutine copy_cohort( currentCohort,copyc ) n%dhdt = o%dhdt n%ddbhdt = o%ddbhdt - ! FIRE + ! FIRE + n%sh = o%sh n%fraction_crown_burned = o%fraction_crown_burned n%fire_mort = o%fire_mort n%crownfire_mort = o%crownfire_mort From 37b478e241555ec384671d2b3dfdcb0b6e272a7f Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Fri, 8 Nov 2019 14:27:35 -0700 Subject: [PATCH 53/72] Remove SH from fuse patch EDPatchDynamics --- biogeochem/EDPatchDynamicsMod.F90 | 2 -- 1 file changed, 2 deletions(-) diff --git a/biogeochem/EDPatchDynamicsMod.F90 b/biogeochem/EDPatchDynamicsMod.F90 index 1f6501aa44..f3aee3ca67 100644 --- a/biogeochem/EDPatchDynamicsMod.F90 +++ b/biogeochem/EDPatchDynamicsMod.F90 @@ -1980,7 +1980,6 @@ subroutine zero_patch(cp_p) 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%sh = 0.0_r8 ! average scorch height for the patch(m) currentPatch%frac_burnt = 0.0_r8 ! fraction burnt daily currentPatch%burnt_frac_litter(:) = 0.0_r8 currentPatch%btran_ft(:) = 0.0_r8 @@ -2296,7 +2295,6 @@ subroutine fuse_2_patches(csite, dp, rp) rp%fi = (dp%fi*dp%area + rp%fi*rp%area) * inv_sum_area rp%fd = (dp%fd*dp%area + rp%fd*rp%area) * inv_sum_area rp%ros_back = (dp%ros_back*dp%area + rp%ros_back*rp%area) * inv_sum_area - rp%sh = (dp%sh*dp%area + rp%sh*rp%area) * inv_sum_area rp%frac_burnt = (dp%frac_burnt*dp%area + rp%frac_burnt*rp%area) * inv_sum_area rp%burnt_frac_litter(:) = (dp%burnt_frac_litter(:)*dp%area + rp%burnt_frac_litter(:)*rp%area) * inv_sum_area rp%btran_ft(:) = (dp%btran_ft(:)*dp%area + rp%btran_ft(:)*rp%area) * inv_sum_area From 8bc50d2b27674805af5185549883775054746a91 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Fri, 8 Nov 2019 21:40:04 -0700 Subject: [PATCH 54/72] include frac_burnt in FI calculation --- fire/SFMainMod.F90 | 132 +++++++++++++++++++++------------------------ 1 file changed, 62 insertions(+), 70 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 190112895c..f51041429e 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -50,9 +50,8 @@ module SFMainMod public :: charecteristics_of_fuel public :: rate_of_spread public :: ground_fuel_consumption - public :: fire_intensity public :: wind_effect - public :: area_burnt + public :: area_burnt_intensity public :: crown_scorching public :: crown_damage public :: cambial_damage_kill @@ -97,8 +96,7 @@ subroutine fire_model( currentSite, bc_in) call charecteristics_of_fuel(currentSite) call rate_of_spread(currentSite) call ground_fuel_consumption(currentSite) - call fire_intensity(currentSite) - call area_burnt(currentSite) + call area_burnt_intensity(currentSite) call crown_scorching(currentSite) call crown_damage(currentSite) call cambial_damage_kill(currentSite) @@ -651,9 +649,11 @@ subroutine ground_fuel_consumption ( currentSite ) end subroutine ground_fuel_consumption + !***************************************************************** - subroutine fire_intensity ( currentSite ) - !***************************************************************** + subroutine area_burnt_intensity ( currentSite ) + !***************************************************************** + !returns the updated currentPatch%FI value for each patch. !currentPatch%FI avg fire intensity of flaming front during day. Backward ROS plays no role here. kJ/m/s or kW/m. @@ -667,88 +667,55 @@ subroutine fire_intensity ( currentSite ) use FatesConstantsMod, only : years_per_day 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 type(ed_patch_type), pointer :: currentPatch real(r8) ROS !m/s real(r8) W !kgBiomass/m2 + 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 + 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 real(r8),parameter :: CG_strikes = .20_r8 !cloud to ground lightning strikes !Latham and Williams (2001) + ! ---initialize site parameters to zero--- + currentSite%frac_burnt = 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 + currentSite%FDI = 1.0_r8 - exp(-SF_val_fdi_alpha*currentSite%acc_NI) + !NF = number of lighting strikes per day per km2 currentSite%NF = ED_val_nignitions * years_per_day * CG_strikes ! If there are 15 lightning strikes per year, per km2. (approx from NASA product for S.A.) - ! then there are 15 * 1/365 strikes/km2 each day + ! then there are 15 * 1/365 strikes/km2 each day + currentPatch => currentSite%oldest_patch; - do while(associated(currentPatch)) - ROS = currentPatch%ROS_front / 60.0_r8 !m/min to m/sec - W = currentPatch%TFC_ROS / 0.45_r8 !kgC/m2 to kgbiomass/m2 - - !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 - endif - !'decide_fire' subroutine shortened and put in here... - if (currentPatch%FI >= fire_threshold .and. currentSite%NF > 0._r8) then !50kW/m threshold for 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) + ! ---initialize patch parameters to zero--- + currentPatch%frac_burnt = 0.0_r8 + + if (currentSite%NF > 0) then ! Equation 14 in Thonicke et al. 2010 ! fire duration in minutes - currentPatch%FD = (SF_val_max_durat+1.0_r8) / (1.0_r8 + SF_val_max_durat * & exp(SF_val_durat_slope*currentSite%FDI)) - if(write_SF == itrue)then if ( hlm_masterproc == itrue ) write(fates_log(),*) 'fire duration minutes',currentPatch%fd endif !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 - else - currentPatch%fire = 0 ! No fire... :-/ - currentPatch%FD = 0.0_r8 - endif - - currentPatch => currentPatch%younger; - enddo !end patch loop - - end subroutine fire_intensity - + !currentPatch%FD = 60.0_r8 * 24.0_r8 !no minutes in a day - !***************************************************************** - subroutine area_burnt ( currentSite ) - !***************************************************************** - - type(ed_site_type), intent(inout), target :: currentSite - type(ed_patch_type), pointer :: currentPatch - - 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 - 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 - - - ! ---initialize site parameters to zero--- - currentSite%frac_burnt = 0.0_r8 - - currentPatch => currentSite%oldest_patch; - do while(associated(currentPatch)) - ! ---initialize patch parameters to zero--- - currentPatch%frac_burnt = 0.0_r8 - - if (currentPatch%fire == 1) then + ! 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. @@ -784,23 +751,48 @@ subroutine area_burnt ( currentSite ) !AB = m2 per km2 per day AB = size_of_fire * currentSite%NF * currentSite%FDI - !frac_burnt in units of m2 here. - currentPatch%frac_burnt = min(0.99_r8, AB / km2_to_m2) + !frac_burnt + currentPatch%frac_burnt = (min(0.99_r8, AB / km2_to_m2)) * currentPatch%area/area if(write_SF == itrue)then if ( hlm_masterproc == itrue ) write(fates_log(),*) 'frac_burnt',currentPatch%frac_burnt endif - endif - endif! fire - ! convert frac_burnt to % prior to accumulating at site level - currentSite%frac_burnt = currentSite%frac_burnt + currentPatch%frac_burnt * currentPatch%area/area + 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 + + !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 + + if(write_sf == itrue)then + if( hlm_masterproc == itrue ) write(fates_log(),*) 'fire_intensity',currentPatch%fi,W,currentPatch%ROS_front + endif + + !'decide_fire' subroutine + if (currentPatch%FI >= fire_threshold) then !50kW/m threshold for self-sustaining fire + currentPatch%fire = 1 ! Fire... :D + + else + currentPatch%fire = 0 ! No fire... :-/ + currentPatch%FD = 0.0_r8 + currentPatch%frac_burnt = 0.0_r8 + 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 - end subroutine area_burnt + end subroutine area_burnt_intensity + + !***************************************************************** subroutine crown_scorching ( currentSite ) From 9747a7622f4ce5ab4e5c0eb9c4c07b787a3b85b1 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Mon, 18 Nov 2019 11:12:35 -0700 Subject: [PATCH 55/72] Update NF conditional to real SFMain --- fire/SFMainMod.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index f51041429e..2f62d26979 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -680,7 +680,7 @@ subroutine area_burnt_intensity ( currentSite ) real(r8) size_of_fire !in m2 real(r8),parameter :: km2_to_m2 = 1000000.0_r8 !area conversion for square km to square m - real(r8),parameter :: CG_strikes = .20_r8 !cloud to ground lightning strikes + real(r8),parameter :: CG_strikes = 0.20_r8 !cloud to ground lightning strikes !Latham and Williams (2001) ! ---initialize site parameters to zero--- @@ -703,7 +703,7 @@ subroutine area_burnt_intensity ( currentSite ) ! ---initialize patch parameters to zero--- currentPatch%frac_burnt = 0.0_r8 - if (currentSite%NF > 0) then + if (currentSite%NF > 0.0_r8) then ! Equation 14 in Thonicke et al. 2010 ! fire duration in minutes From 27748c56df1ac7f7b76f5afb888fc601f6f9f564 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Tue, 19 Nov 2019 09:50:56 -0700 Subject: [PATCH 56/72] Track all fires. Update energy threshold for fires to zero --- fire/SFMainMod.F90 | 2 +- main/EDTypesMod.F90 | 3 ++- 2 files changed, 3 insertions(+), 2 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 2f62d26979..fa3c2907b9 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -771,7 +771,7 @@ subroutine area_burnt_intensity ( currentSite ) endif !'decide_fire' subroutine - if (currentPatch%FI >= fire_threshold) then !50kW/m threshold for self-sustaining fire + if (currentPatch%FI > fire_threshold) then !track fires greater than kW/m2 energy threshold currentPatch%fire = 1 ! Fire... :D else diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index 4a74971c45..b0fa8b1d88 100644 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -145,7 +145,8 @@ module EDTypesMod integer, parameter, public :: dl_sf = 5 ! array index of dead leaf pool for spitfire (dead grass and dead leaves) integer, parameter, public :: lg_sf = 6 ! array index of live grass pool for spitfire - real(r8), parameter, public :: fire_threshold = 50.0_r8 ! threshold for fires that spread or go out. KWm-2 (Pyne 1986) + !real(r8), parameter, public :: fire_threshold = 50.0_r8 ! threshold for fires that spread or go out. KWm-2 (Pyne 1986) + real(r8), parameter, public :: fire_threshold = 0.0_r8 ! track all fires with energy greater than 0 kW/m2 ! PATCH FUSION real(r8), parameter, public :: force_patchfuse_min_biomass = 0.005_r8 ! min biomass (kg / m2 patch area) below which to force-fuse patches From fe7be5d6951bfe4a28aacbde64bda7b94da0b9b8 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Tue, 19 Nov 2019 14:02:56 -0700 Subject: [PATCH 57/72] Fire threshold of 0.01 --- main/EDTypesMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index b0fa8b1d88..5775d2e541 100644 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -146,7 +146,7 @@ module EDTypesMod integer, parameter, public :: lg_sf = 6 ! array index of live grass pool for spitfire !real(r8), parameter, public :: fire_threshold = 50.0_r8 ! threshold for fires that spread or go out. KWm-2 (Pyne 1986) - real(r8), parameter, public :: fire_threshold = 0.0_r8 ! track all fires with energy greater than 0 kW/m2 + real(r8), parameter, public :: fire_threshold = 0.01_r8 ! track all fires with energy greater than .01 kW/m2 ! PATCH FUSION real(r8), parameter, public :: force_patchfuse_min_biomass = 0.005_r8 ! min biomass (kg / m2 patch area) below which to force-fuse patches From b715f76236a289149a83967505ac6559ad6e0b24 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Thu, 21 Nov 2019 14:22:12 -0700 Subject: [PATCH 58/72] Add calculation for canopy base height for crown fire --- fire/SFMainMod.F90 | 94 +++++++++++++++++++++++++++++++++++++++++----- 1 file changed, 85 insertions(+), 9 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index a47e93f903..e1e7869985 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -891,6 +891,7 @@ subroutine crown_damage ( currentSite ) !currentCohort%fraction_crown_burned is the proportion of crown affected by fire + use SFParamsMod, only : SF_VAL_CWD_FRAC type(ed_site_type), intent(in), target :: currentSite @@ -900,15 +901,95 @@ subroutine crown_damage ( currentSite ) 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) - real(r8) :: ignite_crown ! ratio for ignition of passive crown fire,EQ 14 Bessie & Johnson 1995 - + real(r8) :: ignite_passive_crown ! ratio for ignition of passive crown fire,EQ 14 Bessie & Johnson 1995 + real(r8) :: tree_sap_struct_biomass ! above-ground tree struct and sap biomass in current cohort. kgC/m2 + real(r8) :: leaf_c ! leaf carbon (kg) + real(r8) :: twig_sapw_c ! above-ground twig sap biomass in current cohort. kgC/m2 + real(r8) :: twig_struct_c ! above-ground twig struct biomass in current cohort. kgC/m2 + real(r8) :: crown_fuel_biomass ! biomass of 1 hr fuels (leaves,twigs) in current Patch. kgC/m2 + real(r8) :: crown_fuel_density ! density of crown fuel in current Patch. kgC/m3 + real(r8) :: crown_density_per_m ! density crown fuel per 1m section. + real(r8) :: crown_density ! density crown fuel to find min height for passive crown fire + + real(r8),parameter :: min_passive_crown_fuel = 0.11_r8 !minimum crown fuel density to ignite passive crown fire + currentPatch => currentSite%oldest_patch - do while(associated(currentPatch)) + do while(associated(currentPatch)) + + passive_crown_FI = 0.0_r8 + ignite_passive_crown = 0.0_r8 + tree_sap_struct_biomass = 0.0_r8 + leaf_c = 0.0_r8 ! zero here or in cohort loop? + twig_sapw_c = 0.0_r8 + twig_struct_c = 0.0_r8 + canopy_fuel_biomass = 0.0_r8 + canopy_fuel_density = 0.0_r8 + canopy_density_per_m = 0.0_r8 + if (currentPatch%fire == 1) then + canopy_fuel_density = 0.0_r8 + canopy_fuel_biomass = 0.0_r8 + min_height_fuel = 0.0_r8 + max_height_fuel = 0.0_r8 + currentCohort=>currentPatch%tallest + do while(associated(currentCohort)) + + ! Calculate crown 1hr fuel biomass (leaf, twig sapwood, twig structural biomass) + if (EDPftvarcon_inst%woody(currentCohort%pft) == 1) then !trees only + + crown_depth = currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft) + height_cbb = currentCohort%hite - crown_depth + + if (height_cbb < min_height) then + min_height = height_cbb + endif ! find min height for canopy fuels + + if (currentCohort%hite > max_height) then + max_height = currentCohort%hite + endif ! find max height for canopy fuels + + leaf_c = currentCohort%prt%GetState(leaf_organ, all_carbon_elements) + sapw_c = currentCohort%prt%GetState(sapw_organ, all_carbon_elements) + struct_c = currentCohort%prt%GetState(struct_organ, all_carbon_elements) + + tree_sap_struct_biomass = currentCohort%n * & + (EDPftvarcon_inst%allom_agb_frac(currentCohort%pft)*(sapw_c + struct_c)) + + twig_sapw_c = tree_sap_struct_biomass * SF_VAL_CWD_frac(1) ! find sapw_c for twigs + twig_struct_c = tree_sap_struct_biomass * SF_VAL_CWD_frac(1) ! find struct_c for twigs + + canopy_fuel_biomass = canopy_fuel_biomass = & + (currentCohort%n * leaf_c) + twig_sapw_c + twig_struct_c ! canopy fuel biomass. Patch + + endif !trees only + + currentCohort => currentCohort%shorter; + + enddo !end cohort loop + canopy_fuel_density = canopy_fuel_biomass / currentPatch%area ! kgC/m3 in canopy for patch + + canopy_depth = max_height_fuel - min_height_fuel ! depth across patch + + canopy_density_per_m = canopy_fuel_density / canopy_depth ! density per 1m section of canopy + + ! find height where canopy fuel density in 1 m section is >= to 0.11 kg/m3 + ! loop from min_height to max_height + ! canopy_density = 0.0_r8 + ! canopy_denisty = canopy_density + canopy_density_per_m + ! if canopy_density >= min_passive_crown_fuel + ! then height_base_canopy = this height + ! replace height_cbb with height_base_canopy in EQ passive_crown_FI + ! this prevents small trees without much fuel from starting crown fires in patch + + + + + currentCohort=>currentPatch%tallest + do while(associated(currentCohort)) currentCohort%fraction_crown_burned = 0.0_r8 @@ -919,13 +1000,8 @@ subroutine crown_damage ( currentSite ) ! inst%crown = crown_depth_frac (PFT) crown_depth = currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft) height_cbb = currentCohort%hite - crown_depth - passive_crown_FI = 0.0_r8 - ignite_crown = 0.0_r8 - currentCohort%passive_crown_fire_flg = 0 - - ! Evaluate for passive crown fire ignition - if (EDPftvarcon_inst%crown_fire(currentCohort%pft) > 0.0_r8) then + currentCohort%passive_crown_fire_flg = 0 ! Note: crown_ignition_energy to be calculated based on foliar moisture content from FATES-Hydro ! EQ 3 Van Wagner 1977 From 380666ce8463b334b7dad980853aed142a6e0f12 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Fri, 22 Nov 2019 11:38:52 -0700 Subject: [PATCH 59/72] Update canopy_base_height to patch level for crown fire --- fire/SFMainMod.F90 | 75 ++++++++++++++++++++++++---------------------- 1 file changed, 39 insertions(+), 36 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index e1e7869985..6794d167b8 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -911,7 +911,8 @@ subroutine crown_damage ( currentSite ) real(r8) :: crown_density_per_m ! density crown fuel per 1m section. real(r8) :: crown_density ! density crown fuel to find min height for passive crown fire - real(r8),parameter :: min_passive_crown_fuel = 0.11_r8 !minimum crown fuel density to ignite passive crown fire + real(r8),parameter :: min_density_crown_fuel = 0.11_r8 !min crown fuel density to ignite passive crown fire + real(r8),parameter :: crown_ignite_energy = 3060_r8 ! crown ignition energy currentPatch => currentSite%oldest_patch @@ -931,61 +932,71 @@ subroutine crown_damage ( currentSite ) canopy_fuel_density = 0.0_r8 canopy_fuel_biomass = 0.0_r8 - min_height_fuel = 0.0_r8 max_height_fuel = 0.0_r8 currentCohort=>currentPatch%tallest do while(associated(currentCohort)) - ! Calculate crown 1hr fuel biomass (leaf, twig sapwood, twig structural biomass) + ! Calculate crown 1hr fuel biomass (leaf, twig sapwood, twig structural biomass) if (EDPftvarcon_inst%woody(currentCohort%pft) == 1) then !trees only - crown_depth = currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft) - height_cbb = currentCohort%hite - crown_depth - - if (height_cbb < min_height) then - min_height = height_cbb - endif ! find min height for canopy fuels + crown_depth = currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft) + height_cbb = currentCohort%hite - crown_depth if (currentCohort%hite > max_height) then - max_height = currentCohort%hite + max_height = currentCohort%hite !max height on Patch endif ! find max height for canopy fuels leaf_c = currentCohort%prt%GetState(leaf_organ, all_carbon_elements) sapw_c = currentCohort%prt%GetState(sapw_organ, all_carbon_elements) struct_c = currentCohort%prt%GetState(struct_organ, all_carbon_elements) - tree_sap_struct_biomass = currentCohort%n * & + tree_sap_struct_c = currentCohort%n * & (EDPftvarcon_inst%allom_agb_frac(currentCohort%pft)*(sapw_c + struct_c)) - twig_sapw_c = tree_sap_struct_biomass * SF_VAL_CWD_frac(1) ! find sapw_c for twigs - twig_struct_c = tree_sap_struct_biomass * SF_VAL_CWD_frac(1) ! find struct_c for twigs + twig_sapw_struct_c = tree_sap_struct_c * SF_VAL_CWD_frac(1) ! only 1hr fuel + + canopy_fuel_c = (currentCohort%n * leaf_c) + twig_sapw_struct_c ! canopy fuel kg C. Patch - canopy_fuel_biomass = canopy_fuel_biomass = & - (currentCohort%n * leaf_c) + twig_sapw_c + twig_struct_c ! canopy fuel biomass. Patch + canopy_fuel_biomass = canopy_fuel_c / 0.45_r8 ! canopy fuel biomass + + canopy_fuel_per_m = canopy_fuel_biomass / crown_depth ! biomass per m + + do ih = int(height_cbb), int(currentCohort%hite) + biom_matrix(ih) = biom_matrix(ih) + canopy_fuel_per_m + end do endif !trees only currentCohort => currentCohort%shorter; - enddo !end cohort loop + enddo !end cohort loop - canopy_fuel_density = canopy_fuel_biomass / currentPatch%area ! kgC/m3 in canopy for patch + canopy_fuel_density = canopy_fuel_biomass / (currentPatch%area * canopy_depth) ! kgbiomass/m3 - canopy_depth = max_height_fuel - min_height_fuel ! depth across patch + biom_matrix(:) = biom_matrix(:) / currentPatch%area !kg biomass/m3 - canopy_density_per_m = canopy_fuel_density / canopy_depth ! density per 1m section of canopy + do ih=1,70 !loop from 1m up to 70m to find section where density = 0.11kg/m3 + if (biom_matrix(ih) > min_density_crown_fuel) then + height_base_canopy = float(ih) + exit + end if + end do - ! find height where canopy fuel density in 1 m section is >= to 0.11 kg/m3 - ! loop from min_height to max_height - ! canopy_density = 0.0_r8 - ! canopy_denisty = canopy_density + canopy_density_per_m - ! if canopy_density >= min_passive_crown_fuel - ! then height_base_canopy = this height - ! replace height_cbb with height_base_canopy in EQ passive_crown_FI - ! this prevents small trees without much fuel from starting crown fires in patch + !canopy_bulk_denisty (kg/m3) for Patch + canopy_bulk_denisty = sum(biom_matrix) / (max_height - height_base_canopy) + ! Note: crown_ignition_energy to be calculated based on PFT foliar moisture content from FATES-Hydro + ! Use EDPftvarcon_inst%crown_ignite_energy(currentCohort%pft) and complete as weighted average + ! in place of crown_ignite_energy parameter + ! EQ 3 Van Wagner 1977 + ! h = crown_ignite_energy (kJ/kg), m = foliar moisture content based on dry fuel (%) + ! crown_ignite_energy = 460 + 26 * m + ! Crown fuel ignition potential, EQ 8 Bessie and Johnson 1995 + ! Van Wagner EQ4 FI = (Czh)**3/2 where z=crown base height,h=heat crown ignite energy, FI=fire intensity + ! 0.01 = C from Van Wagner 1977 EQ4 for crown of base height 6m, 100% FMC, and FI 2500kW/m + passive_crown_FI = (0.01_r8 * height_base_canopy * crown_ignite_energy**1.5_r8 currentCohort=>currentPatch%tallest @@ -1003,15 +1014,7 @@ subroutine crown_damage ( currentSite ) currentCohort%passive_crown_fire_flg = 0 - ! Note: crown_ignition_energy to be calculated based on foliar moisture content from FATES-Hydro - ! EQ 3 Van Wagner 1977 - ! h = crown_ignite_energy (kJ/kg), m = foliar moisture content based on dry fuel (%) - ! crown_ignite_energy = 460 + 26 * m - - ! Crown fuel ignition potential, EQ 8 Bessie and Johnson 1995 - ! Van Wagner EQ4 FI = (Czh)**3/2 where z=crown base height,h=heat crown ignite energy, FI=fire intensity - ! 0.01 = C from Van Wagner 1977 EQ4 for crown of base height 6m, 100% FMC, and FI 2500kW/m - passive_crown_FI = (0.01_r8 * height_cbb *EDPftvarcon_inst%crown_ignite_energy(currentCohort%pft))**1.5_r8 + if (currentPatch%FI >= crown_fire_threshold) then ! 200 kW/m = threshold for crown fire potential From 7dfced0c6e9ccd8ea78de69202b7885d15c446ad Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Fri, 22 Nov 2019 17:41:30 -0700 Subject: [PATCH 60/72] Clean up crown scorch SFMain --- fire/SFMainMod.F90 | 20 +++----------------- 1 file changed, 3 insertions(+), 17 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index fa3c2907b9..0d66bbcf1d 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -691,7 +691,7 @@ subroutine area_burnt_intensity ( currentSite ) ! 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) - !NF = number of lighting strikes per day per km2 + !NF = number of lighting strikes per day per km2 scaled by cloud to ground strikes currentSite%NF = ED_val_nignitions * years_per_day * CG_strikes ! If there are 15 lightning strikes per year, per km2. (approx from NASA product for S.A.) @@ -806,7 +806,6 @@ subroutine crown_scorching ( currentSite ) 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 real(r8) :: tree_ag_biomass ! total amount of above-ground tree biomass in patch. kgC/m2 real(r8) :: leaf_c ! leaf carbon [kg] real(r8) :: sapw_c ! sapwood carbon [kg] @@ -817,7 +816,6 @@ subroutine crown_scorching ( currentSite ) do while(associated(currentPatch)) tree_ag_biomass = 0.0_r8 - f_ag_bmass = 0.0_r8 if (currentPatch%fire == 1) then currentCohort => currentPatch%tallest; do while(associated(currentCohort)) @@ -834,11 +832,7 @@ subroutine crown_scorching ( currentSite ) currentCohort=>currentCohort%shorter; - enddo !end cohort loop - - !This loop weights the scorch height for the contribution of each cohort to the overall biomass. - - ! does this do anything? I think it might be redundant? RF. + enddo !end cohort loop currentCohort => currentPatch%tallest; do while(associated(currentCohort)) @@ -846,17 +840,9 @@ subroutine crown_scorching ( currentSite ) if (EDPftvarcon_inst%woody(currentCohort%pft) == 1 & .and. (tree_ag_biomass > 0.0_r8)) then !trees only - leaf_c = currentCohort%prt%GetState(leaf_organ, all_carbon_elements) - sapw_c = currentCohort%prt%GetState(sapw_organ, all_carbon_elements) - struct_c = currentCohort%prt%GetState(struct_organ, all_carbon_elements) - - f_ag_bmass = currentCohort%n * (leaf_c + & - EDPftvarcon_inst%allom_agb_frac(currentCohort%pft)*(sapw_c + struct_c)) & - / tree_ag_biomass - !equation 16 in Thonicke et al. 2010 if(write_SF == itrue)then - if ( hlm_masterproc == itrue ) write(fates_log(),*) 'currentCohort%SH',currentCohort%SH,f_ag_bmass + if ( hlm_masterproc == itrue ) write(fates_log(),*) 'currentCohort%SH',currentCohort%SH endif !2/3 Byram (1959) currentCohort%SH = EDPftvarcon_inst%fire_alpha_SH(currentCohort%pft) * (currentPatch%FI**0.667_r8) From adca1fbfe18f090e054780b89dc12818462d210c Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Sat, 23 Nov 2019 10:45:36 -0700 Subject: [PATCH 61/72] remove extra lines crown scorch SFMain --- fire/SFMainMod.F90 | 17 +++++------------ 1 file changed, 5 insertions(+), 12 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 0d66bbcf1d..f16eb4a4c9 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -828,24 +828,17 @@ subroutine crown_scorching ( currentSite ) tree_ag_biomass = tree_ag_biomass + & currentCohort%n * (leaf_c + & EDPftvarcon_inst%allom_agb_frac(currentCohort%pft)*(sapw_c + struct_c)) - endif !trees only - - currentCohort=>currentCohort%shorter; - - enddo !end cohort loop + - currentCohort => currentPatch%tallest; - do while(associated(currentCohort)) currentCohort%SH = 0.0_r8 - if (EDPftvarcon_inst%woody(currentCohort%pft) == 1 & - .and. (tree_ag_biomass > 0.0_r8)) then !trees only + if (tree_ag_biomass > 0.0_r8)) then - !equation 16 in Thonicke et al. 2010 + !Equation 16 in Thonicke et al. 2010 !2/3 Byram (1959) + currentCohort%SH = EDPftvarcon_inst%fire_alpha_SH(currentCohort%pft) * (currentPatch%FI**0.667_r8) + if(write_SF == itrue)then if ( hlm_masterproc == itrue ) write(fates_log(),*) 'currentCohort%SH',currentCohort%SH endif - !2/3 Byram (1959) - currentCohort%SH = EDPftvarcon_inst%fire_alpha_SH(currentCohort%pft) * (currentPatch%FI**0.667_r8) endif !trees only currentCohort=>currentCohort%shorter; From 6fbba177597fa6b87367a71e5dc3d815964a0adf Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Sat, 23 Nov 2019 10:55:45 -0700 Subject: [PATCH 62/72] Remove extra paren crown scorch --- fire/SFMainMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index f16eb4a4c9..8b9c2d634f 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -831,7 +831,7 @@ subroutine crown_scorching ( currentSite ) currentCohort%SH = 0.0_r8 - if (tree_ag_biomass > 0.0_r8)) then + if (tree_ag_biomass > 0.0_r8) then !Equation 16 in Thonicke et al. 2010 !2/3 Byram (1959) currentCohort%SH = EDPftvarcon_inst%fire_alpha_SH(currentCohort%pft) * (currentPatch%FI**0.667_r8) From bd502508935413015a36898a57e46ac9443befaa Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Sat, 23 Nov 2019 11:08:02 -0700 Subject: [PATCH 63/72] Add end if crown scorch SFMain --- fire/SFMainMod.F90 | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 8b9c2d634f..bad8190910 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -831,14 +831,15 @@ subroutine crown_scorching ( currentSite ) currentCohort%SH = 0.0_r8 - if (tree_ag_biomass > 0.0_r8) then + if (tree_ag_biomass > 0.0_r8) then !Equation 16 in Thonicke et al. 2010 !2/3 Byram (1959) currentCohort%SH = EDPftvarcon_inst%fire_alpha_SH(currentCohort%pft) * (currentPatch%FI**0.667_r8) - if(write_SF == itrue)then - if ( hlm_masterproc == itrue ) write(fates_log(),*) 'currentCohort%SH',currentCohort%SH - endif + if(write_SF == itrue)then + if ( hlm_masterproc == itrue ) write(fates_log(),*) 'currentCohort%SH',currentCohort%SH + endif + endif ! tree biomass endif !trees only currentCohort=>currentCohort%shorter; From 299cbf550f3063216cf626f5d3fdb43f5f1f419c Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Sat, 23 Nov 2019 16:22:33 -0700 Subject: [PATCH 64/72] Clean up crown damage SFMain --- fire/SFMainMod.F90 | 79 +++++++++++++++++++++++----------------------- 1 file changed, 40 insertions(+), 39 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 6794d167b8..46bf20ec2d 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -909,10 +909,14 @@ subroutine crown_damage ( currentSite ) real(r8) :: crown_fuel_biomass ! biomass of 1 hr fuels (leaves,twigs) in current Patch. kgC/m2 real(r8) :: crown_fuel_density ! density of crown fuel in current Patch. kgC/m3 real(r8) :: crown_density_per_m ! density crown fuel per 1m section. - real(r8) :: crown_density ! density crown fuel to find min height for passive crown fire + real(r8) :: crown_density ! density crown fuel to find min height for passive crown fire + + real, dimension(70):: biom_matrix ! matrix to track biomass from bottom to 70m - real(r8),parameter :: min_density_crown_fuel = 0.11_r8 !min crown fuel density to ignite passive crown fire - real(r8),parameter :: crown_ignite_energy = 3060_r8 ! crown ignition energy + real(r8),parameter :: min_density_canopy_fuel = 0.011_r8 !min canopy fuel density (kg/m3) sufficient to + !propogate fire vertically through canopy + !Scott and Reinhardt 2001 RMRS-RP-29 + real(r8),parameter :: crown_ignite_energy = 3060_r8 !crown ignition energy (kJ/kg) Van Wagner 1977 currentPatch => currentSite%oldest_patch @@ -924,16 +928,14 @@ subroutine crown_damage ( currentSite ) leaf_c = 0.0_r8 ! zero here or in cohort loop? twig_sapw_c = 0.0_r8 twig_struct_c = 0.0_r8 - canopy_fuel_biomass = 0.0_r8 - canopy_fuel_density = 0.0_r8 - canopy_density_per_m = 0.0_r8 + canopy_fuel_biomass = 0.0_r8 + canopy_fuel_density = 0.0_r8 + canopy_density_per_m = 0.0_r8 + biom_matrix = 0.0_r8 + max_height_fuel = 0.0_r8 if (currentPatch%fire == 1) then - canopy_fuel_density = 0.0_r8 - canopy_fuel_biomass = 0.0_r8 - max_height_fuel = 0.0_r8 - currentCohort=>currentPatch%tallest do while(associated(currentCohort)) @@ -954,11 +956,11 @@ subroutine crown_damage ( currentSite ) tree_sap_struct_c = currentCohort%n * & (EDPftvarcon_inst%allom_agb_frac(currentCohort%pft)*(sapw_c + struct_c)) - twig_sapw_struct_c = tree_sap_struct_c * SF_VAL_CWD_frac(1) ! only 1hr fuel + twig_sapw_struct_c = tree_sap_struct_c * SF_VAL_CWD_frac(1) !only 1hr fuel - canopy_fuel_c = (currentCohort%n * leaf_c) + twig_sapw_struct_c ! canopy fuel kg C. Patch + canopy_fuel_c = (currentCohort%n * leaf_c) + twig_sapw_struct_c !canopy fuel (kgC) Patch - canopy_fuel_biomass = canopy_fuel_c / 0.45_r8 ! canopy fuel biomass + canopy_fuel_biomass = canopy_fuel_c / 0.45_r8 ! canopy fuel (kg biomass) canopy_fuel_per_m = canopy_fuel_biomass / crown_depth ! biomass per m @@ -972,11 +974,11 @@ subroutine crown_damage ( currentSite ) enddo !end cohort loop - canopy_fuel_density = canopy_fuel_biomass / (currentPatch%area * canopy_depth) ! kgbiomass/m3 + canopy_fuel_density = canopy_fuel_biomass / (currentPatch%area * canopy_depth) !kg biomass/m3 biom_matrix(:) = biom_matrix(:) / currentPatch%area !kg biomass/m3 - do ih=1,70 !loop from 1m up to 70m to find section where density = 0.11kg/m3 + do ih=1,70 !loop from 1m to 70m to find bin with density = 0.011 kg/m3 if (biom_matrix(ih) > min_density_crown_fuel) then height_base_canopy = float(ih) exit @@ -987,60 +989,59 @@ subroutine crown_damage ( currentSite ) canopy_bulk_denisty = sum(biom_matrix) / (max_height - height_base_canopy) ! Note: crown_ignition_energy to be calculated based on PFT foliar moisture content from FATES-Hydro - ! Use EDPftvarcon_inst%crown_ignite_energy(currentCohort%pft) and complete as weighted average + ! Use EDPftvarcon_inst%crown_ignite_energy(currentCohort%pft) and compute weighted average ! in place of crown_ignite_energy parameter + ! EQ 3 Van Wagner 1977 ! h = crown_ignite_energy (kJ/kg), m = foliar moisture content based on dry fuel (%) ! crown_ignite_energy = 460 + 26 * m - ! Crown fuel ignition potential, EQ 8 Bessie and Johnson 1995 - ! Van Wagner EQ4 FI = (Czh)**3/2 where z=crown base height,h=heat crown ignite energy, FI=fire intensity - ! 0.01 = C from Van Wagner 1977 EQ4 for crown of base height 6m, 100% FMC, and FI 2500kW/m - passive_crown_FI = (0.01_r8 * height_base_canopy * crown_ignite_energy**1.5_r8 + ! Crown fuel ignition potential, EQ 8 Bessie and Johnson 1995, EQ 4 Van Wagner 1977 + ! FI = (Czh)**3/2 where z=canopy base height,h=heat crown ignite energy, FI=fire intensity + ! 0.01 = C from Van Wagner 1977 EQ4 for canopy base height 6m, 100% FMC, and FI 2500kW/m + passive_crown_FI = (0.01_r8 * height_base_canopy * crown_ignite_energy)**1.5_r8 currentCohort=>currentPatch%tallest 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) + ! height_cbb = clear branch bole height at base of crown (m) ! inst%crown = crown_depth_frac (PFT) - crown_depth = currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft) - height_cbb = currentCohort%hite - crown_depth + crown_depth = currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft) + height_cbb = currentCohort%hite - crown_depth currentCohort%passive_crown_fire_flg = 0 - - - + if (currentPatch%FI >= crown_fire_threshold) then ! 200 kW/m = threshold for crown fire potential ! Initiation of passive crown fire, EQ 14 Bessie and Johnson 1995 - ignite_crown = currentPatch%FI/passive_crown_FI + ignite_passive_crown = currentPatch%FI/passive_crown_FI - 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" - ! Johnson, E.A. 1992 Fire and Veg Dynamics: North American boreal forest. Cambridge Press - currentCohort%fraction_crown_burned = EDPftvarcon_inst%crown_fire(currentCohort%pft) + if (ignite_passive_crown >= 1.0_r8) then + currentCohort%passive_crown_fire_flg = 1 !passive canopy fuels ignited for vertical spread + + !evaluate active crown fire conditions + ! else ! evaluate crown damage based on scorch height - endif ! ignite passive crown fire + endif ! ignite crown fire ! else no crown fire today endif ! crown fire intensity ! else ! not crown fire plant endif ! evaluate passive crown fire ! For surface fires, are flames in the canopy? - ! height_cbb is clear branch bole height or height of bottom of 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 + if (currentCohort%hite > 0.0_r8 .and. currentPatch%SH > height_cbb) + + !add active flag !.and. currentCohort%active_crown_fire_flg == 0) then - currentCohort%fraction_crown_burned = min(1.0_r8, ((currentPatch%SH - height_cbb)/crown_depth)) + currentCohort%fraction_crown_burned = min(1.0_r8, ((currentPatch%SH - height_cbb)/crown_depth)) endif !SH frac crown burnt calculation ! Check for strange values. currentCohort%fraction_crown_burned = min(1.0_r8, max(0.0_r8,currentCohort%fraction_crown_burned)) From f2442792fea4ea9fce47443e2c1a09f862fee9a3 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Sat, 23 Nov 2019 17:51:20 -0700 Subject: [PATCH 65/72] Update passive fire fuels in crown damage --- fire/SFMainMod.F90 | 66 ++++++++++++++++++++++++---------------------- 1 file changed, 35 insertions(+), 31 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 46bf20ec2d..192a40db05 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -900,16 +900,20 @@ subroutine crown_damage ( currentSite ) 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) - real(r8) :: ignite_passive_crown ! ratio for ignition of passive crown fire,EQ 14 Bessie & Johnson 1995 - real(r8) :: tree_sap_struct_biomass ! above-ground tree struct and sap biomass in current cohort. kgC/m2 - real(r8) :: leaf_c ! leaf carbon (kg) - real(r8) :: twig_sapw_c ! above-ground twig sap biomass in current cohort. kgC/m2 - real(r8) :: twig_struct_c ! above-ground twig struct biomass in current cohort. kgC/m2 - real(r8) :: crown_fuel_biomass ! biomass of 1 hr fuels (leaves,twigs) in current Patch. kgC/m2 - real(r8) :: crown_fuel_density ! density of crown fuel in current Patch. kgC/m3 - real(r8) :: crown_density_per_m ! density crown fuel per 1m section. - real(r8) :: crown_density ! density crown fuel to find min height for passive crown fire + real(r8) :: max_height ! max cohort on patch (m) + real(r8) :: passive_crown_FI ! fire intensity for ignition from passive canopy fuel (kW/m), EQ 8 + real(r8) :: ignite_passive_crown ! ratio for ignition from passive canopy fuel,EQ 14 Bessie & Johnson 1995 + real(r8) :: tree_sapw_struct_c ! above-ground tree struct and sap biomass in cohort (kgC) + real(r8) :: leaf_c ! leaf carbon (kgC) + real(r8) :: sapw_c ! sapwood carbon (kgC) + real(r8) :: struct_c ! structure carbon (kgC) + real(r8) :: twig_sapw_struct_c ! above-ground twig sap and struct in cohort (kgC) + real(r8) :: crown_fuel_c ! biomass of 1 hr fuels (leaves,twigs) in cohort (kg C) + real(r8) :: crown_fuel_biomass ! biomass of crown fuel in cohort (kg biomass) + real(r8) :: crown_fuel_per_m ! crown fuel per 1m section in cohort + real(r8) :: canopy_bulk_density ! density of canopy fuel on patch + real(r8) :: height_base_canopy ! lowest height of fuels to carry fire in crown + integer :: ih ! counter real, dimension(70):: biom_matrix ! matrix to track biomass from bottom to 70m @@ -917,6 +921,7 @@ subroutine crown_damage ( currentSite ) !propogate fire vertically through canopy !Scott and Reinhardt 2001 RMRS-RP-29 real(r8),parameter :: crown_ignite_energy = 3060_r8 !crown ignition energy (kJ/kg) Van Wagner 1977 + currentPatch => currentSite%oldest_patch @@ -924,15 +929,18 @@ subroutine crown_damage ( currentSite ) passive_crown_FI = 0.0_r8 ignite_passive_crown = 0.0_r8 - tree_sap_struct_biomass = 0.0_r8 + tree_sapw_struct_c = 0.0_r8 leaf_c = 0.0_r8 ! zero here or in cohort loop? - twig_sapw_c = 0.0_r8 - twig_struct_c = 0.0_r8 - canopy_fuel_biomass = 0.0_r8 - canopy_fuel_density = 0.0_r8 - canopy_density_per_m = 0.0_r8 + sapw_c = 0.0_r8 + struct_c = 0.0_r8 + twig_sapw_struct_c = 0.0_r8 + crown_fuel_c = 0.0_r8 + crown_fuel_biomass = 0.0_r8 + crown_fuel_per_m = 0.0_r8 biom_matrix = 0.0_r8 - max_height_fuel = 0.0_r8 + canopy_bulk_density = 0.0_r8 + max_height = 0.0_r8 + height_base_canopy = 0.0_r8 if (currentPatch%fire == 1) then @@ -953,19 +961,19 @@ subroutine crown_damage ( currentSite ) sapw_c = currentCohort%prt%GetState(sapw_organ, all_carbon_elements) struct_c = currentCohort%prt%GetState(struct_organ, all_carbon_elements) - tree_sap_struct_c = currentCohort%n * & + tree_sapw_struct_c = currentCohort%n * & (EDPftvarcon_inst%allom_agb_frac(currentCohort%pft)*(sapw_c + struct_c)) - twig_sapw_struct_c = tree_sap_struct_c * SF_VAL_CWD_frac(1) !only 1hr fuel + twig_sapw_struct_c = tree_sapw_struct_c * SF_VAL_CWD_frac(1) !only 1hr fuel - canopy_fuel_c = (currentCohort%n * leaf_c) + twig_sapw_struct_c !canopy fuel (kgC) Patch + crown_fuel_c = (currentCohort%n * leaf_c) + twig_sapw_struct_c !crown fuel (kgC) - canopy_fuel_biomass = canopy_fuel_c / 0.45_r8 ! canopy fuel (kg biomass) + crown_fuel_biomass = crown_fuel_c / 0.45_r8 ! crown fuel (kg biomass) - canopy_fuel_per_m = canopy_fuel_biomass / crown_depth ! biomass per m + crown_fuel_per_m = crown_fuel_biomass / crown_depth ! kg biomass per m do ih = int(height_cbb), int(currentCohort%hite) - biom_matrix(ih) = biom_matrix(ih) + canopy_fuel_per_m + biom_matrix(ih) = biom_matrix(ih) + crown_fuel_per_m end do endif !trees only @@ -974,19 +982,17 @@ subroutine crown_damage ( currentSite ) enddo !end cohort loop - canopy_fuel_density = canopy_fuel_biomass / (currentPatch%area * canopy_depth) !kg biomass/m3 - biom_matrix(:) = biom_matrix(:) / currentPatch%area !kg biomass/m3 do ih=1,70 !loop from 1m to 70m to find bin with density = 0.011 kg/m3 - if (biom_matrix(ih) > min_density_crown_fuel) then + if (biom_matrix(ih) > min_density_canopy_fuel) then height_base_canopy = float(ih) exit end if end do !canopy_bulk_denisty (kg/m3) for Patch - canopy_bulk_denisty = sum(biom_matrix) / (max_height - height_base_canopy) + canopy_bulk_density = sum(biom_matrix) / (max_height - height_base_canopy) ! Note: crown_ignition_energy to be calculated based on PFT foliar moisture content from FATES-Hydro ! Use EDPftvarcon_inst%crown_ignite_energy(currentCohort%pft) and compute weighted average @@ -997,7 +1003,7 @@ subroutine crown_damage ( currentSite ) ! crown_ignite_energy = 460 + 26 * m ! Crown fuel ignition potential, EQ 8 Bessie and Johnson 1995, EQ 4 Van Wagner 1977 - ! FI = (Czh)**3/2 where z=canopy base height,h=heat crown ignite energy, FI=fire intensity + ! FI = (Czh)**3/2 where z=canopy base height,h=heat of crown ignite energy, FI=fire intensity ! 0.01 = C from Van Wagner 1977 EQ4 for canopy base height 6m, 100% FMC, and FI 2500kW/m passive_crown_FI = (0.01_r8 * height_base_canopy * crown_ignite_energy)**1.5_r8 @@ -1030,14 +1036,12 @@ subroutine crown_damage ( currentSite ) endif ! ignite crown fire ! else no crown fire today endif ! crown fire intensity - ! else ! not crown fire plant - endif ! evaluate passive crown fire ! 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) + if (currentCohort%hite > 0.0_r8 .and. currentPatch%SH > height_cbb) then !add active flag !.and. currentCohort%active_crown_fire_flg == 0) then From 9a648ea5b92e4944622f6ad85c9a8713c9ac96d3 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Mon, 25 Nov 2019 10:32:38 -0700 Subject: [PATCH 66/72] Maintain 50kW/m threshold for tracking fires --- main/EDTypesMod.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/main/EDTypesMod.F90 b/main/EDTypesMod.F90 index 5775d2e541..53cebaee96 100644 --- a/main/EDTypesMod.F90 +++ b/main/EDTypesMod.F90 @@ -145,8 +145,8 @@ module EDTypesMod integer, parameter, public :: dl_sf = 5 ! array index of dead leaf pool for spitfire (dead grass and dead leaves) integer, parameter, public :: lg_sf = 6 ! array index of live grass pool for spitfire - !real(r8), parameter, public :: fire_threshold = 50.0_r8 ! threshold for fires that spread or go out. KWm-2 (Pyne 1986) - real(r8), parameter, public :: fire_threshold = 0.01_r8 ! track all fires with energy greater than .01 kW/m2 + real(r8), parameter, public :: fire_threshold = 50.0_r8 ! threshold for fires that spread or go out. KWm-2 (Pyne 1986) + !real(r8), parameter, public :: fire_threshold = 0.01_r8 ! track all fires with energy greater than .01 kW/m2 ! PATCH FUSION real(r8), parameter, public :: force_patchfuse_min_biomass = 0.005_r8 ! min biomass (kg / m2 patch area) below which to force-fuse patches From 2ace55f5f3767d5a28a45682f05085c3df8f257f Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Mon, 25 Nov 2019 10:46:38 -0700 Subject: [PATCH 67/72] Add Van Wagner ref crown scrorch SFMain --- fire/SFMainMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index bad8190910..64e3c6cd16 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -833,7 +833,7 @@ subroutine crown_scorching ( currentSite ) currentCohort%SH = 0.0_r8 if (tree_ag_biomass > 0.0_r8) then - !Equation 16 in Thonicke et al. 2010 !2/3 Byram (1959) + !Equation 16 in Thonicke et al. 2010 !Van Wagner 1973 EQ8 !2/3 Byram (1959) currentCohort%SH = EDPftvarcon_inst%fire_alpha_SH(currentCohort%pft) * (currentPatch%FI**0.667_r8) if(write_SF == itrue)then From 38388ade4804148840d3c1694ac2488d12383b83 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Mon, 25 Nov 2019 11:16:44 -0700 Subject: [PATCH 68/72] Add ref and zero patch vars in area_burnt_intensity --- fire/SFMainMod.F90 | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 64e3c6cd16..e8dcc83260 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -701,7 +701,10 @@ subroutine area_burnt_intensity ( currentSite ) currentPatch => currentSite%oldest_patch; do while(associated(currentPatch)) ! ---initialize patch parameters to zero--- + currentPatch%fire = 0 + currentPatch%FD = 0.0_r8 currentPatch%frac_burnt = 0.0_r8 + if (currentSite%NF > 0.0_r8) then @@ -763,6 +766,7 @@ subroutine area_burnt_intensity ( currentSite ) ROS = currentPatch%ROS_front / 60.0_r8 !m/min to m/sec 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 From ea6ecd974a42dacd3903707f0f86b1b55975531a Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Mon, 25 Nov 2019 14:47:38 -0700 Subject: [PATCH 69/72] Update fire parameter names foliar_moisture, crown_fire_flag --- fire/SFMainMod.F90 | 5 ++--- main/EDPftvarcon.F90 | 12 ++++++------ parameter_files/fates_params_default.cdl | 16 ++++++++-------- 3 files changed, 16 insertions(+), 17 deletions(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 50d6669612..644efde37a 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -971,7 +971,7 @@ subroutine crown_damage ( currentSite ) canopy_bulk_density = sum(biom_matrix) / (max_height - height_base_canopy) ! Note: crown_ignition_energy to be calculated based on PFT foliar moisture content from FATES-Hydro - ! Use EDPftvarcon_inst%crown_ignite_energy(currentCohort%pft) and compute weighted average + ! Use EDPftvarcon_inst%foliar_moisture(currentCohort%pft) and compute weighted PFT average ! in place of crown_ignite_energy parameter ! EQ 3 Van Wagner 1977 @@ -996,7 +996,7 @@ subroutine crown_damage ( currentSite ) crown_depth = currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft) height_cbb = currentCohort%hite - crown_depth - currentCohort%passive_crown_fire_flg = 0 ! does tand have canopy fuels for active crown fire? + currentCohort%passive_crown_fire_flg = 0 !does patch have canopy fuels for active crown fire? if (currentPatch%FI >= crown_fire_threshold) then ! 200 kW/m = threshold for crown fire potential @@ -1008,7 +1008,6 @@ subroutine crown_damage ( currentSite ) !evaluate active crown fire conditions - ! else ! evaluate crown damage based on scorch height endif ! ignite crown fire ! else no crown fire today endif ! crown fire intensity diff --git a/main/EDPftvarcon.F90 b/main/EDPftvarcon.F90 index 359c66f5ce..36143624f4 100644 --- a/main/EDPftvarcon.F90 +++ b/main/EDPftvarcon.F90 @@ -49,8 +49,8 @@ module EDPftvarcon ! that is occupied by crown. For fire model. real(r8), allocatable :: bark_scaler(:) ! scaler from dbh to bark thickness. For fire model. real(r8), allocatable :: crown_kill(:) ! crown resistance to fire. (1 = none) For fire model. - real(r8), allocatable :: crown_fire(:) ! Is plant susceptible to passive crown fire?(1=yes,0=no) - real(r8), allocatable :: crown_ignite_energy(:) ! heat of crown foliage ignition [kJ/kg] + real(r8), allocatable :: active_crown_fire(:) ! Is plant susceptible to active crown fire?(1=yes,0=no) + real(r8), allocatable :: foliar_moisture(:) ! foliar moisture content from dry fuel [%] real(r8), allocatable :: initd(:) ! initial seedling density real(r8), allocatable :: seed_suppl(:) ! seeds that come from outside the gridbox. real(r8), allocatable :: BB_slope(:) ! ball berry slope parameter @@ -383,11 +383,11 @@ subroutine Register_PFT(this, fates_params) call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) - name = 'fates_fire_crown_fire' + name = 'fates_fire_active_crown_fire' call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) - name = 'fates_fire_crown_ignite_energy' + name = 'fates_fire_foliar_moisture' call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) @@ -823,11 +823,11 @@ subroutine Receive_PFT(this, fates_params) call fates_params%RetreiveParameterAllocate(name=name, & data=this%crown_kill) - name = 'fates_fire_crown_fire' + name = 'fates_fire_active_crown_fire' call fates_params%RetreiveParameterAllocate(name=name, & data=this%crown_fire) - name = 'fates_fire_crown_ignite_energy' + name = 'fates_fire_foliar_moisture' call fates_params%RetreiveParameterAllocate(name=name, & data=this%crown_ignite_energy) diff --git a/parameter_files/fates_params_default.cdl b/parameter_files/fates_params_default.cdl index 7789cd2986..c29c3c60dd 100644 --- a/parameter_files/fates_params_default.cdl +++ b/parameter_files/fates_params_default.cdl @@ -132,12 +132,12 @@ variables: double fates_fire_crown_kill(fates_pft) ; fates_fire_crown_kill:units = "NA" ; fates_fire_crown_kill:long_name = "crown resistance to fire (1 = none), EQ 22 in Thonicke et al 2010" ; - double fates_fire_crown_fire(fates_pft) ; - fates_fire_crown_fire:units = "binary flag" ; - fates_fire_crown_fire:long_name = "binary flag for plants susceptible to passive crown fire (1=yes)"; - double fates_fire_crown_ignite_energy(fates_pft) ; - fates_fire_crown_ignite_energy:units = "kJ/kg" ; - fates_fire_crown_ignite_energy:long_name = "spitfire parameter, heat of crown foliage ignition" ; + double fates_fire_active_crown_fire(fates_pft) ; + fates_fire_active_crown_fire:units = "logical flag" ; + fates_fire_active_crown_fire:long_name = "Binary flag for active crown fire (1=yes)"; + double fates_fire_foliar_moisture(fates_pft) ; + fates_fire_foliar_moisture:units = "%" ; + fates_fire_foliar_moisture:long_name = "spitfire parameter, foliar moisture content based on dry fuel (%)" ; double fates_fr_fcel(fates_pft) ; fates_fr_fcel:units = "fraction" ; fates_fr_fcel:long_name = "Fine root litter cellulose fraction" ; @@ -724,9 +724,9 @@ data: fates_fire_crown_kill = 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775, 0.775 ; - fates_fire_passive_crown_fire = 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 ; + fates_fire_active_crown_fire = 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 ; - fates_fire_crown_ignite_energy = 3060, 3060, 3060, 3060, 3060, 3060, 3060, 3060, 3060, 3060, 3060, 3060 ; + fates_fire_foliar_moisture = 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100 ; fates_fr_fcel = 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5 ; From 8b9565e65945d989b294902a5825c1dec0159831 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Mon, 25 Nov 2019 14:56:42 -0700 Subject: [PATCH 70/72] Update var namaes for crown fire and foliar moisture --- main/EDPftvarcon.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/main/EDPftvarcon.F90 b/main/EDPftvarcon.F90 index 36143624f4..0425338043 100644 --- a/main/EDPftvarcon.F90 +++ b/main/EDPftvarcon.F90 @@ -825,11 +825,11 @@ subroutine Receive_PFT(this, fates_params) name = 'fates_fire_active_crown_fire' call fates_params%RetreiveParameterAllocate(name=name, & - data=this%crown_fire) + data=this%active_crown_fire) name = 'fates_fire_foliar_moisture' call fates_params%RetreiveParameterAllocate(name=name, & - data=this%crown_ignite_energy) + data=this%foliar_moisture) name = 'fates_recruit_initd' call fates_params%RetreiveParameterAllocate(name=name, & From 23ccb3f4a7b4c48115b238b9792ad5374edccc77 Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Mon, 25 Nov 2019 15:04:17 -0700 Subject: [PATCH 71/72] SH as cohort level in crown damage/SFMain --- fire/SFMainMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 644efde37a..410c62bb22 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -1021,7 +1021,7 @@ subroutine crown_damage ( currentSite ) !add active flag !.and. currentCohort%active_crown_fire_flg == 0) then - currentCohort%fraction_crown_burned = min(1.0_r8, ((currentPatch%SH - height_cbb)/crown_depth)) + currentCohort%fraction_crown_burned = min(1.0_r8, ((currentCohort%SH - height_cbb)/crown_depth)) endif !SH frac crown burnt calculation ! Check for strange values. currentCohort%fraction_crown_burned = min(1.0_r8, max(0.0_r8,currentCohort%fraction_crown_burned)) From 339aca2172839cc3aed5427db8af6982e8170d9a Mon Sep 17 00:00:00 2001 From: Jacquelyn Shuman Date: Tue, 26 Nov 2019 12:13:08 -0700 Subject: [PATCH 72/72] Clean up patch level crown fire calculations --- biogeochem/EDCohortDynamicsMod.F90 | 1 - fire/SFMainMod.F90 | 85 +++++++++++++----------- main/EDPftvarcon.F90 | 14 +--- parameter_files/fates_params_default.cdl | 9 +-- 4 files changed, 49 insertions(+), 60 deletions(-) diff --git a/biogeochem/EDCohortDynamicsMod.F90 b/biogeochem/EDCohortDynamicsMod.F90 index 4cecf9879c..a12406c907 100644 --- a/biogeochem/EDCohortDynamicsMod.F90 +++ b/biogeochem/EDCohortDynamicsMod.F90 @@ -535,7 +535,6 @@ subroutine nan_cohort(cc_p) ! FIRE - currentCohort%passive_crown_fire_flg = nan ! flag to indicate passive crown fire ignition currentCohort%sh = nan ! scorch height affecting crown (m) currentCohort%fraction_crown_burned = nan ! proportion of crown affected by fire currentCohort%cambial_mort = nan ! probability that trees dies due to cambial char P&R (1986) diff --git a/fire/SFMainMod.F90 b/fire/SFMainMod.F90 index 410c62bb22..17401a0f85 100644 --- a/fire/SFMainMod.F90 +++ b/fire/SFMainMod.F90 @@ -841,7 +841,7 @@ subroutine crown_scorching ( currentSite ) currentCohort%SH = 0.0_r8 if (tree_ag_biomass > 0.0_r8) then - !Equation 16 in Thonicke et al. 2010 !Van Wagner 1973 EQ8 !2/3 Byram (1959) + !Equation 16 in Thonicke et al. 2010 !Van Wagner 1973 EQ8 based on 2/3 Byram (1959) currentCohort%SH = EDPftvarcon_inst%fire_alpha_SH(currentCohort%pft) * (currentPatch%FI**0.667_r8) if(write_SF == itrue)then @@ -889,39 +889,44 @@ subroutine crown_damage ( currentSite ) real(r8) :: crown_fuel_per_m ! crown fuel per 1m section in cohort real(r8) :: canopy_bulk_density ! density of canopy fuel on patch real(r8) :: height_base_canopy ! lowest height of fuels to carry fire in crown - integer :: ih ! counter + integer :: ih ! counter + integer :: passive_canopy_fuel_flg ! flag if canopy fuel true for vertical spread real, dimension(70):: biom_matrix ! matrix to track biomass from bottom to 70m real(r8),parameter :: min_density_canopy_fuel = 0.011_r8 !min canopy fuel density (kg/m3) sufficient to !propogate fire vertically through canopy !Scott and Reinhardt 2001 RMRS-RP-29 + real(r8),parameter :: crown_ignite_energy = 3060_r8 !crown ignition energy (kJ/kg) Van Wagner 1977 currentPatch => currentSite%oldest_patch do while(associated(currentPatch)) - + + !zero Patch level variables passive_crown_FI = 0.0_r8 - ignite_passive_crown = 0.0_r8 - tree_sapw_struct_c = 0.0_r8 - leaf_c = 0.0_r8 ! zero here or in cohort loop? - sapw_c = 0.0_r8 - struct_c = 0.0_r8 - twig_sapw_struct_c = 0.0_r8 - crown_fuel_c = 0.0_r8 - crown_fuel_biomass = 0.0_r8 - crown_fuel_per_m = 0.0_r8 + ignite_passive_crown = 0.0_r8 biom_matrix = 0.0_r8 canopy_bulk_density = 0.0_r8 max_height = 0.0_r8 height_base_canopy = 0.0_r8 - + if (currentPatch%fire == 1) then currentCohort=>currentPatch%tallest - do while(associated(currentCohort)) + do while(associated(currentCohort)) + + !zero cohort level variables + tree_sapw_struct_c = 0.0_r8 + leaf_c = 0.0_r8 + sapw_c = 0.0_r8 + struct_c = 0.0_r8 + twig_sapw_struct_c = 0.0_r8 + crown_fuel_c = 0.0_r8 + crown_fuel_biomass = 0.0_r8 + crown_fuel_per_m = 0.0_r8 ! Calculate crown 1hr fuel biomass (leaf, twig sapwood, twig structural biomass) if (EDPftvarcon_inst%woody(currentCohort%pft) == 1) then !trees only @@ -929,9 +934,10 @@ subroutine crown_damage ( currentSite ) crown_depth = currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft) height_cbb = currentCohort%hite - crown_depth + !find patch max height for stand canopy fuel if (currentCohort%hite > max_height) then - max_height = currentCohort%hite !max height on Patch - endif ! find max height for canopy fuels + max_height = currentCohort%hite + endif leaf_c = currentCohort%prt%GetState(leaf_organ, all_carbon_elements) sapw_c = currentCohort%prt%GetState(sapw_organ, all_carbon_elements) @@ -944,10 +950,12 @@ subroutine crown_damage ( currentSite ) crown_fuel_c = (currentCohort%n * leaf_c) + twig_sapw_struct_c !crown fuel (kgC) - crown_fuel_biomass = crown_fuel_c / 0.45_r8 ! crown fuel (kg biomass) + crown_fuel_biomass = crown_fuel_c / 0.45_r8 ! crown fuel (kg biomass) - crown_fuel_per_m = crown_fuel_biomass / crown_depth ! kg biomass per m + crown_fuel_per_m = crown_fuel_biomass / crown_depth ! kg biomass per m + !sort crown fuel into bins from bottom to top of crown + !accumulate across cohorts to find density within canopy 1m sections do ih = int(height_cbb), int(currentCohort%hite) biom_matrix(ih) = biom_matrix(ih) + crown_fuel_per_m end do @@ -960,7 +968,9 @@ subroutine crown_damage ( currentSite ) biom_matrix(:) = biom_matrix(:) / currentPatch%area !kg biomass/m3 - do ih=1,70 !loop from 1m to 70m to find bin with density = 0.011 kg/m3 + !loop from 1m to 70m to find bin with total density = 0.011 kg/m3 + !min canopy fuel density to propogate fire vertically in canopy across patch + do ih=1,70 if (biom_matrix(ih) > min_density_canopy_fuel) then height_base_canopy = float(ih) exit @@ -971,7 +981,8 @@ subroutine crown_damage ( currentSite ) canopy_bulk_density = sum(biom_matrix) / (max_height - height_base_canopy) ! Note: crown_ignition_energy to be calculated based on PFT foliar moisture content from FATES-Hydro - ! Use EDPftvarcon_inst%foliar_moisture(currentCohort%pft) and compute weighted PFT average + ! or create foliar_moisture based on BTRAN + ! Use foliar_moisture(currentCohort%pft) and compute weighted PFT average with EQ3 Van Wagner 1977 ! in place of crown_ignite_energy parameter ! EQ 3 Van Wagner 1977 @@ -982,7 +993,18 @@ subroutine crown_damage ( currentSite ) ! FI = (Czh)**3/2 where z=canopy base height,h=heat of crown ignite energy, FI=fire intensity ! 0.01 = C from Van Wagner 1977 EQ4 for canopy base height 6m, 100% FMC, and FI 2500kW/m passive_crown_FI = (0.01_r8 * height_base_canopy * crown_ignite_energy)**1.5_r8 + + passive_canopy_fuel_flg = 0 !does patch have canopy fuels for vertical spread? + + ! Initiation of passive crown fire, EQ 14a Bessie and Johnson 1995 + ! Are the canopy fuels in the stand large enough to support vertical spread of fire? + ignite_passive_crown = currentPatch%FI/passive_crown_FI + + if (ignite_passive_crown >= 1.0_r8) then + passive_canopy_fuel_flg = 1 !enough passive canopy fuels for vertical spread + endif + !evaluate active crown fire conditions currentCohort=>currentPatch%tallest @@ -996,32 +1018,15 @@ subroutine crown_damage ( currentSite ) crown_depth = currentCohort%hite*EDPftvarcon_inst%crown(currentCohort%pft) height_cbb = currentCohort%hite - crown_depth - currentCohort%passive_crown_fire_flg = 0 !does patch have canopy fuels for active crown fire? - - if (currentPatch%FI >= crown_fire_threshold) then ! 200 kW/m = threshold for crown fire potential - - ! Initiation of passive crown fire, EQ 14 Bessie and Johnson 1995 - ignite_passive_crown = currentPatch%FI/passive_crown_FI - - if (ignite_passive_crown >= 1.0_r8) then - currentCohort%passive_crown_fire_flg = 1 !passive canopy fuels ignited for vertical spread - - !evaluate active crown fire conditions - - endif ! ignite crown fire - ! else no crown fire today - endif ! crown fire intensity - - ! 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 ! flames over bottom of canopy, and potentially over top of canopy if (currentCohort%hite > 0.0_r8 .and. currentCohort%SH >= height_cbb) then - - !add active flag !.and. currentCohort%active_crown_fire_flg == 0) then + !add active flag !.and. currentCohort%active_crown_fire_flg == 0) then currentCohort%fraction_crown_burned = min(1.0_r8, ((currentCohort%SH - height_cbb)/crown_depth)) + !add active flag !.and. currentCohort%active_crown_fire_flg == 1) then + !currentCohort%fraction_crown_burned = 1.0_r8 !1 for active fire = true endif !SH frac crown burnt calculation ! Check for strange values. currentCohort%fraction_crown_burned = min(1.0_r8, max(0.0_r8,currentCohort%fraction_crown_burned)) diff --git a/main/EDPftvarcon.F90 b/main/EDPftvarcon.F90 index 0425338043..807747a252 100644 --- a/main/EDPftvarcon.F90 +++ b/main/EDPftvarcon.F90 @@ -49,8 +49,7 @@ module EDPftvarcon ! that is occupied by crown. For fire model. real(r8), allocatable :: bark_scaler(:) ! scaler from dbh to bark thickness. For fire model. real(r8), allocatable :: crown_kill(:) ! crown resistance to fire. (1 = none) For fire model. - real(r8), allocatable :: active_crown_fire(:) ! Is plant susceptible to active crown fire?(1=yes,0=no) - real(r8), allocatable :: foliar_moisture(:) ! foliar moisture content from dry fuel [%] + real(r8), allocatable :: active_crown_fire(:) ! Is plant susceptible to active crown fire? real(r8), allocatable :: initd(:) ! initial seedling density real(r8), allocatable :: seed_suppl(:) ! seeds that come from outside the gridbox. real(r8), allocatable :: BB_slope(:) ! ball berry slope parameter @@ -387,10 +386,6 @@ subroutine Register_PFT(this, fates_params) call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) - name = 'fates_fire_foliar_moisture' - call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & - dimension_names=dim_names, lower_bounds=dim_lower_bound) - name = 'fates_recruit_initd' call fates_params%RegisterParameter(name=name, dimension_shape=dimension_shape_1d, & dimension_names=dim_names, lower_bounds=dim_lower_bound) @@ -827,10 +822,6 @@ subroutine Receive_PFT(this, fates_params) call fates_params%RetreiveParameterAllocate(name=name, & data=this%active_crown_fire) - name = 'fates_fire_foliar_moisture' - call fates_params%RetreiveParameterAllocate(name=name, & - data=this%foliar_moisture) - name = 'fates_recruit_initd' call fates_params%RetreiveParameterAllocate(name=name, & data=this%initd) @@ -1740,8 +1731,7 @@ subroutine FatesReportPFTParams(is_master) write(fates_log(),fmt0) 'crown = ',EDPftvarcon_inst%crown write(fates_log(),fmt0) 'bark_scaler = ',EDPftvarcon_inst%bark_scaler write(fates_log(),fmt0) 'crown_kill = ',EDPftvarcon_inst%crown_kill - write(fates_log(),fmt0) 'crown_fire = ',EDPftvarcon_inst%crown_fire - write(fates_log(),fmt0) 'crown_ignite_energy = ',EDPftvarcon_inst%crown_ignite_energy + write(fates_log(),fmt0) 'active_crown_fire = ',EDPftvarcon_inst%active_crown_fire write(fates_log(),fmt0) 'initd = ',EDPftvarcon_inst%initd write(fates_log(),fmt0) 'seed_suppl = ',EDPftvarcon_inst%seed_suppl write(fates_log(),fmt0) 'BB_slope = ',EDPftvarcon_inst%BB_slope diff --git a/parameter_files/fates_params_default.cdl b/parameter_files/fates_params_default.cdl index c29c3c60dd..87e1659863 100644 --- a/parameter_files/fates_params_default.cdl +++ b/parameter_files/fates_params_default.cdl @@ -133,11 +133,8 @@ variables: fates_fire_crown_kill:units = "NA" ; fates_fire_crown_kill:long_name = "crown resistance to fire (1 = none), EQ 22 in Thonicke et al 2010" ; double fates_fire_active_crown_fire(fates_pft) ; - fates_fire_active_crown_fire:units = "logical flag" ; - fates_fire_active_crown_fire:long_name = "Binary flag for active crown fire (1=yes)"; - double fates_fire_foliar_moisture(fates_pft) ; - fates_fire_foliar_moisture:units = "%" ; - fates_fire_foliar_moisture:long_name = "spitfire parameter, foliar moisture content based on dry fuel (%)" ; + fates_fire_active_crown_fire:units = "fraction" ; + fates_fire_active_crown_fire:long_name = "resistance to active crown fire (1=none)"; double fates_fr_fcel(fates_pft) ; fates_fr_fcel:units = "fraction" ; fates_fr_fcel:long_name = "Fine root litter cellulose fraction" ; @@ -726,8 +723,6 @@ data: fates_fire_active_crown_fire = 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 ; - fates_fire_foliar_moisture = 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100 ; - fates_fr_fcel = 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5, 0.5 ; fates_fr_flab = 0.25, 0.25, 0.25, 0.25, 0.25, 0.25, 0.25, 0.25, 0.25, 0.25,