Merge pull request #13140 from ericvmueller/master

add VEG_LSET_ROS_FIXED option for static level set spread rate
firemodels · Jul 8, 2024 · be369fb · be369fb
2 parents d8d6479 + f5e690d
commit be369fb
Show file tree

Hide file tree

Showing 4 changed files with 55 additions and 40 deletions.
diff --git a/Manuals/FDS_User_Guide/FDS_User_Guide.tex b/Manuals/FDS_User_Guide/FDS_User_Guide.tex
@@ -7456,6 +7456,7 @@ \section{Level Set Model for Wildland Fire Spread}
 \ee
 Here, $m_{\rm f}''$ is the dry fuel loading in units of kg/m$^2$ ({\ct VEG\_LSET\_SURF\_LOAD}, default 0.3~kg/m$^2$); $\nu_{\rm char}$ is the char fraction ({\ct VEG\_LSET\_CHAR\_FRACTION}, default 0.2); and $\delta t$ is the duration of the fire at a given location. The expression for the burn duration is given by Albini~\cite{Albini:1976}, but you can enter your own burn duration using {\ct VEG\_LSET\_FIREBASE\_TIME} (s).
 
+It is possible to set a custom fuel for which the spread rate does not vary with slope or wind speed. This is done by specifying {\ct VEG\_LSET\_ROS\_FIXED=T} on the {\ct SURF} line. The spread rate is then fixed at the no-wind, no-slope value. The main purpose of this feature is for reproduction of manual ignition patterns (e.g. ignition with a drip torch at a fixed walking speed).
 
 
 \begin{sidewaystable}[p]
@@ -13162,6 +13163,7 @@ \section{\texorpdfstring{{\tt SURF}}{SURF} (Surface Properties)}
 {\ct VEG\_LSET\_ROS\_BACK}                                & Real            & Section~\ref{info:level_set}              & m/s                 & 0.                      \\ \hline
 {\ct VEG\_LSET\_ROS\_FLANK}                               & Real            & Section~\ref{info:level_set}              & m/s                 & 0.                      \\ \hline
 {\ct VEG\_LSET\_ROS\_HEAD}                                & Real            & Section~\ref{info:level_set}              & m/s                 & 0.                      \\ \hline
+{\ct VEG\_LSET\_ROS\_FIXED}                               & Logical         & Section~\ref{info:level_set}              &                     & {\ct F}                        \\ \hline
 {\ct VEG\_LSET\_SIGMA}                                    & Real            & Section~\ref{info:level_set}              & 1/m                 & 0.                      \\ \hline
 {\ct VEG\_LSET\_SURF\_LOAD}                               & Real            & Section~\ref{info:level_set}              & kg/m$^2$            & 0.3                     \\ \hline
 {\ct VEG\_LSET\_TAN2}                                     & Real            & Section~\ref{info:level_set}              &                     &                         \\ \hline

diff --git a/Source/read.f90 b/Source/read.f90
@@ -7512,7 +7512,7 @@ SUBROUTINE READ_SURF(QUICK_READ)
             VEG_LSET_M1,VEG_LSET_M10,VEG_LSET_M100,VEG_LSET_MLW,VEG_LSET_MLH,VEG_LSET_SURF_LOAD,VEG_LSET_FIREBASE_TIME,&
             VEG_LSET_CHAR_FRACTION,VEL_PART,INIT_PER_AREA,TIME_STEP_FACTOR
 LOGICAL :: DEFAULT,VEG_LSET_SPREAD,VEG_LSET_TAN2,TGA_ANALYSIS,COMPUTE_EMISSIVITY,&
-           COMPUTE_EMISSIVITY_BACK,VARIABLE_THICKNESS,HT3D,THERM_THICK
+           COMPUTE_EMISSIVITY_BACK,VARIABLE_THICKNESS,HT3D,THERM_THICK,VEG_LSET_ROS_FIXED
 LOGICAL, ALLOCATABLE, DIMENSION(:) :: DUPLICATE
 ! Ember generating variables
 REAL(EB) :: EMBER_GENERATION_HEIGHT(2),EMBER_IGNITION_POWER_MEAN,EMBER_IGNITION_POWER_SIGMA,EMBER_TRACKING_RATIO,EMBER_YIELD
