Skip to content

Commit

Permalink
Merge pull request #249 from jkshuman/jkshuman-AreaBurnt
Browse files Browse the repository at this point in the history
Update AreaBurnt, add DEBUG, fix indexing
  • Loading branch information
rgknox authored Jul 25, 2017
2 parents 03186af + 4cf5aa9 commit 692a8b2
Show file tree
Hide file tree
Showing 4 changed files with 29 additions and 26 deletions.
2 changes: 1 addition & 1 deletion biogeochem/EDCanopyStructureMod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -1233,7 +1233,7 @@ subroutine leaf_area_profile( currentSite , snow_depth_si, frac_sno_eff_si)
write(fates_log(), *) 'ED: canopy-area-profile wrong', &
sum(currentPatch%canopy_area_profile(L,1:numpft_ed,1)), &
currentPatch%patchno, L
write(fates_log(), *) 'ED: areas',currentPatch%canopy_area_profile(L,1:2,1),currentPatch%patchno
write(fates_log(), *) 'ED: areas',currentPatch%canopy_area_profile(L,1:numpft_ed,1),currentPatch%patchno

currentCohort => currentPatch%shortest

Expand Down
2 changes: 1 addition & 1 deletion biogeochem/EDPhysiologyMod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -960,7 +960,7 @@ subroutine Growth_Derivatives( currentSite, currentCohort, bc_in)
!FIX(RF,032414) - to fix high bl's. needed to prevent numerical errors without the ODEINT.
if (currentCohort%balive > target_balive*1.1_r8)then
va = 0.0_r8; vs = 1._r8
write(fates_log(),*) 'using high bl cap',target_balive,currentCohort%balive
if (DEBUG) write(fates_log(),*) 'using high bl cap',target_balive,currentCohort%balive
endif

else
Expand Down
31 changes: 16 additions & 15 deletions biogeophys/EDSurfaceAlbedoMod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -830,11 +830,12 @@ subroutine ED_Norman_Radiation (nsites, sites, bc_in, bc_out )
endif
enddo
enddo

if (lai_change(1,2,1).gt.0.0.and.lai_change(1,2,2).gt.0.0)then
! write(fates_log(),*) 'lai_change(1,2,12)',lai_change(1,2,1:4)
! write(fates_log(),*) 'lai_change(1,2,12)',lai_change(1,2,1:4)
endif
if (lai_change(1,2,2).gt.0.0.and.lai_change(1,2,3).gt.0.0)then
! write(fates_log(),*) ' lai_change (1,2,23)',lai_change(1,2,1:4)
! write(fates_log(),*) ' lai_change (1,2,23)',lai_change(1,2,1:4)
endif
if (lai_change(1,1,3).gt.0.0.and.lai_change(1,1,2).gt.0.0)then
! NO-OP
Expand All @@ -846,7 +847,7 @@ subroutine ED_Norman_Radiation (nsites, sites, bc_in, bc_out )
endif
if (lai_change(1,1,4).gt.0.0.and.lai_change(1,1,5).gt.0.0)then
! NO-OP
! write(fates_log(),*) 'first layer of lai_change 4 5',lai_change(1,1,1:5)
! write(fates_log(),*) 'first layer of lai_change 4 5',lai_change(1,1,1:5)
endif

