Ibtissam Hallaci*
| Fares Innal
© 2026 The authors. This article is published by IIETA and is licensed under the CC BY 4.0 license (http://creativecommons.org/licenses/by/4.0/).
OPEN ACCESS
Industrial development, while raising the standard of living, has led to large-scale installations that pose significant risks of major accidents. These facilities handle substantial quantities of hazardous materials under severe operating conditions, such as high temperature and pressure. Consequently, rigorous risk assessment methodologies are essential, requiring the quantification of accident scenarios in terms of both frequency and severity. However, such quantifications are inherently subject to uncertainties in the underlying parameters, which can profoundly impact the accuracy of the results and, ultimately, risk-informed decision-making. Therefore, a structured treatment of these uncertainties is critical for providing a robust foundation for safety measures. This paper has two primary objectives. First, it provides a thorough description of a physical model for quantifying thermal radiation from a pool fire using simplified engineering equations. Second, it demonstrates a methodology for handling uncertainty in the quantification of pool fire safety distances for 3 kW/m2, 5 kW/m2 and 8 kW/m2 thermal radiation. This is achieved by applying a Monte Carlo sampling approach, where input uncertainties are characterized by probability distributions. The proposed framework is demonstrated through a case study of a gasoline pool fire. The results reveal a substantial difference between the average predicted safety distances and the more conservative estimates derived from upper percentiles (e.g., the 95th). Furthermore, a sensitivity analysis is conducted to identify the parameters that contribute most significantly to the variability in these safety distances.
pool fire, thermal radiation thresholds, quantitative risk assessment, uncertainty, Monte Carlo sampling, sensitivity evaluation
Industrial sites are, by their very nature, prone to risks. Serious incidents may occur, including the discharge of hazardous substances and forms of energy, such as thermal radiation, blast waves, and flying debris resulting from fires and explosions. Such events have the potential to cause catastrophic harm, leading to numerous casualties, extensive property destruction, and severe environmental consequences, alongside considerable economic losses. Consequently, there is a growing need to ensure that risks to people, property, and the environment from the activities of potentially hazardous industries are properly assessed and managed. This concern is always high on the agenda in most governments, which establish regulatory documents in order to help site operators demonstrate risk control over potential unwanted events or accidents. For instance, we could mention the Control of Major Accident Hazards (COMAH) regulation applicable in the United Kingdom [1], the Occupational Safety and Health Administration Process Safety Management (OSHA-PSM) [2] and the Environmental Protection Agency Risk Management Plan (EPA-RMP) applicable in the United States [3], and the Seveso Directive which concerns all European Union member states [4]. The two main objectives of such regulations are the provision of all required safeguards to prevent major incidents, and to mitigate their impact should they occur. The achievement of these objectives requires the use of Quantitative Risk Assessment (QRA), an analytical framework designed to evaluate both the likelihood of hazardous events and their tangible impacts, encompassing health, economic and environmental dimensions. Subsequently, the associated risk level could be generally identified as a function of the frequency-severity pair through an appropriate frequency/severity risk matrix. More formally, risk may be understood as the possibility of harm arising from contact with a hazard. The process of risk analysis involves addressing the following key questions [5, 6]:
•What kinds of events might happen (accident scenarios: Si)?
•How often might each one occur (scenario frequency: Fi)?
•If one does occur, what would the effects be (Ci)?
Thus, risk (R) can be stated in quantitative terms as the set of triplets shown below: R = (Si, Fi, Ci) | i = 1, 2, …, n.
Addressing the first question demands technical expertise regarding the potential events or root causes that may give rise to adverse outcomes (i.e., scenarios). For that purpose, qualitative approaches such as HAZard and OPerability study (HAZOP) and Failure Mode and Effect Analysis (FMEA) are generally used. The second question is tackled through the application of probabilistic techniques. In practice, this quantification commonly relies on Boolean logic-based tools such as Event Tree Analysis (ETA), Fault Tree Analysis (FTA), Bow-tie, and Layer of Protection Analysis (LOPA). The third question is answered through deterministic models that describe the physical impacts of the studied accident scenario (heat radiation, blast overpressure, toxic concentration, …). These models are generally incorporated within user-friendly software such as ALOHA [7], PHAST [8], and EFFECTS [9]. Note that more sophisticated models could be developed by the analysts using CFD (computational fluid dynamics) tools.
Within the QRA framework, a major concern is the suitability of the models used and the reliability of the input data, which come into play in the quantification of both frequency and severity. Indeed, uncertainties are inherently present when conducting risk analysis [10]. These arise from incomplete knowledge regarding the studied system, which necessitates the adoption of simplifying assumptions to address this shortcoming. The uncertainty induced by the lack of knowledge is referred to as “epistemic uncertainty” and can be diminished as new information becomes available. Generally, epistemic uncertainty can be traced to three primary causes [11, 12]: completeness uncertainty, model uncertainty, and parameter uncertainty. Completeness uncertainty concerns factors that have not been adequately incorporated in the study (e.g., overlooking certain human errors in frequency determination, or disregarding certain accident scenarios, either deliberately or inadvertently). Model uncertainty is associated with the simplification of the reality that the model is intended to represent (e.g., using simplified formulas or dedicated software to model fires, explosions, and toxic release, or using simplified formulas or dedicated methods to evaluate the accident frequency). Parameter uncertainty arises from a lack of accurate data on the parameter values used in the analysis (for instance, failure rate of a safety barrier, repair duration, ignition probability, wind speed, etc.). This paper is primarily concerned with parameter uncertainty.
More precisely, this study focuses on uncertainty handling in heat radiation threat zones induced by pool fire scenario, where uncertainty is mainly related to weather conditions and fuel properties. It is worth noting that deterministic approaches produce single-value risk indicators that fail to capture this variability, potentially leading to overconfident predictions and suboptimal resource allocation. Uncertainty quantification addresses these limitations by characterizing the distribution of possible outcomes rather than a single prediction. Different techniques are in use when it comes to uncertainty quantification, including interval analysis, probability bound analysis, Monte Carlo sampling, fuzzy sets, and evidence theory. More insights into these approaches are given in Section 3. Within the framework of fire safety engineering, Monte Carlo technique is the most common approach. For instance, the propagation of uncertainty associated with NO2 atmospheric dispersion following a crude oil tank fire has been investigated using an appropriate CFD model [13]. Uncertainties were defined a priori as probability distributions for wind speed, NO2 emission rate, viscosity, and diffusivity coefficients. An LNG release scenario has also been studied by propagating uncertainties in source terms and dispersion models to calculate hazard zones [14]. It is worth noting that Monte Carlo sampling for uncertainty propagation in complex models, such as CFD models, is computationally expensive. Surrogate modeling and Latin Hypercube Sampling (LHS) offer an effective solution to this limitation, thereby enabling fast Monte Carlo sampling [15, 16]. Besides the computational expense of Monte Carlo sampling, representing parameter uncertainty as statistical distributions may be subject to criticism, particularly in the absence of sufficient statistical data. The Bayesian approach offers a complementary or alternative framework that allows new data or observations (new evidence) to be integrated with prior knowledge to produce updated posterior probability distributions, thereby reducing uncertainty iteratively as data become available. Bayesian updating and inference have been applied in different fire safety studies [17, 18].
For parameter uncertainties that cannot be adequately represented using statistical distributions, for example in the case of subjective expert judgments or incomplete data, alternative approaches including fuzzy sets [19] and evidence theory [20] are more suitable. An overview of these and other uncertainty quantification alternatives in the context of fire safety engineering [21].
The main purpose of this study is the deployment of Monte Carlo sampling applied to a gasoline pool fire scenario physical impacts model. The quantification of severity, in terms of heat radiation threat zones, relies on well-established and widely recognized empirical correlations. Accordingly, this paper aims first to provide a thorough description of the physical model, incorporating new literature-based correlations for some parameters for the case of gasoline instead of typical values. More specifically, these new correlations relate to the fraction of the flame surface covered by soot (instead of using a typical value of 0.8) and the radiation fraction (traditionally, 0.3 or 0.35 values are used). Moreover, a correction to an existing view factor formula is provided. Unlike the cylindrical and planar models, the simpler model adopted here accounts for target location (not just flame front distance), enabling uncertainty quantification for the view factor. The second aim of this paper is uncertainty quantification (propagation and sensitivity analysis) for thermal radiation threshold distances using Monte Carlo sampling. Several sensitivity indices are computed, including the newly developed Delta measure. To our knowledge, no study to date has evaluated the impact of uncertainties on thermal radiation threshold distances.
Due to the physical model’s low complexity, advanced sampling methods (e.g., LHS) or surrogate modeling are unnecessary for reducing the computational load. Although statistical distributions are unsuitable for some uncertain parameters, available information does not justify alternatives such as Bayesian, fuzzy sets, or evidence theory. Moreover, sensitivity analysis using these alternatives remains immature.
This paper is structured as follows. Section 2 details the model for quantifying thermal radiation from a pool fire. Section 3 briefly describes the Monte Carlo technique for uncertainty analysis. Section 4 demonstrates the application of this technique to assess uncertainties in estimating thermal radiation threshold distances, using a pool fire scenario caused by the overfilling or leakage of a gasoline storage tank as a case study. The analysis carefully respects the actual variability in meteorological conditions characterizing the case study, specifically temperature, wind speed, and relative humidity. Conclusions and future work are presented in Section 5.
The standard formula for thermal radiation (heat flux) released by a pool fire and received by a target situated x meters away from the flame front is:
$\emptyset(x)=\emptyset_0 \cdot \tau(x) \cdot F(x)$ (1)
where,
$\emptyset(x)$: heat flux at a certain distance x (received heat flux) (kW/m2).
$\varnothing_0$: flame surface emissive power (kW/m2).
$\tau(x)$: atmospheric transmissivity to thermal radiation.
$F(x)$: view factor.
2.1 Flame surface emissive power
The flame surface emissive power ($\varnothing_0$) depends on the size of the fire and the type of fuel. There are several models for calculating this quantity. In this paper, we adopt the correlation provided by the TNO in the Yellow Book [22] that constitutes the base document for the PHAST software.
$\phi_0=\phi_{\max} \cdot(1-\xi)+\phi_{\text {soot}} \cdot \xi$ (2)
where,
$\emptyset_{\max}$: maximum surface emissive power of a flame without soot production $\left(\mathrm{kW} / \mathrm{m}^2\right)$.
$\emptyset_{\text {soot}}$: surface emissive power of soot (black smoke) $\left(\mathrm{kW} / \mathrm{m}^2\right)$.
$\xi$: fraction of the flame surface covered by soot.
Note that the adopted correlation, including its subsequent relations, is the most comprehensive one since it considers the flame length, surface area, and the product type.
It has been observed that $\emptyset_{\text {soot}}$ for pool fires with a diameter of 10 m is $20 \mathrm{~kW} / \mathrm{m}^2$ when the temperature is approximately 800 K [23, 24]. More recently, a value of $40 \mathrm{~kW} / \mathrm{m}^2$ was recorded during experiments involving gasoline and diesel oil [25]. In addition, a value of 0.8 is generally used as a representative figure for oil products for $\xi$. However, the following expression was established for the derivation of $\xi$ in the case of gasoline [25]:
$\xi=\left\{\begin{array}{c}0.55, \text { if } D_{e q}<5 \mathrm{~m} \\ 1-\left(1.80 D_{e q}{ }^{-0.377}-0.533\right), \text { if } 5 \mathrm{~m} \leq D_{e q}<20 \mathrm{~m} \\ 0.95, \text { if } D_{e q} \geq 20 \mathrm{~m}\end{array}\right.$ (3)
$D_{e q}$ refers to the equivalent pool diameter (m) of a confined pool. It can be determined according to two different geometric definitions for the equivalent diameter of a non-circular pool:
$\begin{gathered}D_{e q}=\sqrt{4 . S / \pi} \\ D_{e q}=4 . S / P\end{gathered}$ (4)
where,
S: surface area of the pool (m2). The estimated surface could be taken equal to the whole bund surface (conservative approach), or a net surface is calculated by subtracting the surface of tanks settled within the same bund.
P: fire (pool) perimeter (m).
The quantity $\sqrt{4 . S / \pi}$ calculates the diameter of a circle whose area is equal to the pool surface area $S$ (areaequivalent diameter), whereas $4 . S / P$ represents the hydraulic diameter, which is the characteristic length based on perimeter.
$\emptyset_{\max }$ is calculated using the Moorhouse and Pritchard approach [26]:
$\emptyset_{\max}=f_{r a d} \cdot \ddot{m} \cdot \Delta H_C \cdot \frac{1}{1+4 L_F / D_{e q}}$ (5)
where,
$f_{\text {rad}}$: radiation fraction, i.e., the fraction of the combustion energy radiated from the flame surface.
$\ddot{m}$: mass burning rate per unit area at still weather conditions $\left(\mathrm{kg} /\left(\mathrm{m}^2 . \mathrm{s}\right)\right)$.
$\Delta H_C$: heat of combustion of the flammable material at its boiling point (kJ/kg).
$L_F$: flame length (m).
A value of 0.35 is considered conservative for $f_{\text {rad}}$ [27]. This value is used within the RAST tool [28]. Moreover, ALOHA software [7] approximates $f_{\text {rad}}$ as a constant 0.3 , as suggested by Roberts [29]. However, these values are adequate for pool fires with small diameters. In fact, there is a significant dependence between the value of $f_{\text {rad}}$ and pool diameter [30]. In general, experiments show that $f_{\text {rad}}$ associated with large pool fires is relatively low (less than 0.1) [31]. In this context, Muñoz et al. [25] developed the following correlation for gasoline and diesel oil:
$f_{r a d}=\left\{\begin{array}{l}0.158 D_{e q}^{0.15} \text { if } D_{e q} \leq 5 \mathrm{~m} \\ 0.436 D_{e q}^{-0.58} \text { if } D_{e q}>5 \mathrm{~m}\end{array}\right.$ (6)
The mass burning rate per unit area $\ddot{m}$ depends on the material and the pool diameter, and is determined by the subsequent correlation [32]:
$\ddot{m}=\ddot{m}_{\infty} .\left(1-e^{-k \beta \cdot D_{e q}}\right)$ (7)
$\ddot{m}_{\infty}$ represents the mass burning rate for an infinitely large pool in $\mathrm{kg} /\left(\mathrm{m}^2 \mathrm{~s}\right), k$ denotes the absorption extinction coefficient of the flame $\left(\mathrm{m}^{-1}\right)$ and $\beta$ is the mean beam length corrector. For gasoline, $m_{\infty}=0.055 \pm 0.002 \mathrm{~kg} /\left(\mathrm{m}^2 \mathrm{~s}\right)$ and $k \beta=2.1 \pm 0.3 \mathrm{~m}^{-1}[32]$.
The flame length $\left(L_F\right)$ is determined using the well-established Thomas correlations [33, 34]:
$L_F=\left\{\begin{array}{l}42 \cdot D_{e q} \cdot\left(\frac{\ddot{m}}{\rho_{a i r} \cdot \sqrt{g \cdot D_{e q}}}\right)^{0.61} \text { if } u_w<1 . \\ 55 \cdot D_{e q} \cdot\left(\frac{\ddot{m}}{\rho_{a i r} \cdot \sqrt{g \cdot D_{e q}}}\right)^{0.67} \cdot\left(u^*\right)^{-0.21}, \\ \text { otherwise.} \\ \text { With: } u^*=\max \left\{1, \frac{u_w}{\left(\frac{g \times \ddot{m} \times D_{e q}}{\rho_{\text {air }}}\right)^{\frac{1}{3}}}\right\}\end{array}\right.$ (8)
where,
$u_w$: wind speed $(\mathrm{m} / \mathrm{s})$.
$g$: gravitational acceleration $\left(9.80665 \mathrm{~m} / \mathrm{s}^2\right)$.
$u^*$: scaled wind velocity.
$\rho_{\text {air}}$: air density $\left(\mathrm{kg} / \mathrm{m}^3\right)$. It is calculated in this paper using the subsequent equation:
$\rho_{\text {air}}=\frac{1}{287.06(\vartheta+273.15)}\left(p-230.617 . R H . e^{\left(\frac{17.5043 \vartheta}{241.2+\vartheta}\right)}\right)$ (9)
where,
$\vartheta$: ambient temperature $\left({ }^{\circ} \mathrm{C}\right)$.
$R H$: relative humidity (%).
$p$: barometric pressure ( Pa ). Its useful mathematical form with respect to altitude $h(\mathrm{m})$ is given as follows [35]:
$p(h)=p_0 e^{-\frac{M g h}{R(\vartheta+273.15)}}$ (10)
where, $p_0$ is the atmospheric pressure at sea level (101325 Pa), R is the gas constant (8314.462 j/(kmol·K)), and M is the molar mass of Earthly air (28.9644 kg/kmol).
2.2 Atmospheric transmissivity
The atmospheric transmissivity $\tau(x)$ reflects the fact that the emitted radiation is partially absorbed by the air between the radiating surface (flame front) and the target. Different correlations have been established for the estimation of this factor. The Bagster and Pitblado correlation [36] is adopted in this paper, as it accounts for relative humidity and ambient temperature in a straightforward manner and is the correlation embedded within the PHAST software.
$\tau(x)=2.02\left(P_w \cdot x\right)^{-0.09}=2.02\left(P_s \cdot R H \cdot x\right)^{-0.09}$ (11)
$P_w$ refers to the water vapor partial pressure. $P_s$ is the saturation water vapor pressure (Pa), which is calculated here using the following approximation (in Pa) [37]:
$P_s=\left\{\begin{array}{l}610.5 e^{\left(\frac{17.269 \vartheta}{237.3+\vartheta}\right)} \text {if} \vartheta \geq 0^{\circ} \mathrm{C} \\ 610.5 e^{\left(\frac{21.875 \vartheta}{265.5+\vartheta}\right)} \text { if } \vartheta<0^{\circ} \mathrm{C}\end{array}\right.$ (12)
2.3 View factor
The view factor $F(x)$ characterizes the geometric arrangement between the flame and the target surfaces. It represents the fraction of the target’s field of view that is occupied by the flame [22] (Figure 1). Thus, it is purely geometric and depends only on the relative arrangement of the two surfaces and their respective geometry. The view factor translates the fraction of the energy emitted by a fire and which is received by a target.
Figure 1. Angular configuration between flame and target
The available models are the cylindrical and planar view factors depending on the pool geometry. For instance, PHAST software implements cylindrical correlation regardless of the pool fire geometry, whereas ALOHA software employs a more accurate and complex approach based on numerical integration. In this paper, for a rectangular pool, the simplified formula given by Eq. (13) is retained [38] (Figure 2), in which the flame is treated as a point source, such that the radiative heat flux at a target located at distance $r$ from the flame center and angle $\theta_1$ measured from the flame surface normal (with $\theta_2=\pi / 2-\theta_1$ ) decays as $1 / r^2$. The angle $\theta_1$ varies from 0 (target aligned with the midpoint of the bund length) to $\pi / 2$ (target aligned with the midpoint of the bund width). All other target locations are covered by this range, with the four corners of the rectangle being symmetric.
Figure 2. Pool fire dimensions for the derivation of the view factor
The rationale for adopting this model lies in its capacity to consider target location (rather than merely flame front distance), which is essential for quantifying view factor uncertainties. By contrast, alternative correlations only consider target locations that maximize the view factor.
$F(x)=\frac{\cos \theta_1 \cdot L_L \cdot L_F}{\pi r^2}+\frac{\cos \theta_2 \cdot L_W \cdot L_F}{\pi r^2}=\frac{\cos \theta_1 \cdot L_L \cdot L_F}{\pi r^2}+\frac{\sin \theta_1 \cdot L_W \cdot L_F}{\pi r^2}=\frac{L_F}{\pi r^2}\left(\cos \theta_1 \cdot L_L+\sin \theta_1 \cdot L_W\right)$ (13)
where,
$L_F$: flame length (m).
$L_L$: bund length (m).
$L_W$: bund width (m).
r:distance between the center of the pool and the target.
$0 \leq \theta_1 \leq \frac{\pi}{2} ; \theta_2=\frac{\pi}{2}-\theta_1$.
In the reference [38], an ambiguity exits relating to the parameter r. In fact, it is used as the distance between the center of the pool and the target in the estimation of the view factor and as the distance from the flame front to the target when calculating the atmospheric transmissivity. To remedy this confusion, based on Figure 2, we write r as a function of x (the distance between the flame front and the target, and x being contained in segment r.) under the following form:
$r=\left\{\begin{array}{l}\frac{L_W}{2 \cos \theta_1}+x \text { if } 0 \leq \theta_1 \leq \tan ^{-1}\left(\frac{L_L}{L_W}\right) \\ \frac{L_L}{2 \sin \theta_1}+x \text { if } \tan ^{-1}\left(\frac{L_L}{L_W}\right)<\theta_1 \leq \frac{\pi}{2}\end{array}\right.$ (14)
Note that $r$ should not take values lower than those given by Eq. (14) where $x=0$, i.e., the target is located outside the fire. If lower values are used, the view factor would be greater than the unity. Moreover, for a given $r$, the maximum view factor is obtained for $\theta_1=\tan ^{-1}\left(\frac{L_W}{L_L}\right)\left(\theta_1\right)$ that maximizes $\cos \theta_1 \cdot L_L+\sin \theta_1 \cdot L_W$ in Eq. (13).
From the description of the adopted model, it can be seen that most of the parameters used are uncertain. For instance, weather conditions (ambient temperature, relative humidity, and wind speed) depend on the moment the accident occurs. In addition, the literature shows variability in certain product characteristics, such as heat of combustion, soot surface emissive power, fraction of the flame surface covered by soot, and radiation fraction. The position of the target relative to the pool fire also affects the received heat radiation. Therefore, these uncertainties must be accounted for to obtain reliable results. The purpose of the next Section is to present the Monte Carlo method, which is one of the most widely used techniques in the field of uncertainty quantification.
3.1 Introduction
Uncertainty quantification enables model users to better assess the confidence that can be attributed to model outputs [39]. The Monte Carlo method has become the industry standard for propagating uncertainties [40], as it offers an efficient and straightforward approach based on sampling from probability density functions (PDFs) assigned to uncertain inputs. The main advantages of this technique can be summarized as follows [41]:
•Adaptable, as it handles model complexity and nonlinearity without resorting to simplifying approximations.
•Simple to implement and communicate.
•Able to incorporate information on correlations among variables.
However, MC technique requires substantial empirical data to define credible probability distributions. When such data are lacking, the analyst must make subjective assumptions about distributions and dependencies, which can lead to risk misestimation. Additionally, it is computationally intensive due to repeated random sampling and model evaluation. As noted in Section 1, this computational burden can be substantially reduced through advanced sampling techniques and surrogate modeling.
It is worth mentioning that several alternative approaches for incorporating uncertainty into risk assessment exist. These are briefly summarized in terms of their strengths and drawbacks [42].
•Interval analysis: although simple to implement thanks to interval arithmetic, this approach cannot account for detailed empirical information when available. Moreover, it may yield overly conservative results.
•Probability bounds analysis (PBA): generalizing both interval analysis and probability theory, PBA is computationally faster than Monte Carlo simulation. However, like interval analysis, it may overestimate uncertainty. Furthermore, probability bounds do not indicate the most likely values, as they provide no second-order information.
•Evidence theory: this method is suitable for representing situations where the available knowledge is more detailed than an interval but insufficient to define a unique probability distribution. It does not lead to overestimated output results. Nevertheless, as with PBA, it does not identify the most likely values. Further research is still required to make evidence theory more widely accepted and operational within the safety engineering framework.
•Possibility theory: based on fuzzy numbers (possibility distributions), this approach produces a full distribution instead of a simple range. Fuzzy numbers generalize and refine intervals, whose bounds depend on a chosen confidence level. Additionally, possibility theory addresses non-statistical (non-probabilistic) uncertainty, making it versatile for different uncertainty contexts. However, it may be overly conservative and implicitly assumes total dependence between uncertain variables.
Besides the above-mentioned drawbacks, particularly the inability to identify the most likely values, one should also note the considerable difficulty of implementing some of these approaches on the presented pool fire model. Furthermore, sensitivity analysis (which seeks to apportion variations in model output to changes in input parameters) remains, so far, completely immature within these alternative frameworks.
3.2 Uncertainty quantification
The whole uncertainty quantification process involves two portions: uncertainty propagation and sensitivity analysis, as illustrated in Figure 3. They are briefly presented below based on the MC technique.
Figure 3. Parameter uncertainty framework
3.2.1 Uncertainty propagation
This part addresses the influence of potential inaccuracies in model inputs on the model outputs (e.g., uncertainty in wind speed uw). The key stages of uncertainty propagation for a given model are described below:
•Define a PDF for each uncertain input parameter using available information. Suitable distributions include uniform, triangular, normal, lognormal, Chi-square, beta, and gamma, depending on the extent of information available for each parameter.
•Generate one input parameter set by sampling from the assigned PDFs using random numbers uniformly distributed in [0, 1], as illustrated in Figure 4.
•Evaluate the output function using the sampled random values. In this study, the output functions correspond to the distances from the flame front to the target for the regulatory thermal radiation thresholds of 3 KW/m2, 5 KW/m2 and 8 KW/m2.
•Repeat steps 2 and 3 for n iterations (typically 10,000 or more) to produce n independent realizations for each output function. The resulting sample approximates the probability distribution of the output variable.
•Compute summary statistics from the obtained sample for each output, including the mean $(\bar{X})$, standard deviation $\sigma$, and confidence intervals (percentiles).
Figure 4. Sampling of input values process
3.2.2 Sensitivity analysis
Sensitivity analysis examines how model output variability can be attributed to input parameter variations [43]. This allows users to identify which parameters most significantly affect model predictions. The main sensitivity analysis approaches are screening, local, and global methods [43].
Screening tools can be rapidly applied to identify the most influential input parameters prior to conducting a more detailed analysis. They are especially valuable for computationally intensive models. A simple yet effective visual screening tool consists of scatter plots showing the relationship between the output and each individual input.
Local sensitivity analysis involves choosing a base case set of input values and perturbing each input one-at-a-time (OAT) by a fixed percentage while the remaining inputs are held unchanged. The results are strongly influenced by the selected base values.
Global sensitivity analysis accounts for the full variation range of input parameters through their probability distributions. Among the various global sensitivity indices are the Partial Rank Correlation Coefficient (PRCC), which uses rank transformations, the Fourier Amplitude Sensitivity Test (FAST), and Sobol’s method that is based on variance decomposition.
4.1 Description
The aim of this case study is to evaluate thermal radiation threat zones associated with a gasoline storage tank bund fire, while considering parameter uncertainty. It is worth noting that there are other likely accidents relating to this type of storage, including: tank explosion, vapor cloud explosion (UVCE), flash fire, and toxic gas dispersion. The gasoline tank is located at the RA1K refinery in Skikda (Algeria), situated 30 meters above sea level. The analyzed scenario may originate from an uncontrolled filling operation leading to the overfilling of the tank, or from a tank wall leakage. Subsequent ignition of the released product escalates into a large-scale bund fire, as illustrated in Figure 5. The bund is square with a side length of 90 m.
Figure 5. Illustration of bund fire development
The thermal radiation thresholds for human impacts are gathered in Table 1. The adopted reference threshold values applicable to classified facilities, for which regulatory bodies require the derivation of threat zones, are set at 3 kW/m² (irreversible effects), 5 kW/m² (first lethal effects), and 8 kW/m² (significant lethal effects).
Table 1. Thermal radiation thresholds for human effects
|
Heat Flux (KW/m²) |
Effect |
Typical Consequence |
|
1.5 |
Pain threshold (~60 s) |
Possible escape without injury |
|
3 |
Threshold of irreversible effects |
Serious burns after short exposure; difficult escape without protection |
|
5 |
Threshold of first lethal effects |
Fatalities possible for exposure > 20–30 s |
|
8 |
Threshold of significant lethal effects |
High probability of fatalities for short exposure < 20 s |
|
10 |
Certain lethality within seconds |
100% fatality zone; immediate lethal effect |
4.2 Baseline configuration
Prior to accounting for parameter uncertainty, a baseline simulation is conducted using constant parameter values. The resulting predictions are compared against those derived from the ALOHA software. This preliminary evaluation is only intended to provide additional context for the results under identical input conditions. The following parameters are used: $\vartheta=20^{\circ} \mathrm{C}, \mathrm{RH}=70 \%, u_w=5 \mathrm{~m} / \mathrm{s}, \rho_{\text {air }}=1.161 \mathrm{~kg} / \mathrm{m}^3, \xi=0.0$ (no soot); $D_{e q}=90 \mathrm{~m}, f_{r a d}=0.3, \ddot{m}=0.083 \mathrm{~kg} /\left(\mathrm{m}^2 \mathrm{~s}\right), \Delta H_C=47900 \mathrm{~kJ} / \mathrm{kg}, \theta_1=\frac{\pi}{2}$ (the target is facing the midpoint of a given side of the pool). The parameter values for $\ddot{m}$ and $\Delta H_C$ are based on N-octane, as it is the chemical available in ALOHA whose properties most closely match those of gasoline. Gasoline itself is not a selectable chemical in the ALOHA software. Table 2 presents the results obtained. Figure 6 shows the thermal radiation threat zones generated by the ALOHA software, while Figure 7 illustrates the decay in radiation intensity with distance, as calculated using the adopted pool fire model. Note that the distances shown in Figure 6 are measured from the bund center. Therefore, a reduction of 45 m $\left(D_{e q} / 2\right)$ is required to obtain the distance from the flame front. The distances obtained from the implemented model are in excellent agreement with those generated by the ALOHA software. While this finding provides a useful cross-check that supports our results, it should not be interpreted as a formal validation.
Table 2. Thermal radiation threshold distances (m)
|
Thermal Radiations |
Implemented Model |
ALOHA |
|
3 kW/m2 |
308.45 |
309 |
|
5 kW/m2 |
232.30 |
233 |
|
8 kW/m2 |
176.93 |
178 |
Figure 6. Thermal radiation threat zones (ALOHA)
Figure 7. Thermal radiation intensity profile
4.3 Uncertainty consideration
This subsection considers potential uncertainties associated with the model input parameters. The different parameters used in the model and their related uncertainties are detailed in Table 3. The parameter bounds are based on values reported in the literature and on those commonly adopted in practice. Table 4 presents the numerical values of average ambient temperature $(\vartheta)$, average relative humidity $(R H)$ and wind speed $\left(u_w\right)$. The $\vartheta$ and $R H$ data are based on historical records from the Skikda local weather station [44], whereas $u_w$ is derived from Meteoblue simulation/reanalysis platform [45] and are not direct measurements. The Meteoblue data are based on the ERA5 model, a global reanalysis dataset produced by the European Centre for Medium-Range Weather Forecasts (ECMWF) [46]. Since these parameters are correlated, they are provided on a monthly basis rather than as value ranges. Figure 8 illustrates the variation of these three parameters throughout the year [44].
Table 3. Input parameters data with uncertainty
|
Input Parameter |
Definition |
Value |
Comments |
|
Lw |
Bund width (m) |
90 |
Constant value depending on the bund dimensions |
|
Lt |
Bund length (m) |
90 |
Constant value depending on the bund dimensions |
|
h |
Altitude (m) |
30 |
Constant value depending on the geographic location |
|
g |
Gravitational acceleration (m/s2). |
9.80665 |
Constant value |
|
R |
Gas constant (j/(kmol·K)) |
8314.462 |
Constant value |
|
p0 |
Atmospheric pressure at sea level (Pa) |
101325 |
Constant value |
|
Deq |
Equivalent confined pool diameter (m) |
Deq ∼ Discrete {(101.554, 0.5), (90, 0.5)} |
Discrete {(v1, p), (v2, 1-p)}, means that the random variable = v1 with probability p and = v2 with probability 1-p. The values v1 and v2 are obtained according to Eq. (4). This discrete distribution allows equivalent selection between the two possible modeling choices given by Eq. (4) (Deq model uncertainty) |
|
ΔHC |
Flammable material combustion at its boiling point (kJ/kg) |
ΔHC ∼ U ( 42500 , 47300) |
U (a, b) denotes a uniform distribution with lower bound a and upper bound b. Gasoline is a complex mixture of various hydrocarbons, and their specific blend affects the overall combustion heat. The interval [42.5, 47.3] MJ/Kg covers the entire range of possible gasoline mixtures |
|
$\emptyset_{ {soot}}$ |
Surface emissive power of soot (black smoke) (kW/m2) |
$\emptyset_{ {soot}}$ ∼LN (3.343, 0.211) |
LN (µ, σ) stands for the log-normal distribution, where µ and σ represent the mean and standard deviation of the variable’s natural logarithm, respectively. This distribution is common in reliability and safety studies for uncertainty modeling. It is skewed to the right and therefore favors higher values. A heat radiation of 20 kW/m² is generally used [23, 24], whereas 40 kW/m² was recently registered during experiments for gasoline and diesel oil [25]. Hence, we consider the range of $\emptyset_{ {soot}}$to be [20, 40] kW/m². By choosing the LN distribution, we give greater weight to the higher values within this interval. The corresponding µ and $\sigma$ are derived using the following relations: $\left\{\begin{array}{l}40=e^{(\mu+1.645 \sigma)} \\ 20=e^{(\mu-1.645 \sigma)}\end{array}\right.$ |
|
$\xi$ |
Fraction of the flame surface covered by soot |
$\xi$ ∼ U (0.8, 0.95) |
0.8 is the representative figure usually used for oil products, whereas 0.95 is the value given in [25] for gasoline when $D_{e q} \geq 20 \mathrm{~m}$ (see Eq. (3)) |
|
$f_{\text {rad}}$ |
Radiation fraction |
$f_{\text {rad}}$ ∼ U ($f_{\min}$ , 0.1) |
$f_{\min}$ is derived from Eq. (6), depending on the value of Deq , while 0.1 is taken based on the study [31] |
|
$\ddot{m}$ |
Mass burning rate (kg/(m2 s)) |
$\ddot{m}$ ∼ U (0.053, 0.057) |
Given Eq. (7) and the value of $D_{e q}$: $\ddot{n}=\ddot{m}_{\infty}=0.055 \pm 0.002 \mathrm{~kg} /\left(\mathrm{m}^2 \mathrm{~s}\right)[32]$ |
|
$u_w$ |
wind speed (m/s) |
$u_w \sim\left\{U\left(a_k, b_k\right), p_{i k}\right\}$ |
The wind speed varies throughout the year. For each month i, its distribution is represented by a five-class mixture of uniform distributions $U\left(a_k, b_k\right), k=1, \ldots, 5$, with associated probabilities $p_{i k}$ (Table 4 [45]) |
|
$\vartheta$ |
Ambient temperature (℃) |
$\vartheta \sim T\left(\vartheta_{a i}, \vartheta_{c i}, \vartheta_{b i}\right)$ |
$T(a, \mathrm{c}, b)$ stands for Triangular distribution with lower bound $a$, the mode $c$ and upper bound $b$. The values for $\vartheta_{a i}, \vartheta_{c i}$ and $\vartheta_{b i}$ depend on the $i$-th month and are gathered in Table 4 [44]. Note that the monthly average temperature $\vartheta_{c i}$ was used as an approximation of the mode, which is a common simplification when the true mode is unavailable |
|
$R H$ |
Relative humidity (%) |
$R H \sim$ Discrete $\left\{\left(R H_i, i\right)\right\}$ |
RH depends on the month of the year. $R H_i$ denotes the relative humidity of the i-th month. Table 4 summarizes these different values [44] |
|
$\theta_1$ |
See Figure 2 |
$\theta_1 \sim U\left(0, \frac{\pi}{2}\right)$ |
Depends on the position of the target |
Table 4. Ambient temperature ($\vartheta$), relative humidity (RH), and wind speed (uw) variations
|
Weather Parameters |
Jan |
Feb |
Mar |
Apr |
May |
Jun |
Jul |
Aug |
Sep |
Oct |
Nov |
Dec |
|
Average high temperature (℃) |
16.60 |
16.85 |
19.11 |
20.66 |
23.72 |
27.06 |
29.50 |
30.07 |
27.56 |
25.72 |
21.00 |
17.75 |
|
Average temperature (℃) |
12.50 |
12.80 |
14.80 |
16.70 |
19.80 |
22.80 |
25.50 |
25.30 |
24.00 |
22.00 |
17.00 |
13.50 |
|
Average low temperature (℃) |
8.82 |
8.87 |
10.86 |
12.83 |
15.84 |
19.27 |
22.36 |
22.96 |
20.31 |
17.77 |
13.15 |
8.82 |
|
Average relative humidity (%) |
77.33 |
74.33 |
74.33 |
74.16 |
72.83 |
73.50 |
73.83 |
71.83 |
73 |
75.33 |
72.33 |
75.66 |
|
Wind speed (km/h) |
|
|
|
|
|
|
|
|
|
|
|
|
|
5–10 |
0.135 |
0.067 |
0.042 |
0.010 |
0.010 |
0.004 |
0.000 |
0.000 |
0.013 |
0.061 |
0.117 |
0.148 |
|
10–20 |
0.548 |
0.628 |
0.674 |
0.757 |
0.823 |
0.890 |
0.903 |
0.916 |
0.870 |
0.806 |
0.603 |
0.542 |
|
20–30 |
0.252 |
0.230 |
0.229 |
0.203 |
0.152 |
0.103 |
0.097 |
0.084 |
0.113 |
0.119 |
0.233 |
0.248 |
|
30–40 |
0.061 |
0.071 |
0.055 |
0.030 |
0.015 |
0.003 |
0.000 |
0.000 |
0.004 |
0.014 |
0.040 |
0.058 |
|
40–50 |
0.004 |
0.004 |
0.000 |
0.000 |
0.000 |
0.000 |
0.000 |
0.000 |
0.000 |
0.000 |
0.007 |
0.004 |
Figure 8. Monthly variation of average high, mean, and low temperatures, relative humidity, and wind speed
Note that in Table 4, for wind speed, the numerical values correspond to the probabilities of belonging to the respective wind speed classes. To capture seasonal dependencies among weather parameters, the sampling process starts by randomly choosing a month. Thereafter, each weather parameter is sampled separately according to the selected month and the corresponding probability distribution. For instance, for January (month 1): RH = 77.33%, $\vartheta$ ∼T (8.82, 12.50, 16.60), and wind speed follows a mixture of five classes with probabilities 0.135, 0.548, 0.252, 0.061, and 0.004. A uniform random number r $\in$ [0, 1] selects the wind class as follows: [0, 0.135] → U (5, 10); [0.135, 0.683] → U (10, 20); …;[0.966, 1] → U (40, 50). Once the class is selected, a value is sampled uniformly within it and converted to m/s before model evaluation.
The uniform distribution is predominantly adopted in this study, as it is the logical choice when only bounds are known, assigning equal probability across the interval and reflecting maximum uncertainty.
4.4 Obtained results
The distances to the thermal radiation thresholds are computed using the pool fire model presented in Section 2. The Monte Carlo technique is employed to propagate the uncertainties associated with the input parameters and to analyze their effects on the output variables. A total of 104 trials have been performed. A Python code was developed to perform all required calculations and generate the corresponding graphics. The variation in thermal radiation as a function of distance x is shown in Figure 9, along with 5th–95th percentile uncertainty band (blue shaded area). It is clearly shown that the radiation profile $\emptyset(x)$ is strictly monotonically decreasing with distance.
Figure 9. Mean radiation intensity profile with confidence interval
The safety distances corresponding to thermal radiation thresholds (3, 5, and 8 kW/m²) are computed numerically by finding the roots of the radiation intensity profile using Brent’s method, selected for its robustness and super-linear convergence. The resulted PDF and cumulative distribution function (CDF) are depicted in Figures 10–12. The derived means, standard deviations (std) and various percentiles values (5th, 25th, 50th, 75th and 95th) are gathered in Table 5. These measures help describe the distribution and variability of the distance values.
Figure 10. Probability density function (PDF) and cumulative distribution function (CDF) for 8 kW/m2
Figure 11. Probability density function (PDF) and cumulative distribution function (CDF) for 5 kW/m2
Figure 12. Probability density function (PDF) and cumulative distribution function (CDF) for 3 kW/m2
Table 5. Statistics on radiation thresholds distances (m)
|
Radiations |
Mean |
Std |
5th |
25th |
50th |
75th |
95th |
|
3 kW/m2 |
77.89 |
11.64 |
59.65 |
69.84 |
77.43 |
85.19 |
97.74 |
|
5 kW/m2 |
50.89 |
9.12 |
36.52 |
44.57 |
50.62 |
56.60 |
66.48 |
|
8 kW/m2 |
31.45 |
7.34 |
19.80 |
26.47 |
31.23 |
36.09 |
43.96 |
Note that the obtained mean distances are significantly lower than those from the baseline case. Besides the variability of input parameters, this reduction is primarily attributable to the high radiation attenuation effect of soot $\left(\xi \in[0.8,0.95]\right.$ and $\left.\emptyset_{\text {soot}}=[20,40] \mathrm{kW} / \mathrm{m}^2\right)$, combined with the low radiation fraction $\left(f_{\text {rad}} \in[0.03,0.1]\right)$, see Eq. (2). This leads to a substantial drop in the flame surface emissive power ( $\emptyset_0$ ), whose maximum attainable value is only $66.42 \mathrm{~kW} / \mathrm{m}^2$. By contrast, under the baseline configuration ( $\xi=0$ and $\left.f_{\text {rad }}=0.3\right), \emptyset_0$ is much higher, with $\emptyset_0=\emptyset_{\text {max}}=254.25 \mathrm{~kW} / \mathrm{m}^2$. It should be noted, however, that the two configurations should not be directly compared, as they rely on different assumptions.
When accounting for uncertainties, analysts should consider distance values that provide greater assurance of safety, as relying solely on mean values can be misleading. It is therefore advisable to use higher percentiles, such as the 95th percentile, which indicates that 95% of observed distances fall below it. This provides a more conservative estimate that better accounts for data variability and uncertainty, offering a more robust basis for decision-making. For example, at a thermal radiation level of 3 kW/m², the mean distance is 77.89 m, while the 95th percentile distance is 97.74 m (a substantial difference of 20 m).
Using the 95th percentile distance ensures a more cautious and safer design margin. This approach has direct practical implications for safety management and regulatory compliance. Specifically, the 95th percentile distance can be directly integrated into safety cases by specifying it as a design-basis consequence distance in bow-tie analysis and ALARP demonstrations, as well as for establishing emergency response zones. For land-use planning, regulators could adopt the 95th percentile (rather than the mean) for defining exclusion zones and setback requirements, ensuring risk contours account for input variability. Moreover, higher values significantly raise the likelihood of a domino effect [47]. This uncertainty-aware distance approach provides a transparent and conservative benchmark that supports regulatory validation and strengthens the evidence base for safety decisions.
As a preliminary assessment of the sensitivity of distances associated with radiation thresholds to various input parameters, scatter plots are presented in Figure 13. These plots provide a clear and direct visual representation of the sensitivity relationships. Furthermore, the contribution of each input parameter to the uncertainty of the studied outputs is quantified using several traditional sensitivity indices: SPEA (Spearman Rank Correlation Coefficient), PEAR (Pearson Product-Moment Correlation Coefficient), SRC (Standardized Regression Coefficient), SRRC (Standardized Rank Regression Coefficient), and PRCC (Partial Rank Correlation Coefficient). In addition to these correlation- and regression-based indices, a global sensitivity measure is incorporated: Delta Moment-Independent Measure [48, 49]. This index evaluates the impact of an input variable on the output variance without requiring specific variance decomposition. It is important to note that traditional global sensitivity indices such as Sobol and FAST are not applicable in this study due to the existing correlations among ambient temperature (ϑ), relative humidity (RH), and wind speed (uw). Table 6 summarizes all sensitivity indices for the 3 kW/m² threshold distance. The same ranking pattern of Table 7 was observed for the 5 kW/m² and 8 kW/m² thresholds, indicating consistent parameter importance across radiation levels.
Table 6. Sensitivity indices for 3 kW/m² threshold distance
|
Sensitivity Indices |
∆HC |
$\emptyset_{\text {soot}}$ |
ξ |
$\ddot{m}$ |
θ1 |
frad |
uw |
ϑ |
RH |
|
SPEA |
3.042E-02 |
8.776E-01 |
-1.165E-01 |
8.925E-02 |
1.341E-02 |
8.346E-02 |
-2.206E-01 |
-4.940E-02 |
1.706E-02 |
|
PEAR |
3.068E-02 |
8.895E-01 |
-1.061E-01 |
8.842E-02 |
1.269E-02 |
8.202E-02 |
-2.311E-01 |
-5.160E-02 |
2.498E-02 |
|
SRC |
3.279E-02 |
8.940E-01 |
-1.126E-01 |
8.944E-02 |
-2.348E-03 |
8.453E-02 |
-2.532E-01 |
-8.733E-02 |
-1.040E-02 |
|
SRRC |
3.347E-02 |
8.824E-01 |
-1.199E-01 |
9.492E-02 |
-4.535E-03 |
8.937E-02 |
-2.384E-01 |
-7.159E-02 |
-1.323E-02 |
|
PRCC |
1.561E-01 |
9.264E-01 |
-2.451E-01 |
3.090E-01 |
5.879E-02 |
2.969E-01 |
-4.944E-01 |
4.130E-03 |
1.437E-01 |
|
Delta |
4.463E-02 |
4.880E-01 |
7.334E-02 |
6.008E-02 |
1.058E-01 |
6.921E-02 |
9.280E-02 |
5.613E-02 |
5.229E-02 |
Table 7. Parameter importance ranking for 3 kW/m² threshold distance
|
Sensitivity Indices |
∆HC |
$\emptyset_{\text {soot}}$ |
ξ |
$\ddot{m}$ |
θ1 |
frad |
uw |
ϑ |
RH |
|
SPEA |
7 |
1 |
3 |
4 |
9 |
5 |
2 |
6 |
8 |
|
PEAR |
7 |
1 |
3 |
4 |
9 |
5 |
2 |
6 |
8 |
|
SRC |
7 |
1 |
3 |
4 |
9 |
6 |
2 |
5 |
8 |
|
SRRC |
7 |
1 |
3 |
4 |
9 |
5 |
2 |
6 |
8 |
|
PRCC |
6 |
1 |
5 |
3 |
8 |
4 |
2 |
9 |
7 |
|
Delta |
9 |
1 |
4 |
6 |
2 |
5 |
3 |
7 |
8 |
At first glance, from scatter plots of Figure 13, the parameter $\emptyset_{\text {soot}}$ is the most important one in influencing the three outputs as it imparts more “shape” on them. This initial visual assessment is formally confirmed with the derived ranking of Table 7 based on absolute values of Table 6.
Note that the negative values in Table 6 mean that the increase of the input parameter leads to a decrease in the output measure. The inspection of the ranking table shows that $\emptyset_{\text {soot }}$ is overwhelmingly dominant (rank 1 across all used indices), indicating that soot radiation is by far the most influential parameter controlling the thermal radiation hazard distance. $u_w$ is ranked 2 in five of the six methods (classical indices), and 3 in Delta method. This consistency confirms that wind speed is the second most critical parameter, as it affects the flame length. $\xi, \ddot{m}$ and $f_{\text {rad}}$ typically rank between 3 and 6 , showing moderate influence that varies somewhat by method. $R H, \vartheta, \theta_1$ and $\Delta H_C$ are in bottom ranks, except for $\theta_1$ that is ranked 2 by the Delta method but ranks 8 or 9 by all other methods. This discrepancy is attributable to the non-linear and non-monotonic behavior of $\theta_1$ (Figure 13), which is captured by the moment-independent Delta method but overlooked by linear correlation-based approaches. In fact, PEAR and SRC assume smooth and polynomial inputoutput relationships. In addition, SPEA, SRRC and PRCC capture monotonic trends and miss non-monotonic behavior. However, the Delta method essentially answers, "how much the entire output distribution changes when we fix an input parameter?". It can capture non-linear, non-monotonic and discontinuous effects (the case of $\theta_1$ where the view factor switches at critical angle). Therefore, the Delta method provides complementary information that should not be overlooked. This finding highlights the value of using multiple sensitivity measures.
The sensitivity ranking provides crucial guidance for experimental design by focusing on the most significant parameters. In our case, it indicates that measurement efforts should prioritize accurate characterization of soot properties ( $\phi_{\text {soot}}$ and $\xi$ ) and careful control of wind conditions, as these account for the majority of outputs variability within the specified uncertainty ranges.
This paper addressed uncertainty in estimating pool fire thermal radiation regulatory distances using Monte Carlo sampling. First, the pool fire thermal radiation physical model was introduced, including recent findings related to gasoline. This model accounts for various parameters that contribute to evaluating flame surface emissive power, atmospheric transmissivity, and view factor. Notably, a correction of an existing view factor model was carried out to align the distance reference with other factors (flame surface emissive power and atmospheric transmissivity): all distances now start from the flame front rather than the bund center. Subsequently, uncertainty treatment through Monte Carlo technique was introduced in terms of uncertainty propagation and sensitivity analysis. A gasoline bund fire was then considered to illustrate the challenge of uncertainty consideration in QRA studies.
The variability of input parameters was defined by suitable probability distributions, with their representative parameters reflecting common practice and new findings. A comprehensive Python code was developed to derive a full set of estimates and generate plots: thermal radiation decay versus distance between a hypothetical target and the flame front, distances to threshold thermal radiations (3 kW/m², 5 kW/m² and 8 kW/m²), corresponding percentiles, and several sensitivity indices (classical and advanced).
The results from uncertainty propagation show non-negligible variability for threshold distances. Integrating higher percentile values (e.g., 95th percentile), as opposed to averages, into safety case frameworks and land-use planning constitutes a prudent strategy that not only promotes safer facility design and strengthens operators’ credibility with regulators, but also reduces the risk of encroachment into hazardous areas and provides a more defensible foundation for societal risk assessments. The 95th percentile is proposed here as a precautionary threshold, reflecting a reasonable upper bound commonly adopted in statistical practice. While regulations require the explicit identification and characterization of uncertainty, no precise percentile is universally mandated. We recognize that different land-use or emergency planning contexts may warrant different thresholds. The selected percentile should therefore be scaled according to the stakes and complexity of the regulatory decision, ensuring an appropriate balance between safety and operational feasibility.
The sensitivity study produced parameter importance rankings, indicating that $\emptyset_{\text {soot}}$ and uw are the two most influential parameters. The differing rankings obtained with the Delta method, compared to classical indices for certain parameters, are attributable to its ability to capture non-linear, non-monotonic, and discontinuous effects that correlation-based methods may overlook. This is particularly evident for the threshold behavior in view factor calculations influenced by the angle θ1. This sensitivity ranking could significantly optimize experimental design.
In this study, input parameter variability was exclusively described by probability distributions, which may not always be optimal, particularly for vague data. Alternatives such as fuzzy sets and evidence theory are better suited for such cases. In future work, we plan to extend this methodology to incorporate these alternatives alongside Monte Carlo technique. Moreover, conducting a sensitivity analysis on the choice of input distributions and comparing the effects of alternative assumptions on model outputs (particularly the 95th percentile) represent another promising avenue for future research.
|
$\emptyset(x)$ |
heat flux at a certain distance x (received heat flux) (kW/m2) |
|
$\varnothing_0$ |
flame surface emissive power (kW/m2) |
|
$\tau(x)$ |
atmospheric transmissivity to thermal radiation |
|
$F(x)$ |
view factor |
|
x |
distance between the flame front and the target |
|
r |
distance between the center of the pool and the target |
|
$\emptyset_{\max}$ |
maximum surface emissive power of a flame without soot production (kW/m2) |
|
$\emptyset_{\text {soot}}$ |
surface emissive power of soot (black smoke) (kW/m2) |
|
$\xi$ |
fraction of the flame surface covered by soot. |
|
$D_{e q}$ |
equivalent pool diameter (m) |
|
S |
surface area of the pool (m2) |
|
P |
fire (pool) perimeter (m) |
|
$f_{\text {rad}}$ |
radiation fraction, i.e., the fraction of the combustion energy radiated from the flame surface |
|
$\ddot{m}$ |
mass burning rate per unit area at still weather conditions (kg/(m2.s)) |
|
$m_{\infty}$ |
mass burning rate for a pool of infinite diameter in kg/(m2.s) |
|
$\Delta H_C$ |
heat of combustion of the flammable material at its boiling point (kJ/kg) |
|
$L_F$ |
flame length (m) |
|
k |
absorption extinction coefficient of the flame (m−1) |
|
$\beta$ |
mean beam length corrector |
|
$u_w$ |
wind speed (m/s) |
|
$g$ |
gravitational acceleration (9.80665 m/s2) |
|
u* |
scaled wind velocity |
|
$\rho_{\text {air}}$ |
air density (kg/m3) |
|
$\vartheta$ |
ambient temperature (℃) |
|
RH |
relative humidity (%) |
|
p |
barometric pressure (Pa) |
|
h |
altitude (m) |
|
$p_0$ |
atmospheric pressure at sea level (101325 Pa) |
|
R |
gas constant (8314.462 j/(kmol·K)) |
|
M |
molar mass of Earthly air (28.9644 kg/kmol) |
|
$P_w$ |
partial pressure of water vapor |
|
$P_S$ |
saturation water vapor pressure (Pa) |
|
$L_L$ |
bund length (m) |
|
$L_W$ |
bund width (m) |
|
θ1 |
angle reflecting the location of the target regarding the fire |
[1] Legislation.gov.uk. (2015). The Control of Major Accident Hazards Regulations 2015. http://www.legislation.gov.uk/uksi/2015/483/contents/made.
[2] U.S. Department of Labor, Occupational Safety and Health Administration. (1992). Process safety management of highly hazardous chemicals (Standard No. 1910.119). https://www.osha.gov/laws-regs/regulations/standardnumber/1910/1910.119.
[3] U.S. Environmental Protection Agency. Chemical accident prevention provisions: Risk management program under the Clean Air Act, Section 112(r)(7). https://www.ecfr.gov/current/title-40/chapter-I/subchapter-C/part-68.
[4] European Parliament and Council of the European Union. (2012). Directive 2012/18/EU of the European Parliament and of the Council of 4 July 2012 on the control of major-accident hazards involving dangerous substances, amending and subsequently repealing Council Directive 96/82/EC. Official Journal of the European Union, L 197/1. http://data.europa.eu/eli/dir/2012/18/oj.
[5] Modarres, M., Kaminskiy, M.P., Krivtsov, V. (2016). Reliability Engineering and Risk Analysis: A Practical Guide. CRC Press.
[6] Zio, E. (2007). An Introduction to the Basics of Reliability and Risk Analysis. World Scientific.
[7] Environmental Protection Agency. (2026). Areal Locations of Hazardous Atmospheres (ALOHA) (Version 5.4.7). https://www.epa.gov/cameo/aloha-software.
[8] DNV. (2023). Process Hazard Analysis Software Tool (PHAST) (Version 8.6). https://www.dnv.com/services/phast-software-for process-hazard-analysis.
[9] Gexcon. (2023). EFFECTS (Version 10.3). https://www.gexcon.com/products-services/effects/.
[10] Ferdous, R., Khan, F., Sadiq, R., Amyotte, P., Veitch, B. (2013). Analyzing system safety and risks under uncertainty using a bow-tie diagram: An innovative approach. Process Safety and Environmental Protection, 91(1-2): 1-18. https://doi.org/10.1016/j.psep.2011.08.010
[11] Abrahamsson, M. (2002). Uncertainty in Quantitative Risk Analysis: Characterisation and Methods of Treatment. Licentiate Thesis, Lund University Publications.
[12] U.S. Nuclear Regulatory Commission. (2009). Guidance on the treatment of uncertainties associated with PRAs in risk-informed decision making (NUREG-1855). U.S. Nuclear Regulatory Commission, Washington, DC.
[13] Chettouh, S., Hamzi, R., Innal, F., Haddad, D. (2014). Industrial fire simulation and uncertainty associated with the emission dispersion model. Clean Technologies and Environmental Policy, 16(7): 1265-1273. https://doi.org/10.1007/s10098-014-0792-x
[14] Siuta, D., Markowski, A.S., Mannan, M.S. (2013). Uncertainty techniques in liquefied natural gas (LNG) dispersion calculations. Journal of Loss Prevention in the Process Industries, 26(3): 418-426. https://doi.org/10.1016/j.jlp.2012.07.020
[15] Chaudhary, R.K., Van Coile, R., Gernay, T. (2021). Potential of surrogate modelling for probabilistic fire analysis of structures: RK chaudhary. Fire Technology, 57(6): 3151-3177. https://doi.org/10.1007/s10694-021-01126-w
[16] Xie, Q., Lu, S., Kong, D., Wang, J. (2013). Treatment of evacuation time uncertainty using polynomial chaos expansion. Journal of Fire Protection Engineering, 23(1): 31-49. https://doi.org/10.1177/1042391512470578
[17] Kurzawski, A., Cabrera, J.M., Ezekoye, O.A. (2020). Model considerations for fire scene reconstruction using a Bayesian framework. Fire Technology, 56(2): 445-467. https://doi.org/10.1007/s10694-019-00886-w
[18] Zhang, Y., Chen, L., Shi, Y., Wang, H. (2025). Dynamic modeling and real-time risk assessment of tunnel fires using pressure-state-response (PSR)-enhanced Bayesian network. Results in Engineering, 28: 107716. https://doi.org/10.1016/j.rineng.2025.107716
[19] Zadeh, L.A. (1978). Fuzzy sets as a basis for a theory of possibility. Fuzzy Sets and Systems, 1(1): 3-28. https://doi.org/10.1016/0165-0114(78)90029-5
[20] Shafer, G. (1976). A Mathematical Theory of Evidence. Princeton University Press.
[21] Cadena, J.E., Osorio, A.F., Torero, J.L., Reniers, G., Lange, D. (2020). Uncertainty-based decision-making in fire safety: Analyzing the alternatives. Journal of Loss Prevention in the Process Industries, 68: 104288. https://doi.org/10.1016/j.jlp.2020.104288
[22] Ministry of Housing, Spatial Planning and the Environment. (2005). Methods for the Calculation of Physical Effects: Due to Releases of Hazardous Materials (Liquids And Gases). Third Edition. https://publications.tno.nl/publication/34634119/QIKv78/TNO-2005-yellow.pdf.
[23] Hägglund, B., Persson, L.E. (1976). The heat radiation from petroleum fires. FOA Rapport C 20126-D6(A3), National Defence Research Institute, Stockholm, Sweden.
[24] Ingason, H., Lönnermark, A. (2011). Fire spread between industrial premises. In Fire Safety Science-Proceedings of the Tenth International Symposium, pp. 1305-1318. https://doi.org/10.3801/IAFSS.FSS.10-1305
[25] Muñoz, M., Planas, E., Ferrero, F., Casal, J. (2007). Predicting the emissive power of hydrocarbon pool fires. Journal of Hazardous Materials, 144(3): 725-729. https://doi.org/10.1016/j.jhazmat.2007.01.121
[26] Moorhouse, J., Pritchard, M.J. (1982). Thermal radiation hazards from large pool fires and fireballs: A literature review. The Institution of Chemical Engineers Symposium Series, 71: 397-428. https://www.icheme.org/media/25191/the-assessment-of-major-hazards-icheme-symposium-series-71-1982-25-moorhouse.pdf#4#1.
[27] Duiser, J.A. (1989). Warmteuitstraling (radiation of heat). Methods for the calculation of the physical effects of the escape of dangerous materials (liquids and gases). Report of the Committee for the Prevention of Disasters, 2nd ed. Ministry of Social Affaire, The Netherlands.
[28] Center for Chemical Process Safety. (2025). Risk analysis screening tool (RAST) and chemical hazard engineering fundamentals (CHEF). American Institute of Chemical Engineers. https://www.aiche.org/ccps/resources/tools/risk-analysis-screening-tool-rast-and-chemical-hazard-engineering-fundamentals-chef.
[29] Roberts, A.F. (1981). Thermal radiation hazards from releases of LPG from pressurised storage. Fire Safety Journal, 4(3): 197-212. https://doi.org/10.1016/0379-7112(81)90018-7
[30] Yang, J.C., Hamins, A., Kashiwagi, T. (1994). Estimate of the effect of scale on radiative heat loss fraction and combustion efficiency. Combustion Science and Technology, 96(1-3): 183-188. https://doi.org/10.1080/00102209408935354
[31] Koseki, H. (1989). Combustion properties of large liquid pool fires. Fire Technology, 25(3): 241-255. https://doi.org/10.1007/BF01039781
[32] Babrauskas, V. (1983). Estimating large pool fire burning rates. Fire Technology, 19(4): 251-261. https://doi.org/10.1007/BF02380810
[33] Thomas, P.H. (1963). The size of flames from natural fires. Symposium (International) on Combustion, 9(1): 844-859. https://doi.org/10.1016/S0082-0784(63)80091-0
[34] Thomas, P.H. (1965). Fire spread in wooden cribs: Part III the effect of wind. ire Research Note No. 600, Fire Research Station, Borehamwood, UK. https://publications.iafss.org/publications/frn/600/-1/view/frn_600.pdf.
[35] Lente, G., Ősz, K. (2020). Barometric formulas: Various derivations and comparisons to environmentally relevant observations. ChemTexts, 6(2): 13. https://doi.org/10.1007/s40828-020-0111-6
[36] Bagster, D.F., Pitblado, R.M. (1989). Thermal hazards in the process industry. Chemical Engineering Progress, 85(7): 69-75.
[37] International Organization for Standardization. (2012). Hygrothermal performance of building components and building elements-Internal surface temperature to avoid critical surface humidity and interstitial condensation-Calculation methods (ISO Standard No. 13788). https://www.iso.org/standard/51615.html.
[38] Rew, P.J., Hulbert, W.G., Deaves, D.M. (1997). Modelling of thermal radiation from external hydrocarbon pool fires. Process Safety and Environmental Protection, 75(2): 81-89. https://doi.org/10.1205/095758297528841
[39] United States Environmental Protection Agency. (2009). Guidance on the development, evaluation, and application of environmental models. https://www.epa.gov/cea/guidance-development-evaluation-and-application-environmental-models.
[40] National Aeronautics and Space Administration. (2002). Probabilistic risk assessment procedures guide for NASA managers and practitioners. NASA Office of Safety and Mission Assurance. https://ntrs.nasa.gov/api/citations/20120001369/downloads/20120001369.pdf.
[41] Ferson, S., Ginzburg, L.R. (1996). Different methods are needed to propagate ignorance and variability. Reliability Engineering & System Safety, 54(2-3): 133-144. https://doi.org/10.1016/S0951-8320(96)00071-3
[42] Zio, E., Pedroni, N. (2013). Methods for representing uncertainty: A literature review. Cahiers de la Sécurité Industrielle, No. 2013-03, Fondation pour une Culture de Sécurité Industrielle, Toulouse, France. https://doi.org/10.57071/124ure
[43] Saltelli, A., Chan, K., Scott, E.M. (2000). Sensitivity Analysis. John Wiley & Sons.
[44] Meteoblue. (2026). Climate simulation Skikda. https://www.meteoblue.com/fr/meteo/historyclimate/climatemodelled/skikda_alg%c3%a9rie_2479536.
[45] ONM. (2026). Weather records from local Skikda station. Algerian National Meteorological Office.
[46] Hersbach, H., Bell, B., Berrisford, P., et al. (2020). The ERA5 global reanalysis. Quarterly Journal of the Royal Meteorological Society, 146(730): 1999-2049. https://doi.org/10.1002/qj.3803
[47] Amira, A., Innal, F. (2025). Risk assessment of industrial domino effects using stochastic Petri nets. International Journal of Safety and Security Engineering, 15(8): 1691-1701. https://doi.org/10.18280/ijsse.150814
[48] Borgonovo, E. (2007). A new uncertainty importance measure. Reliability Engineering & System Safety, 92(6): 771-784. https://doi.org/10.1016/j.ress.2006.04.015
[49] Plischke, E., Borgonovo, E., Smith, C.L. (2013). Global sensitivity measures from given data. European Journal of Operational Research, 226(3): 536-550. https://doi.org/10.1016/j.ejor.2012.11.047