@@ -7545,7 +7545,7 @@ SUBROUTINE READ_SURF(QUICK_READ)
                 TMP_BACK,TMP_FRONT,TMP_FRONT_INITIAL,TMP_GAS_BACK,TMP_GAS_FRONT,TMP_INNER,TRANSPARENCY,&
                 VEG_LSET_BETA,VEG_LSET_CHAR_FRACTION,VEG_LSET_FIREBASE_TIME,VEG_LSET_FUEL_INDEX,VEG_LSET_HT,VEG_LSET_IGNITE_TIME,&
                 VEG_LSET_M1,VEG_LSET_M10,VEG_LSET_M100,VEG_LSET_MLW,VEG_LSET_MLH,VEG_LSET_QCON,&
-                VEG_LSET_ROS_00,VEG_LSET_ROS_BACK,VEG_LSET_ROS_FLANK,VEG_LSET_ROS_HEAD,VEG_LSET_SIGMA,&
+                VEG_LSET_ROS_00,VEG_LSET_ROS_BACK,VEG_LSET_ROS_FLANK,VEG_LSET_ROS_HEAD,VEG_LSET_ROS_FIXED,VEG_LSET_SIGMA,&
                 VEG_LSET_SURF_LOAD,VEG_LSET_TAN2,VEG_LSET_WIND_EXP,&
                 VEL,VEL_BULK,VEL_GRAD,VEL_PART,VEL_T,VOLUME_FLOW,WIDTH,XYZ,Z0,Z_0
 
@@ -7823,6 +7823,7 @@ SUBROUTINE READ_SURF(QUICK_READ)
    SF%VEG_LSET_ROS_HEAD     = VEG_LSET_ROS_HEAD     ! head fire rate of spread (m/s), McArthur model
    SF%VEG_LSET_ROS_FLANK    = VEG_LSET_ROS_FLANK    ! flank fire rate of spread, McArthur model
    SF%VEG_LSET_ROS_BACK     = VEG_LSET_ROS_BACK     ! back fire rate of spread, McArthur model
+   SF%VEG_LSET_ROS_FIXED    = VEG_LSET_ROS_FIXED    ! ROS is a fixed value and does not change with slope or wind
    SF%VEG_LSET_WIND_EXP     = VEG_LSET_WIND_EXP     ! exponent on wind cosine in ROS formula
    SF%VEG_LSET_SIGMA        = VEG_LSET_SIGMA * 0.01 ! SAV for Farsite emulation in LSET converted to 1/cm
    SF%VEG_LSET_HT           = VEG_LSET_HT
@@ -8740,6 +8741,7 @@ SUBROUTINE SET_SURF_DEFAULTS
 VEG_LSET_ROS_HEAD       = 0.0_EB
 VEG_LSET_ROS_FLANK      = 0.0_EB
 VEG_LSET_ROS_BACK       = 0.0_EB
+VEG_LSET_ROS_FIXED      = .FALSE.
 VEG_LSET_WIND_EXP       = 1.0_EB
 VEG_LSET_TAN2           = .FALSE.
 VEG_LSET_HT             = 1.0_EB

diff --git a/Source/type.f90 b/Source/type.f90
@@ -628,7 +628,7 @@ MODULE TYPES
    REAL(EB), ALLOCATABLE, DIMENSION(:) :: R50
    REAL(EB) :: OXR                  !< Required oxygen for complete combustion (gm/gm-species)
    REAL(EB) :: OXA                  !< Available oxygen for combustion (gm/gm-species)
-   
+
 
 END TYPE SPECIES_MIXTURE_TYPE
 
@@ -965,7 +965,7 @@ MODULE TYPES
 
    ! Level Set Firespread
 