if (radtype == 1)then
Expand All @@ -865,10 +866,10 @@ subroutine ED_Norman_Radiation (nsites, sites, bc_in, bc_out )
write(fates_log(),*) 'Large Dir Radn consvn error',error ,ifp,ib
write(fates_log(),*) 'diags', bc_out(s)%albd_parb(ifp,ib), bc_out(s)%ftdd_parb(ifp,ib), &
bc_out(s)%ftid_parb(ifp,ib), bc_out(s)%fabd_parb(ifp,ib)
write(fates_log(),*) 'lai_change',lai_change(currentpatch%ncl_p,1:2,1:4)
write(fates_log(),*) 'elai',currentpatch%elai_profile(currentpatch%ncl_p,1:2,1:4)
write(fates_log(),*) 'esai',currentpatch%esai_profile(currentpatch%ncl_p,1:2,1:4)
write(fates_log(),*) 'ftweight',ftweight(1,1:2,1:4)
write(fates_log(),*) 'lai_change',lai_change(currentpatch%ncl_p,1:numpft_ed,1:4)
write(fates_log(),*) 'elai',currentpatch%elai_profile(currentpatch%ncl_p,1:numpft_ed,1:4)
write(fates_log(),*) 'esai',currentpatch%esai_profile(currentpatch%ncl_p,1:numpft_ed,1:4)
write(fates_log(),*) 'ftweight',ftweight(1,1:numpft_ed,1:4)
write(fates_log(),*) 'cp',currentPatch%area, currentPatch%patchno
write(fates_log(),*) 'bc_in(s)%albgr_dir_rb(ib)',bc_in(s)%albgr_dir_rb(ib)

Expand All @@ -884,16 +885,16 @@ subroutine ED_Norman_Radiation (nsites, sites, bc_in, bc_out )
write(fates_log(),*) '>5% Dif Radn consvn error',error ,ifp,ib
write(fates_log(),*) 'diags', bc_out(s)%albi_parb(ifp,ib), bc_out(s)%ftii_parb(ifp,ib), &
bc_out(s)%fabi_parb(ifp,ib)
write(fates_log(),*) 'lai_change',lai_change(currentpatch%ncl_p,1:2,1:4)
write(fates_log(),*) 'elai',currentpatch%elai_profile(currentpatch%ncl_p,1:2,1:4)
write(fates_log(),*) 'esai',currentpatch%esai_profile(currentpatch%ncl_p,1:2,1:4)
write(fates_log(),*) 'ftweight',ftweight(currentpatch%ncl_p,1:2,1:4)
write(fates_log(),*) 'lai_change',lai_change(currentpatch%ncl_p,1:numpft_ed,1:4)
write(fates_log(),*) 'elai',currentpatch%elai_profile(currentpatch%ncl_p,1:numpft_ed,1:4)
write(fates_log(),*) 'esai',currentpatch%esai_profile(currentpatch%ncl_p,1:numpft_ed,1:4)
write(fates_log(),*) 'ftweight',ftweight(currentpatch%ncl_p,1:numpft_ed,1:4)
write(fates_log(),*) 'cp',currentPatch%area, currentPatch%patchno
write(fates_log(),*) 'bc_in(s)%albgr_dif_rb(ib)',bc_in(s)%albgr_dif_rb(ib)
write(fates_log(),*) 'rhol',rhol(1:2,:)
write(fates_log(),*) 'ftw',sum(ftweight(1,:,1)),ftweight(1,1:2,1)
write(fates_log(),*) 'present',currentPatch%present(1,1:2)
write(fates_log(),*) 'CAP',currentPatch%canopy_area_profile(1,1:2,1)
write(fates_log(),*) 'rhol',rhol(1:numpft_ed,:)
write(fates_log(),*) 'ftw',sum(ftweight(1,:,1)),ftweight(1,1:numpft_ed,1)
write(fates_log(),*) 'present',currentPatch%present(1,1:numpft_ed)
write(fates_log(),*) 'CAP',currentPatch%canopy_area_profile(1,1:numpft_ed,1)

bc_out(s)%albi_parb(ifp,ib) = bc_out(s)%albi_parb(ifp,ib) + error
end if
Expand Down
20 changes: 11 additions & 9 deletions fire/SFMainMod.F90
Original file line number Diff line number Diff line change
Expand Up @@ -344,7 +344,7 @@ subroutine wind_effect ( currentSite, bc_in)
currentCohort => currentPatch%tallest

