Merge remote-tracking branch 'firemodels/master' into FireX

Author: cxp484

Manuals/FDS_User_Guide/FDS_User_Guide.tex

@@ -10444,7 +10444,7 @@ \subsection{Flame Height and Flame Tilt}
 The ``flame height,'' and to some extent the ``flame tilt,'' are common quantities used to assess the hazard from a large fire. FDS does not output these quantities directly, but there is a combination of devices and controls that yield these values. This is all best explained with an example. Suppose you are simulating a fire in a domain that is 10~m high with a grid resolution of 0.2~m. The fire is centered about the point $(x,y)=(0,0)$ at $z=0$. The following input lines compute the flame height and tilt:
 \begin{lstlisting}
 &DEVC ID='HRR', Z_ID='z', XBP=0.0,0.0,0.0,0.0,0.1,9.9, QUANTITY='HRRPUV',
-SPATIAL_STATISTIC='VOLUME INTEGRAL', DX=5, DY=5, DZ=0.1, POINTS=100,
+SPATIAL_STATISTIC='VOLUME INTEGRAL', DX=5, DY=5, DZ=0.1, POINTS=50,
 STATISTICS_START=10. /
 &CTRL ID='H', FUNCTION_TYPE='PERCENTILE', INPUT_ID='HRR', PERCENTILE=0.97 /
 &DEVC ID='x_max', XB=-5,5,-5,5,9.0,9.2, QUANTITY='TEMPERATURE',

Source/dump.f90

@@ -3629,6 +3629,7 @@ SUBROUTINE INITIALIZE_DIAGNOSTIC_FILE(DT)
       WRITE(LU_OUTPUT,'(A,ES10.3)') '     Packing Ratio               ', SF%VEG_LSET_BETA
       WRITE(LU_OUTPUT,'(A,ES10.3)') '     Surface Area/Volume (1/m)   ', SF%VEG_LSET_SIGMA*100.  ! Convert from 1/cm to 1/m
       WRITE(LU_OUTPUT,'(A,ES10.3)') '     Fuel Depth (m)              ', SF%VEG_LSET_HT
+      WRITE(LU_OUTPUT,'(A,ES10.3)') '     Fuel Load (kg/m2)           ', SF%VEG_LSET_SURF_LOAD
    ENDIF
 
 ENDDO SURFLOOP

Source/read.f90

@@ -8326,15 +8326,17 @@ SUBROUTINE READ_SURF(QUICK_READ)
    ENDIF
 
    ! Level set vegetation fire spread specific
+   ! VEG_LSET_FUEL_INDEX: -1=undefined, 0=custom Rothermel from layer props, 1-13=Albini fuel models
 
    VEG_LSET_SPREAD = .FALSE.
-   IF (VEG_LSET_IGNITE_TIME < 1.E6_EB .OR. VEG_LSET_FUEL_INDEX>0 .OR. VEG_LSET_ROS_00>0._EB) VEG_LSET_SPREAD = .TRUE.
+   ! SURF_LOAD without ROS_00/index implies custom fuel (index 0)
+   IF (VEG_LSET_FUEL_INDEX<0 .AND. VEG_LSET_ROS_00<=0._EB .AND. VEG_LSET_SURF_LOAD>0._EB) VEG_LSET_FUEL_INDEX = 0
+   IF (VEG_LSET_IGNITE_TIME < 1.E6_EB .OR. VEG_LSET_FUEL_INDEX>=0 .OR. VEG_LSET_ROS_00>0._EB) VEG_LSET_SPREAD = .TRUE.
    IF (VEG_LSET_SPREAD .AND. LEVEL_SET_MODE==0) THEN
       WRITE(MESSAGE,'(A,A,A)') 'ERROR(305): SURF ',TRIM(ID),' indicates a level set simulation, but LEVEL_SET_MODE not set on MISC.'
       CALL SHUTDOWN(MESSAGE) ; RETURN
    ENDIF
-   IF (VEG_LSET_FUEL_INDEX>0 .AND. LEVEL_SET_COUPLED_FIRE) HRRPUA = 1._EB  ! HRRPUA to be set properly later
-   IF (VEG_LSET_ROS_00    >0 .AND. LEVEL_SET_COUPLED_FIRE) HRRPUA = 1._EB
+   IF ((VEG_LSET_FUEL_INDEX>=0 .OR. VEG_LSET_ROS_00>0._EB) .AND. LEVEL_SET_COUPLED_FIRE) HRRPUA = 1._EB  ! set properly later
 
    SF%VEG_LSET_SPREAD       = VEG_LSET_SPREAD
    SF%VEG_LSET_ROS_00       = VEG_LSET_ROS_00       ! no-wind, no-slope RoS (m/s), Rothermel model
@@ -8376,8 +8378,8 @@ SUBROUTINE READ_SURF(QUICK_READ)
       END SELECT
    ENDIF
 
