using variable fuel moisture in fmc_g

(cherry picked from commit 5a6c4b8)
http://github.com/openwfm/WRF-Fire-merge/commit/5a6c4b8c06

phys/module_fr_fire_driver.F

@@ -104,6 +104,7 @@ subroutine fire_driver_em ( grid , config_flags                    &
     fp%fgip => grid%fgip     ! a rate of spread formula coeff 
     fp%ischap => grid%ischap ! a rate of spread formula switch
     fp%iboros => grid%iboros ! Ib divided by ROS
+    fp%fmc_g => grid%fmc_g   ! fuel moisture, ground
 
     ! get ignition data 
     call fire_ignition_convert (config_flags,fire_max_ignitions,fire_ignition_longlat, &
@@ -1553,6 +1554,7 @@ subroutine set_flags(config_flags)
 fire_viscosity_band     = config_flags%fire_viscosity_band 
 fire_viscosity_ngp      = config_flags%fire_viscosity_ngp
 fire_slope_factor       = config_flags%fire_slope_factor 
+fire_fmc_read           = config_flags%fire_fmc_read
 
 
 

phys/module_fr_fire_model.F

@@ -440,7 +440,7 @@ subroutine fire_model (                    &
                    ifts,ifte,jfts,jfte, &
                    fuel_frac_burnt,'model: fuel_frac burned')
 
-    call heat_fluxes(dt,                          &
+    call heat_fluxes(dt,fp,                       &
         ifms,ifme,jfms,jfme,                      &
         ifts,ifte,jfts,jfte,                      &
         ifts,ifte,jfts,jfte,                      &  ! fuel_frac_burned is tile dimensioned

phys/module_fr_fire_phys.F

@@ -20,6 +20,8 @@ module module_fr_fire_phys
 real,pointer,dimension(:,:):: fgip                 ! init mass of surface fuel (kg/m^2)
 real,pointer,dimension(:,:):: ischap               ! if fuel is chaparral and want rate of spread treated differently
 real,pointer,dimension(:,:):: iboros               ! fire instensity over rate of spread -> avialable fuel x heat of combustion
+real,pointer,dimension(:,:):: fuel_time            ! time to burn to 1/e (s) 
+real,pointer,dimension(:,:):: fmc_g                ! fuel moisture contents, ground (1) 
 end type fire_params
 
 ! use as
@@ -82,7 +84,7 @@ module module_fr_fire_phys
 ! scalar fuel coefficients
    REAL, SAVE:: cmbcnst,hfgl,fuelmc_g,fuelmc_c
 ! computed values
-   REAL, SAVE:: bmst,fuelheat
+   REAL, SAVE:: fuelheat
 
 ! defaults, may be changed in init_fuel_cats
    DATA cmbcnst  / 17.433e+06/             ! J/kg dry fuel
@@ -197,7 +199,6 @@ subroutine init_fuel_cats
 
 ! compute derived scalars
 
-bmst     = fuelmc_g/(1+fuelmc_g)
 fuelheat = cmbcnst * 4.30e-04     ! convert J/kg to BTU/lb
 
 ! compute derived fuel category coefficients 
@@ -223,8 +224,6 @@ subroutine init_fuel_cats
 call message(msg)
 write(msg,8)'fuelmc_c   ',fuelmc_c
 call message(msg)
-write(msg,8)'bmst       ',bmst
-call message(msg)
 write(msg,8)'fuelheat   ',fuelheat
 call message(msg)
 write(msg,7)'nfuelcats  ',nfuelcats
@@ -296,7 +295,7 @@ subroutine write_fuels_m(nsteps,maxwind,maxslope)
 !real,pointer,dimension(:,:):: fgip                 ! init mass of surface fuel (kg/m^2)
 !real,pointer,dimension(:,:):: ischap               ! 1 if chapparal
 !end type fire_params
-real, dimension(1:2,1:nsteps), target::vx,vy,zsf,dzdxf,dzdyf,bbb,betafl,phiwc,r_0,fgip,ischap
+real, dimension(1:2,1:nsteps), target::vx,vy,zsf,dzdxf,dzdyf,bbb,betafl,phiwc,r_0,fgip,ischap,fmc_g
 real, dimension(1:2,1:nsteps), target::iboros
 real, dimension(1:2,1:nsteps)::nfuel_cat,fuel_time,ros
 real::ros_base,ros_wind,ros_slope,propx,propy,r
@@ -312,6 +311,7 @@ subroutine write_fuels_m(nsteps,maxwind,maxslope)
 fp%fgip=>fgip
 fp%ischap=>ischap
 fp%iboros=>iboros 
+fp%fmc_g=>fmc_g
 
 iounit = open_text_file('fuels.m','write')
 
@@ -433,6 +433,7 @@ subroutine set_fire_params( &
         qig, epsilon, rhob, wn, betaop, e, c, &
         xifr, etas, etam, a, gammax, gamma, ratio, ir, &
         fuelloadm,fdxinv,fdyinv
+real:: bmst
 integer:: i,j,k
 integer::nerr
 character(len=128)::msg
@@ -556,16 +557,20 @@ subroutine set_fire_params( &
 
         fp%ischap(i,j)=ichap(k)
         fp%fgip(i,j)=fgi(k)
+        if(fire_fmc_read.eq.1)then
+           fp%fmc_g(i,j)=fuelmc_g
+        endif
 
         !     ...Settings of fire spread parameters from Rothermel follows. These
         !        don't need to be recalculated later.
         
+        bmst     = fp%fmc_g(i,j) / (1.+fp%fmc_g(i,j))
         fuelloadm= (1.-bmst) * fgi(k)  !  fuelload without moisture
         fuelload = fuelloadm * (.3048)**2 * 2.205    ! to lb/ft^2
         fueldepth = fueldepthm(k)/0.3048               ! to ft
         fp%betafl(i,j) = fuelload/(fueldepth * fueldens(k))! packing ratio
         betaop = 3.348 * savr(k)**(-0.8189)     ! optimum packing ratio
-        qig = 250. + 1116.*fuelmc_g            ! heat of preignition, btu/lb
+        qig = 250. + 1116.*fp%fmc_g(i,j)            ! heat of preignition, btu/lb
         epsilon = exp(-138./savr(k) )    ! effective heating number
         rhob = fuelload/fueldepth    ! ovendry bulk density, lb/ft^3
 
@@ -581,7 +586,7 @@ subroutine set_fire_params( &
         gamma = gammax *(ratio**a) *exp(a*(1.-ratio)) !optimum rxn vel, 1/min
 
         wn = fuelload/(1 + st(k))       ! net fuel loading, lb/ft^2
-        rtemp1 = fuelmc_g/fuelmce(k)
+        rtemp1 = fp%fmc_g(i,j)/fuelmce(k)
         etam = 1.-2.59*rtemp1 +5.11*rtemp1**2 -3.52*rtemp1**3  !moist damp coef
         etas = 0.174* se(k)**(-0.19)                ! mineral damping coef
         ir = gamma * wn * fuelheat * etam * etas  !rxn intensity,btu/ft^2 min
@@ -611,7 +616,7 @@ end subroutine set_fire_params
 !*******************
 !
 
-subroutine heat_fluxes(dt,                        &
+subroutine heat_fluxes(dt,fp,                     &
         ifms,ifme,jfms,jfme,                      &  ! memory dims
         ifts,ifte,jfts,jfte,                      &  ! tile dims
         iffs,iffe,jffs,jffe,                      &  ! fuel_frac_burnt dims
@@ -623,6 +628,7 @@ subroutine heat_fluxes(dt,                        &
 ! compute the heat fluxes on the fire grid cells
 
 !*** arguments
+type(fire_params), intent(in)::fp
 real, intent(in)::dt          ! dt  the fire time step (the fire model advances time by this)
 integer, intent(in)::ifts,ifte,jfts,jfte,ifms,ifme,jfms,jfme,iffs,iffe,jffs,jffe   ! dimensions                   
 real, intent(in),dimension(ifms:ifme,jfms:jfme):: fgip
@@ -631,14 +637,16 @@ subroutine heat_fluxes(dt,                        &
 
 !*** local
 integer::i,j
-real:: dmass
+real:: dmass,bmst
+logical::latent
 
 !*** executable        
 do j=jfts,jfte
     do i=ifts,ifte
          dmass =                     &     ! surface fuel dry mass burnt this call (kg/m^2)
              fgip(i,j)               &     ! init mass from fuel model no (kg/m^2) = fgi(nfuel_cat(i,j)
              * fuel_frac_burnt(i,j)        ! fraction burned this call    (1)
+         bmst     = fp%fmc_g(i,j)/(1.+fp%fmc_g(i,j))        ! bmst varies by (i,j)
          grnhft(i,j) = (dmass/dt)*(1.-bmst)*cmbcnst         ! surface fire sensible heat flux W/m^2
          grnqft(i,j) = (bmst+(1.-bmst)*.56)*(dmass/dt)*xlv  ! surface fire latent heat flux W/m
          ! xlv is defined in module_model_constants.. Assume 56% of cellulose molecule mass is water.

phys/module_fr_fire_util.F

@@ -28,6 +28,7 @@ module module_fr_fire_util
  fire_test_steps=0,     &  ! 0=no fire, 1=normal, >1 = do specified number of steps and terminate (testing only)
  fire_topo_from_atm=1,  &  ! 0 = expect ZSF set correctly on entry, 1 = populate by interploating from atmosphere
  fire_advection=0,      &  ! 0 = fire spread from normal wind/slope (CAWFE), 1 = full speed projected
+ fire_fmc_read=1,       &  ! fuel moisture: 0 from wrfinput, 1 from namelist.fire, 2 read from file in ideal
  fire_upwinding_reinit=4, &  ! spatial discretization for reinitialization PDE: 1=WENO3, 2=WENO5, 3=hybrid WENO3-ENO1, 4=hybrid WENO5-ENO1 
  fire_lsm_reinit_iter=1,  &  ! number of iterations for the reinitialization PDE 
  fire_lsm_band_ngp=4,     &  ! number of grid points around lfn=0 that the higher-order scheme WENO5/WENO3 is used (ENO1 elsewhere), for ire_upwinding=8,9 and fire_upwinding_reinit=3,4 and