do while(associated(currentCohort))
write(fates_log(),*) 'SF currentCohort%c_area ',currentCohort%c_area
if (DEBUG) write(fates_log(),*) 'SF currentCohort%c_area ',currentCohort%c_area
if(EDPftvarcon_inst%woody(currentCohort%pft) == 1)then
currentPatch%total_tree_area = currentPatch%total_tree_area + currentCohort%c_area
else
Expand Down Expand Up @@ -477,23 +477,23 @@ subroutine rate_of_spread ( currentSite )
! Equation A5 in Thonicke et al. 2010
! phi_wind (unitless)
! convert wind_elev_fire from m/min to ft/min for Rothermel ROS eqn
! wind max per Lasslop et al 2014 to lenearly reduce ROS for high wind speeds
! wind max per Lasslop et al 2014 to linearly reduce ROS for high wind speeds
!OLD! phi_wind = c * ((3.281_r8*currentPatch%effect_wspeed)**b)*(beta_ratio**(-e))
if (currentPatch%effect_wspeed .le. wind_max) then
wind_elev_fire = currentPatch%effect_wspeed
phi_wind = c * ((3.281_r8*wind_elev_fire)**b)*(beta_ratio**(-e))
if (debug_windspeed) write(fates_log(),*) 'SF wind LESS max ', currentPatch%effect_wspeed
if (debug_windspeed) write(fates_log(),*) 'month and day', hlm_current_month, hlm_current_day
else
! max conditional 225 ft/min from Lasslop 2014 converted to 68.577 m/min
!max condition 225 ft/min (FIREMIP Rabin table A10 JSBACH-Spitfire) convert to 68.577 m/min
wind_elev_fire = max(0.0_r8,(68.577-0.5*currentPatch%effect_wspeed))
phi_wind = c * ((3.281_r8*wind_elev_fire)**b)*(beta_ratio**(-e))
if (debug_windspeed) write(fates_log(),*) 'SF wind GREATER max ', currentPatch%effect_wspeed
if (debug_windspeed) write(fates_log(),*) 'month and day', hlm_current_month, hlm_current_day
endif

! ---propagating flux----
! Equation A2 in Thonicke et al.2010
! Equation A2 in Thonicke et al.2010 and Eq. 42 Rothermal 1972
! xi (unitless)
xi = (exp((0.792_r8 + 3.7597_r8 * (currentPatch%fuel_sav**0.5_r8)) * (beta+0.1_r8))) / &
(192_r8+7.9095_r8 * currentPatch%fuel_sav)
Expand Down Expand Up @@ -752,20 +752,22 @@ subroutine area_burnt ( currentSite )
! THIS SHOULD HAVE THE COLUMN AND LU AREA WEIGHT ALSO, NO?

gridarea = km2_to_m2 ! 1M m2 in a km2
currentPatch%NF = ED_val_nignitions * currentPatch%area/area /365 * &
currentSite%FDI
!NF = number of lighting strikes per day per km2
currentPatch%NF = ED_val_nignitions * currentPatch%area/area /365

! If there are 15 lightening strickes per year, per km2. (approx from NASA product)
! then there are 15/365 s/km2 each day.

! Equation 1 in Thonicke et al. 2010
! To Do: Connect here with the Li & Levis GDP fire suppression algorithm.
! Equation 16 in arora and boer model.
! Equation 16 in arora and boer model JGR 2005
!currentPatch%AB = currentPatch%AB *3.0_r8

!size of fire = equation 14 Arora and Boer JGR 2005
size_of_fire = ((3.1416_r8/(4.0_r8*lb))*((df+db)**2.0_r8))

!AB is daily area burnt = size of fires in m2 * number of ignitions
currentPatch%AB = size_of_fire * currentPatch%NF
!AB = daily area burnt = size fires in m2 * num ignitions * prob ignition starts fire
currentPatch%AB = size_of_fire * currentPatch%NF * currentSite%FDI

patch_area_in_m2 = gridarea*currentPatch%area/area

Expand Down

0 comments on commit 692a8b2

Please sign in to comment.