-   ! Set defaults for custom fuel model (for indexed fuels they are either user-specified or calculated later)
-   IF (SF%VEG_LSET_FUEL_INDEX==0._EB) THEN
+   ! Defaults for undefined fuel properties (set elsewhere for indexed fuels)
+   IF (SF%VEG_LSET_FUEL_INDEX<=0) THEN
       IF (SF%VEG_LSET_SIGMA<0._EB) SF%VEG_LSET_SIGMA = 50._EB
       IF (SF%VEG_LSET_HT<0._EB) SF%VEG_LSET_HT = 0.2_EB
       IF (SF%VEG_LSET_SURF_LOAD<0._EB) SF%VEG_LSET_SURF_LOAD = 1.0_EB
@@ -9306,7 +9308,7 @@ SUBROUTINE SET_SURF_DEFAULTS
 VEG_LSET_M100           = 0.05_EB
 VEG_LSET_MLW            = 0.70_EB
 VEG_LSET_MLH            = 0.70_EB
-VEG_LSET_FUEL_INDEX     = 0
+VEG_LSET_FUEL_INDEX     = -1  ! -1=undefined, 0=custom Rothermel, 1-13=Albini models
 VEG_LSET_SURF_LOAD      = -1.0_EB !kg/m^2
 VEG_LSET_FIREBASE_TIME  = -1.0_EB
 VEG_LSET_CHAR_FRACTION  = 0.20_EB

Source/vege.f90

@@ -46,17 +46,18 @@ SUBROUTINE INITIALIZE_LEVEL_SET_FIRESPREAD_1(NM)
 
 M => MESHES(NM)
 
-! Loop through all SURFace types and find level set cases that need a calculated RoS
+! Loop through all SURFace types and find level set cases that need a calculated RoS.
+! VEG_LSET_FUEL_INDEX: -1=undefined (user ROS_00), 0=custom Rothermel, 1-13=Albini models.
 
 DO SURF_INDEX=0,N_SURF
    SF => SURFACE(SURF_INDEX)
-   IF (SF%VEG_LSET_SPREAD .AND. SF%VEG_LSET_FUEL_INDEX>0) THEN
+   IF (.NOT.SF%VEG_LSET_SPREAD) CYCLE
+   IF (SF%VEG_LSET_FUEL_INDEX>=0) THEN
       SF%VEG_LSET_ROS_00 = ROS_NO_WIND_NO_SLOPE(SF%VEG_LSET_FUEL_INDEX,SURF_INDEX)
-   ENDIF
-   IF (SF%VEG_LSET_SPREAD .AND. SF%VEG_LSET_FUEL_INDEX==0) THEN
-     SF%BURN_DURATION = SF%VEG_LSET_FIREBASE_TIME
-     IF (LEVEL_SET_COUPLED_FIRE) SF%MASS_FLUX(REACTION(1)%FUEL_SMIX_INDEX) = &
-       (1._EB-SF%VEG_LSET_CHAR_FRACTION)*SF%VEG_LSET_SURF_LOAD/SF%VEG_LSET_FIREBASE_TIME
+   ELSE
+      SF%BURN_DURATION = SF%VEG_LSET_FIREBASE_TIME
+      IF (LEVEL_SET_COUPLED_FIRE) SF%MASS_FLUX(REACTION(1)%FUEL_SMIX_INDEX) = &
+        (1._EB-SF%VEG_LSET_CHAR_FRACTION)*SF%VEG_LSET_SURF_LOAD/SF%VEG_LSET_FIREBASE_TIME
    ENDIF
 ENDDO
 
@@ -1010,6 +1011,11 @@ REAL(EB) FUNCTION ROS_NO_WIND_NO_SLOPE(ROTHERMEL_FUEL_INDEX,SURF_INDEX)
 mlh = SF%VEG_LSET_MLH
 
 SELECT CASE(ROTHERMEL_FUEL_INDEX)
+   CASE(0)  ! Custom homogeneous fuel from SURF layer properties
+      w0d1=SF%VEG_LSET_SURF_LOAD ; w0d2=0._EB ; w0d3=0._EB ; w0lh=0._EB ; w0lw=0._EB
+      svd1=SF%VEG_LSET_SIGMA*100._EB  ! stored as 1/cm; Rothermel formulas use 1/m
+      svd2=358._EB ; svd3=98._EB ; svlh=4921._EB ; svlw=4921._EB
+      mx=0.15_EB ; depth=SF%VEG_LSET_HT ; rhop=512._EB ; heat=18607._EB ; st=0.0555_EB ; se=0.01_EB
    CASE(1)  ! 'Short Grass'
       w0d1=0.1659     ; w0d2=0.        ; w0d3=0.        ; w0lh=0.        ; w0lw=0.     ! dry mass per unit area (kg/m2)
       svd1=11483.     ; svd2=358.      ; svd3=98.       ; svlh=4921.     ; svlw=4921.  ! surface area to volume (1/m)
@@ -1076,6 +1082,8 @@ REAL(EB) FUNCTION ROS_NO_WIND_NO_SLOPE(ROTHERMEL_FUEL_INDEX,SURF_INDEX)
    w0lh = SF%VEG_LSET_SURF_LOAD/w0*w0lh
    w0lw = SF%VEG_LSET_SURF_LOAD/w0*w0lw
    w0 = (w0d1 + w0d2 + w0d3 + w0lh + w0lw)
+ELSE
+   SF%VEG_LSET_SURF_LOAD = w0
 ENDIF
 
 ! Auxiliary functions