FireX: Merge with firemodels/master

Manuals/Copy_Firebot_FiguresGH.sh

@@ -1,7 +1,12 @@
 #!/bin/bash
+LATEST=$1
+if [ "$LATEST" != "" ]; then
+  LATEST=_latest
+fi
+
 DOWNLOADFIGURES ()
 {
-  FILE=$1_figures
+  FILE=$1_figures$LATEST
   echo ***Downloading $FILE.tar.gz
   gh release download FDS_TEST -p $FILE.tar.gz -R github.com/firemodels/test_bundles --clobber
   if [ ! -e $FILE.tar.gz ]; then
@@ -31,25 +36,25 @@ FBVG=$BASEDIR/FIGS/FDS_VERG
 FBVAL=$BASEDIR/FIGS/FDS_VALG
 
 # Copy Tech Guide Figures
-if [ -e $BASEDIR/FDS_TG_figures.tar.gz ]; then
+if [ -e $BASEDIR/FDS_TG_figures$LATEST.tar.gz ]; then
   cd $FBTG
-  tar xf $BASEDIR/FDS_TG_figures.tar.gz
+  tar xf $BASEDIR/FDS_TG_figures$LATEST.tar.gz
   cp * $BASEDIR/FDS_Technical_Reference_Guide/SCRIPT_FIGURES/
   echo FDS Technical Guide figures copied to FDS_Technical_Reference_Guide/SCRIPT_FIGURES
 fi
 
 # Copy User's Guide Figures
-if [ -e $BASEDIR/FDS_UG_figures.tar.gz ]; then
+if [ -e $BASEDIR/FDS_UG_figures$LATEST.tar.gz ]; then
   cd $FBUG
-  tar xf $BASEDIR/FDS_UG_figures.tar.gz
+  tar xf $BASEDIR/FDS_UG_figures$LATEST.tar.gz
   cp * $BASEDIR/FDS_User_Guide/SCRIPT_FIGURES/
   echo FDS Users Guide figures copied to FDS_User_Reference_Guide/SCRIPT_FIGURES
 fi
 
 # Copy Verification Guide Figures
-if [ -e $BASEDIR/FDS_VERG_figures.tar.gz ]; then
+if [ -e $BASEDIR/FDS_VERG_figures$LATEST.tar.gz ]; then
   cd $FBVG
-  tar xf $BASEDIR/FDS_VERG_figures.tar.gz
+  tar xf $BASEDIR/FDS_VERG_figures$LATEST.tar.gz
   cp *.pdf $BASEDIR/FDS_Verification_Guide/SCRIPT_FIGURES/.
   cp *.png $BASEDIR/FDS_Verification_Guide/SCRIPT_FIGURES/.
   cp *.tex $BASEDIR/FDS_Verification_Guide/SCRIPT_FIGURES/.
@@ -58,9 +63,9 @@ if [ -e $BASEDIR/FDS_VERG_figures.tar.gz ]; then
 fi
 
 # Copy Validation Guide Figures
-if [ -e $BASEDIR/FDS_VALG_figures.tar.gz ]; then
+if [ -e $BASEDIR/FDS_VALG_figures$LATEST.tar.gz ]; then
   cd $FBVAL
-  tar xf $BASEDIR/FDS_VALG_figures.tar.gz
+  tar xf $BASEDIR/FDS_VALG_figures$LATEST.tar.gz
   cp -R * $BASEDIR/FDS_Validation_Guide/SCRIPT_FIGURES
   echo FDS Validation Guide Figures copied to FDS_Validation_Guide_Reference_Guide/SCRIPT_FIGURES
 fi

Manuals/FDS_User_Guide/FDS_User_Guide.tex

@@ -4621,7 +4621,7 @@ \subsubsection{Enthalpy}
 \be
    c_{p,\alpha} = \frac{\gamma}{\gamma-1} \frac{\R}{W_\alpha}
 \ee
-The ratio of specific heats, \ct{GAMMA}, is 1.4 by default and can be changed on the \ct{MISC} line. If you want all the gas specific heats to follow this relation, set \ct{CONSTANT_SPECIFIC_HEAT_RATIO=T} on the \ct{MISC} line. For high molecular weight species, use of the default gamma will result in very low values of the specific heat which can cause issues with the default extinction model. In this case it is recommended that the specific heat be defined. If the FDS is determining the specific heat using this relation, then it is recommended that you check the \ct{CHID.out} and verify that reasonable specific heat values have been created.
+The ratio of specific heats, \ct{GAMMA}, is 1.4 by default and can be changed on the \ct{MISC} line. If you want all the gas specific heats to follow this relation, set \ct{CONSTANT_SPECIFIC_HEAT_RATIO=T} on the \ct{MISC} line. For high molecular weight species, use of the default gamma will result in very low values of the specific heat. This can cause issues with temperature in regions with high fuel mass fractions and can cause issues with the \ct{EXTINCTION 2} extinction model. For high molecular weight species it is recommended that the specific heat be defined. If the FDS is determining the specific heat using \ct{GAMMA}, then it is recommended that you check the \ct{CHID.out} file and verify that reasonable specific heat values have been created.
 
 Note that species used in chemical reactions must either be predefined or have an explicitly defined $h(T_{\rm ref})$. The exception is the species defined as \ct{FUEL} on the \ct{REAC} input as long as a \ct{HEAT_OF_COMBUSTION} is defined or the reaction is a simple chemistry reaction where \ct{EPUMO2} applies. Specifying either \ct{SPECIFIC_HEAT} or \ct{RAMP_CP} is considered to have explicitly defined $h(T_{\rm ref})$. If you still wish to use the default \ct{EPUMO2} or 13,100~kJ/kg, then you will need to explicitly define it on the \ct{REAC} input.
 
@@ -6888,12 +6888,12 @@ \subsection{Basic Equations}
 
 \begin{figure}[p]
    \begin{tabular*}{\textwidth}{l@{\extracolsep{\fill}}r}
-      \includegraphics[width=3.2in]{FIGURES/vel_L=-100}  &
-      \includegraphics[width=3.2in]{FIGURES/tmp_L=-100}  \\
-      \includegraphics[width=3.2in]{FIGURES/vel_L=100}  &
-      \includegraphics[width=3.2in]{FIGURES/tmp_L=100}  \\
-      \includegraphics[width=3.2in]{FIGURES/vel_L=-350}  &
-      \includegraphics[width=3.2in]{FIGURES/tmp_L=-350}
+      \includegraphics[width=3.2in]{FIGURES/Monin_Obukhov_vel_L_-100}  &
+      \includegraphics[width=3.2in]{FIGURES/Monin_Obukhov_tmp_L_-100}  \\
+      \includegraphics[width=3.2in]{FIGURES/Monin_Obukhov_vel_L_100}  &
+      \includegraphics[width=3.2in]{FIGURES/Monin_Obukhov_tmp_L_100}  \\
+      \includegraphics[width=3.2in]{FIGURES/Monin_Obukhov_vel_L_-350}  &
+      \includegraphics[width=3.2in]{FIGURES/Monin_Obukhov_tmp_L_-350}
    \end{tabular*}
    \caption[Sample vertical wind and temperature profiles]{Sample vertical wind and temperature profiles.}
    \label{wind_profiles}
@@ -12267,7 +12267,7 @@ \section{\texorpdfstring{{\tt DEVC}}{DEVC} (Device Parameters)}
 \ct{SMOOTHING_TIME}       & Real            & Section~\ref{info:basic_control}                                & s     &               \\ \hline
 \ct{SPATIAL_STATISTIC}    & Character       & Section~\ref{info:statistics}                                   &       &               \\ \hline
 \ct{SPEC_ID}              & Character       & Section~\ref{info:gasoutputquantities}                          &       &               \\ \hline
-\ct{STATISTICS_END}       & Real            & Section~\ref{info:rmscovcorr}                                   & s     & \ct{T_BEGIN}\\ \hline
+\ct{STATISTICS_END}       & Real            & Section~\ref{info:rmscovcorr}                                   & s     & \ct{T_END}\\ \hline
 \ct{STATISTICS_START}     & Real            & Section~\ref{info:rmscovcorr}                                   & s     & \ct{T_BEGIN}\\ \hline
 \ct{SURF_ID}              & Character       & Section~\ref{info:statistics}                                   &       &               \\ \hline
 \ct{TEMPORAL_STATISTIC}   & Character       & Section~\ref{info:statistics}                                   &       &               \\ \hline

Manuals/FDS_User_Guide/FIGURES/Monin_Obukhov_Profiles.py

@@ -0,0 +1,83 @@
+
+# Atmospheric Boundary Layer profiles, based on M-O theory as described in
+# "Falcon Series Data Report", GRI-89/0138, June 1990.
+# These plots are used as illustrations in the FDS User's Guide.
+
+import numpy as np
+import matplotlib.pyplot as plt
+import fdsplotlib
+
+plot_style = fdsplotlib.get_plot_style('fds')
+
+rho_0 = 1.2
+g = 9.81
+cp = 1005
+kappa = 0.41
+p_0 = 100000
+
+L   = np.array([-100,100,-350])
+Ls  = np.array(['-100','100','-350'])
+z_0 = np.array([0.01,0.005,0.005])
+z_r = 2.0
+u_r = 5.0
+
+for i in range(3):
+
+    p_r = p_0 - rho_0 * g * (z_r - z_0[i])
+    T_r = 20 + 273.15
+    theta_r = T_r * (p_0 / p_r) ** 0.285
+
+    if L[i] < 0:
+        x_r = (1 - 16 * z_r / L[i]) ** 0.25
+        psi_m_r = 2 * np.log((1 + x_r) / 2) + np.log((1 + x_r**2) / 2) - 2 * np.arctan(x_r) + np.pi / 2
+        psi_h_r = 2 * np.log((1 + x_r**2) / 2)
+    else:
+        psi_m_r = -5 * z_r / L[i]
+        psi_h_r = psi_m_r
+
+    u_star = kappa * u_r / (np.log(z_r / z_0[i]) - psi_m_r)
+    theta_0 = theta_r / (1 + u_star**2 * (np.log(z_r / z_0[i]) - psi_h_r) / (g * kappa**2 * L[i]))
+    theta_star = u_star**2 * theta_0 / (g * kappa * L[i])
+
+    z = np.array([z_0[i], 10*z_0[i], 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 15, 20, 25, 30, 50, 100])
+
+    # Figure 1 - Velocity profile
+
+    if L[i] < 0:
+        x = (1 - 16 * z / L[i]) ** 0.25
+        psi_m = 2 * np.log((1 + x) / 2) + np.log((1 + x**2) / 2) - 2 * np.arctan(x) + np.pi / 2
+        psi_h = 2 * np.log((1 + x**2) / 2)
+    else:
+        psi_m = -5 * z / L[i]
+        psi_h = psi_m
+
+    u = (u_star / kappa) * (np.log(z / z_0[i]) - psi_m)
+    theta = theta_0 + (theta_star / kappa) * (np.log(z / z_0[i]) - psi_h)
+    T = theta * (p_0 / (p_0 - rho_0 * g * (z - z_0[i]))) ** -0.285
+
+    fig = fdsplotlib.plot_to_fig(x_data=u, y_data=z, marker_style='ko-', 
+                                 x_label='Velocity (m/s)', y_label='Height (m)',
+                                 x_min=0, x_max=10, y_min=0, y_max=30)
+
+    ax1 = fig.axes[0]
+    ax1.text(0.05, 0.90, r'$u(2 \; \mathrm{m})=5$ m/s', transform=ax1.transAxes, verticalalignment='top')
+    ax1.text(0.05, 0.80, f'$L={L[i]:4.0f}$ m', transform=ax1.transAxes, verticalalignment='top')
+    ax1.text(0.05, 0.70, f'$z_0={z_0[i]:5.3f}$ m', transform=ax1.transAxes, verticalalignment='top')
+
+    plt.savefig(f'Monin_Obukhov_vel_L_{Ls[i]}.pdf', format='pdf')
+    plt.close()
+
+    # Figure 2 - Temperature profile
+
+    fig = fdsplotlib.plot_to_fig(x_data=T-273.15, y_data=z, marker_style='ko-',
+                                 x_label='Temperature (°C)', y_label='Height (m)',
+                                 x_min=15, x_max=25, y_min=0, y_max=30)
+
+    ax2 = fig.axes[0]
+    ax2.text(0.05, 0.90, r'$T(2 \; \mathrm{m})=20\;^\circ$C', transform=ax2.transAxes, verticalalignment='top')
+    ax2.text(0.05, 0.80, f'$L={L[i]:4.0f}$ m', transform=ax2.transAxes, verticalalignment='top')
+    ax2.text(0.05, 0.70, f'$z_0={z_0[i]:5.3f}$ m', transform=ax2.transAxes, verticalalignment='top')
+
+    plt.savefig(f'Monin_Obukhov_tmp_L_{Ls[i]}.pdf', format='pdf')
+    plt.close()
+

Utilities/Python/scripts/spray_pattern.py → Manuals/FDS_User_Guide/FIGURES/spray_pattern.py

@@ -7,8 +7,6 @@
 
 plot_style = fdsplotlib.get_plot_style('fds')
 
-pltdir = '../../Manuals/FDS_User_Guide/FIGURES/'
-
 b = [0, 5, 10, 100, 1000]
 mu = np.array([0, 10, 20]) * np.pi / 180
 mu_str = ['0', '10', '20']
@@ -34,6 +32,6 @@
         ax.invert_yaxis()
         ax.text(45, 0.1, r'$\mu=$' + str(mu[j]*180/np.pi) + r', $\phi_{\rm max}=$' + str(phi_max[k]*180/np.pi) + ' degrees', fontsize=plot_style['Key_Font_Size'])
 
-        plt.savefig(pltdir + 'spray_pattern_mu_' + mu_str[j] + '_phimax_' + phi_str[k] + '.pdf', format='pdf')
+        plt.savefig('spray_pattern_mu_' + mu_str[j] + '_phimax_' + phi_str[k] + '.pdf', format='pdf')
         plt.close()
 

Source/part.f90

@@ -1013,9 +1013,9 @@ SUBROUTINE INSERT_VOLUMETRIC_PARTICLES
 ! If the INIT volume is outside the current mesh, return
 
 IF (IN%SHAPE/='RING' .AND. IN%SHAPE/='LINE') THEN
-   IF ((IN_X1-XF)>50._EB*TWO_EPSILON_EB .OR. (IN_X2-XS)<-50._EB*TWO_EPSILON_EB .OR. &
-       (IN_Y1-YF)>50._EB*TWO_EPSILON_EB .OR. (IN_Y2-YS)<-50._EB*TWO_EPSILON_EB .OR. &
-       (IN_Z1-ZF)>50._EB*TWO_EPSILON_EB .OR. (IN_Z2-ZS)<-50._EB*TWO_EPSILON_EB) RETURN
+   IF ((IN_X1-XF)>-50._EB*TWO_EPSILON_EB .OR. (IN_X2-XS)<50._EB*TWO_EPSILON_EB .OR. &
+       (IN_Y1-YF)>-50._EB*TWO_EPSILON_EB .OR. (IN_Y2-YS)<50._EB*TWO_EPSILON_EB .OR. &
+       (IN_Z1-ZF)>-50._EB*TWO_EPSILON_EB .OR. (IN_Z2-ZS)<50._EB*TWO_EPSILON_EB) RETURN
 ELSEIF (IN%SHAPE=='RING') THEN
    IF (RING_MESH_INTERSECTION_ARC(NM,IN%X0,IN%Y0,IN%RADIUS)<TWO_EPSILON_EB) RETURN
 ELSE

Source/pres.f90

@@ -3769,7 +3769,7 @@ SUBROUTINE CHECK_UNSUPPORTED_MESH(SUPPORTED_MESH)
    TRN_ME(2)=TRN_ME(1)
    IF (N_MPI_PROCESSES > 1) CALL MPI_ALLREDUCE(TRN_ME(1),TRN_ME(2),1,MPI_INTEGER,MPI_SUM,MPI_COMM_WORLD,IERR)
    IF (TRN_ME(2) > 0) THEN ! Meshes at different refinement levels. Not Unsupported.
-      IF (MY_RANK == 0) WRITE(LU_ERR,*) 'GLMAT Setup Error : Meshes at different refinement levels unsupported.'
+      IF (MY_RANK == 0) WRITE(LU_ERR,*) 'GLMAT Setup Error: Meshes at different refinement levels currently unsupported for &GEOM.'
       SUPPORTED_MESH = .FALSE.
       STOP_STATUS = SETUP_STOP
       RETURN

Source/read.f90

@@ -12999,12 +12999,12 @@ SUBROUTINE READ_INIT
       MESH_LOOP: DO NM=1,NMESHES
          M=>MESHES(NM)
          IF (XYZ(1)>=M%XS .AND. XYZ(1)<=M%XF .AND. XYZ(2)>=M%YS .AND.  XYZ(2)<=M%YF .AND. XYZ(3)>=M%ZS .AND. XYZ(3)<=M%ZF) THEN
-            IF (ABS(XYZ(1)-M%XS)<TWO_EPSILON_EB) XYZ(1) = XYZ(1) + 0.01_EB*M%DXI
-            IF (ABS(XYZ(1)-M%XF)<TWO_EPSILON_EB) XYZ(1) = XYZ(1) - 0.01_EB*M%DXI
-            IF (ABS(XYZ(2)-M%YS)<TWO_EPSILON_EB) XYZ(2) = XYZ(2) + 0.01_EB*M%DETA
-            IF (ABS(XYZ(2)-M%YF)<TWO_EPSILON_EB) XYZ(2) = XYZ(2) - 0.01_EB*M%DETA
-            IF (ABS(XYZ(3)-M%ZS)<TWO_EPSILON_EB) XYZ(3) = XYZ(3) + 0.01_EB*M%DZETA
-            IF (ABS(XYZ(3)-M%ZF)<TWO_EPSILON_EB) XYZ(3) = XYZ(3) - 0.01_EB*M%DZETA
+            IF (ABS(XYZ(1)-M%XS)<50._EB*TWO_EPSILON_EB) XYZ(1) = XYZ(1) + 0.01_EB*M%DXI
+            IF (ABS(XYZ(1)-M%XF)<50._EB*TWO_EPSILON_EB) XYZ(1) = XYZ(1) - 0.01_EB*M%DXI
+            IF (ABS(XYZ(2)-M%YS)<50._EB*TWO_EPSILON_EB) XYZ(2) = XYZ(2) + 0.01_EB*M%DETA
+            IF (ABS(XYZ(2)-M%YF)<50._EB*TWO_EPSILON_EB) XYZ(2) = XYZ(2) - 0.01_EB*M%DETA
+            IF (ABS(XYZ(3)-M%ZS)<50._EB*TWO_EPSILON_EB) XYZ(3) = XYZ(3) + 0.01_EB*M%DZETA
+            IF (ABS(XYZ(3)-M%ZF)<50._EB*TWO_EPSILON_EB) XYZ(3) = XYZ(3) - 0.01_EB*M%DZETA
             EXIT MESH_LOOP
          ENDIF
       ENDDO MESH_LOOP

Utilities/Matlab/FDS_validation_script.m

@@ -54,6 +54,4 @@
 
 % Miscellaneous other scripts for special cases
 
-backward_facing_step
-
 display('validation scripts completed successfully!')