-   LOGICAL :: VEG_LSET_SPREAD,VEG_LSET_TAN2
+   LOGICAL :: VEG_LSET_SPREAD,VEG_LSET_TAN2,VEG_LSET_ROS_FIXED
    REAL(EB) :: VEG_LSET_IGNITE_T,VEG_LSET_ROS_HEAD,VEG_LSET_ROS_00,VEG_LSET_QCON,VEG_LSET_ROS_FLANK,VEG_LSET_ROS_BACK, &
                VEG_LSET_WIND_EXP,VEG_LSET_SIGMA,VEG_LSET_HT,VEG_LSET_BETA,&
                VEG_LSET_M1,VEG_LSET_M10,VEG_LSET_M100,VEG_LSET_MLW,VEG_LSET_MLH,VEG_LSET_SURF_LOAD,VEG_LSET_FIREBASE_TIME, &

diff --git a/Source/vege.f90 b/Source/vege.f90
@@ -769,61 +769,72 @@ SUBROUTINE LEVEL_SET_SPREAD_RATE
 
       IF (LEVEL_SET_ELLIPSE) THEN
 
-         ! Effective wind direction (theta) is clockwise from y-axis (Richards 1990)
-         COS_THETA = COS(THETA_ELPS(I,J)) !V_LS(I,J) / MAG_U
-         SIN_THETA = SIN(THETA_ELPS(I,J)) !U_LS(I,J) / MAG_U
+         ! ROS does not change with wind or slope
+         IF (SURFACE(LS_SURF_INDEX(I,J))%VEG_LSET_ROS_FIXED) THEN
 
-         ROS_TMP = ROS_HEAD(I,J)
+            MAG_SR=SURFACE(LS_SURF_INDEX(I,J))%VEG_LSET_ROS_00
+            SR_X_LS(I,J) = MAG_SR*NORMAL_FIRELINE(1) !spread rate components
+            SR_Y_LS(I,J) = MAG_SR*NORMAL_FIRELINE(2)
 
-         ! Magnitude of wind speed at midflame height must be in units of m/s here
+         ELSE
 
-         UMF_DUM = UMF(I,J)/60.0_EB
+            ! Effective wind direction (theta) is clockwise from y-axis (Richards 1990)
+            COS_THETA = COS(THETA_ELPS(I,J)) !V_LS(I,J) / MAG_U
+            SIN_THETA = SIN(THETA_ELPS(I,J)) !U_LS(I,J) / MAG_U
 
-         ! Length to breadth ratio of ellipse based on effective UMF (Bova et al., Eq. A6)
+            ROS_TMP = ROS_HEAD(I,J)
 
-         LB = 0.936_EB * EXP(0.2566_EB * UMF_DUM) + 0.461_EB * EXP(-0.1548_EB * UMF_DUM) - 0.397_EB
+            ! Magnitude of wind speed at midflame height must be in units of m/s here
 
-         ! Constraint LB max = 8 from Finney 2004
+            UMF_DUM = UMF(I,J)/60.0_EB
 
-         LB = MAX(1.0_EB,MIN(LB,8.0_EB))  ! (Bova et al., Eq. A7)
+            ! Length to breadth ratio of ellipse based on effective UMF (Bova et al., Eq. A6)
 
-         ! Head to back ratio based on LB
+            LB = 0.936_EB * EXP(0.2566_EB * UMF_DUM) + 0.461_EB * EXP(-0.1548_EB * UMF_DUM) - 0.397_EB
 
-         LBD = SQRT(LB**2 - 1.0_EB)
-         HB = (LB + LBD) / (LB - LBD)
+            ! Constraint LB max = 8 from Finney 2004
 
-         ! A_ELPS and B_ELPS notation is consistent with Farsite and Richards
+            LB = MAX(1.0_EB,MIN(LB,8.0_EB))  ! (Bova et al., Eq. A7)
 
-         B_ELPS =  0.5_EB * (ROS_TMP + ROS_TMP/HB)
-         B_ELPS2 = B_ELPS**2
-         A_ELPS =  B_ELPS / LB
-         A_ELPS2=  A_ELPS**2
-         C_ELPS =  B_ELPS - (ROS_TMP/HB)
+            ! Head to back ratio based on LB
 
-         ! Denominator used in spread rate equation from Richards, Intnl. J. Num. Methods Eng. 1990
-         ! and in LS vs Farsite paper, Bova et al., Intnl. J. Wildland Fire, 25(2):229-241, 2015
+            LBD = SQRT(LB**2 - 1.0_EB)
+            HB = (LB + LBD) / (LB - LBD)
 
-         AROS  = XSF*COS_THETA - YSF*SIN_THETA
-         BROS  = XSF*SIN_THETA + YSF*COS_THETA
-         DENOM = A_ELPS2*BROS**2 + B_ELPS2*AROS**2
+            ! A_ELPS and B_ELPS notation is consistent with Farsite and Richards
 
-         IF (DENOM > 0._EB) THEN
-            DENOM = 1._EB / SQRT(DENOM)
-         ELSE
-            DENOM = 0._EB
-         ENDIF
+            B_ELPS =  0.5_EB * (ROS_TMP + ROS_TMP/HB)
+            B_ELPS2 = B_ELPS**2
+            A_ELPS =  B_ELPS / LB
+            A_ELPS2=  A_ELPS**2
+            C_ELPS =  B_ELPS - (ROS_TMP/HB)
 
-         ! This is with A_ELPS2 and B_ELPS2 notation consistent with Finney and Richards and in Bova et al. 2015 IJWF 2015
+            ! Denominator used in spread rate equation from Richards, Intnl. J. Num. Methods Eng. 1990
+            ! and in LS vs Farsite paper, Bova et al., Intnl. J. Wildland Fire, 25(2):229-241, 2015
 
-         SR_X_LS(I,J) = DENOM * ( A_ELPS2*COS_THETA*BROS - B_ELPS2*SIN_THETA*AROS) + C_ELPS*SIN_THETA  ! Bova et al., Eq. A8
-         SR_Y_LS(I,J) = DENOM * (-A_ELPS2*SIN_THETA*BROS - B_ELPS2*COS_THETA*AROS) + C_ELPS*COS_THETA  ! Bova et al., Eq. A9
+            AROS  = XSF*COS_THETA - YSF*SIN_THETA
+            BROS  = XSF*SIN_THETA + YSF*COS_THETA
+            DENOM = A_ELPS2*BROS**2 + B_ELPS2*AROS**2
+
+            IF (DENOM > 0._EB) THEN
+               DENOM = 1._EB / SQRT(DENOM)
+            ELSE
+               DENOM = 0._EB
+            ENDIF
 
-         ! Project spread rates from slope to horizontal plane
+            ! This is with A_ELPS2 and B_ELPS2 notation consistent with Finney and Richards and in Bova et al. 2015 IJWF 2015
 
-         IF (ABS(DZTDX(I,J)) > 0._EB) SR_X_LS(I,J) = SR_X_LS(I,J) * ABS(COS(ATAN(DZTDX(I,J))))
-         IF (ABS(DZTDY(I,J)) > 0._EB) SR_Y_LS(I,J) = SR_Y_LS(I,J) * ABS(COS(ATAN(DZTDY(I,J))))
+            SR_X_LS(I,J) = DENOM * ( A_ELPS2*COS_THETA*BROS - B_ELPS2*SIN_THETA*AROS) + C_ELPS*SIN_THETA  ! Bova et al., Eq. A8
+            SR_Y_LS(I,J) = DENOM * (-A_ELPS2*SIN_THETA*BROS - B_ELPS2*COS_THETA*AROS) + C_ELPS*COS_THETA  ! Bova et al., Eq. A9
 
-         MAG_SR = SQRT(SR_X_LS(I,J)**2 + SR_Y_LS(I,J)**2)
+            ! Project spread rates from slope to horizontal plane
+
+            IF (ABS(DZTDX(I,J)) > 0._EB) SR_X_LS(I,J) = SR_X_LS(I,J) * ABS(COS(ATAN(DZTDX(I,J))))
+            IF (ABS(DZTDY(I,J)) > 0._EB) SR_Y_LS(I,J) = SR_Y_LS(I,J) * ABS(COS(ATAN(DZTDY(I,J))))
+
+            MAG_SR = SQRT(SR_X_LS(I,J)**2 + SR_Y_LS(I,J)**2)
+
+         ENDIF
 
       ELSE ! McArthur Spread Model