Dynamic Simulation on Deflagration of LNG Spill

The deflagration characteristics of premixed LNG vapour-air mixtures with different mixing ratios were quantitatively and qualitatively investigated by using CFD (computational fluid dynamics) method. The CFD model was initially established based on theoretical analysis and then validated by a lab-scale deflagration experiment. The flame propagation behaviour, pressuretime history, and flame speed were compared with the experimental data, upon which a good agreement was achieved. A largescale deflagration fire during LNG bunkering process was conducted using the model to investigate the flame development and overpressure effects. Mesh independence and time scale were tested in order to obtain the suitable grid resolution and time step. Deflagration cases with two different LNG vapour volume fractions, i.e., 10.4% and 15.0%, were simulated and compared. The one with a volume fraction of 10.4% which was around stoichiometric mixing ratio had the highest flame propagating speed. High flame velocity observed in the simulation was coupled with the thin flame front where overpressure occurred. The CFD model could capture the main features of deflagration combustion and account for LNG fire hazard which could provide an in-depth insight when dealing with complicated cases.


Introduction
Due to the fast growing world economy and increasing globalization, the shipping industry is of great importance to the worldwide economic growth.Environmental concerns have forced the industry to reduce the emission of greenhouse and polluting gases, and to investigate the alternative fuels for sustainable shipping.It has been obliged by EU Bunkering Regulation [1] to require the shipping network to facilitate LNG fuelling by 2025.However, the safety issues with respect to the new LNG bunkering facilities concern the general public and authorities importantly.The hazards related to the accidental release during LNG transporting and bunkering include the cryogenic vapour dispersion, fire radiation, and explosion [2].LNG fire could happen when the vapour fire flashes back to the LNG impoundment (known as LNG pool fire) or the disconnected vapour cloud is ignited (known as deflagration).The industrial standard NFPA 59A [3] detailed the safety aspects with regard to LNG vapour dispersion and fire radiation.Increasing attention was obtained to the global discussion on assessing the safety exclusion zone during shipto-ship LNG bunkering.However, the universal consensus was not yet reached due to the lack of analysis of hazardous consequence [4].In this study, the hazards of LNG vapour deflagration during ship-to-ship bunkering were investigated.This paper aimed at studying the flame propagation and the overpressure effect in the process of deflagration by conducting CFD (computational fluid dynamics) simulation.It could provide an in-depth understanding of conducting risk assessment for LNG bunkering process.
In the previous studies, substantial research work on methane/natural gas deflagration, both experimentally and numerically, was performed to deal with flame-obstacles interactions.Lind [5] was among the earliest researchers to perform the larger-field-scale LNG spill test by assuming that a collision occurs between the tanker and another ship, leading to a maximum boil-off gas 26000 kg/sec in which natural gas combustion properties were investigated.Large-scale filed tests, such as Maplin Sands Experiments [6] and Coyote series tests [7], were performed to monitor detailed LNG combustion behaviour.Harrison and Eyre [8] conducted a medium-scale LNG vapour combustion experiment in which 4000 m 3 of premixed natural gas/air was ignited in a wedge-shaped and obstructed enclosure in order to monitor the flame speed and overpressure.However, largescale experiments with built-in solid obstacles are very hard to explore the deflagration mechanism due to the difficulty of conducting detailed measurements [9,10].Using advantaged test technology of laser diagnostic, some lab-scale experiments [11][12][13][14][15][16] were performed to investigate the mechanism of deflagration in complex geometries and short periods.By performing the deflagration test [13] in a vented chamber with solid obstacles varying in numbers and locations, it was concluded that obstacles had significant impacts on the overpressure effect and flame structure.Deflagration experiment on the flame-obstacle interactions was conducted [14] by considering three different obstacle configurations, i.e., side, central, and staggered, respectively.A comprehensive experimental study of natural gas combustion in a stainless steel cylinder vessel (diameter 180 mm and height 210 mm) was conducted by Tang, Zhang et al. [17] to investigate the influence of the initial conditions such as equivalence ratios, pressures, and temperatures, as well as monitoring the peak overpressures to provide fundamental data for natural gas engine timing control.Ajrash, Zanganeh et al. [16] conducted a flame deflagration test of premixed methane-air mixture in a large-scale detonation tube with a diameter of 0.5 m and a total length of 30 m.It was to investigate the length of reactive sections (RS) and relations of flame and pressure wave at different methane concentrations.
Because of major improvement of computational power, numerical research work on premixed combustion was developed dramatically in the past decades.RANS-based (Reynolds Averaged Navier-Stokes) turbulence models (e.g.,  −  model) were widely applied [12,[18][19][20].Due to the fast development of computing capability and requirements for more accurate and precise prediction in combustion, LES (Large Eddy Simulation) model attracted an increasing interest and was employed in many different combustion problems [21].A great number of studies on modelling premixed deflagration and flame propagation demonstrated the validity and accuracy of LES model [22][23][24][25][26][27][28].Compared with the classical RANS approach, LES model could provide an improved explanation of turbulence and flame-obstacle interaction in premixed flame deflagration.Besides the turbulence model, submodel of premixed combustion was an important input to calculate reaction rate and flame speed.Gavelli, Davis et al. [29] used the software FLACS to investigate the consequences of vapour cloud explosions (VCEs) during the process of LNG carrier offloading.Wen, Yu et al. [26] discussed combustion models in different subgrid scales.Xu, Cong et al. [30] applied Zimont model [31] as turbulence flame speed model and modified the flame speed constant in order to better explain deflagration mechanism.
This study aimed at characterising the combustion features in deflagration.The mechanism of deflagration and laminar flame speed was analysed.The LES model of methane-air deflagration in an obstructed domain was established in order to investigate the hazards of LNG vapour fire during bunkering process.Due to the lack of experimental data in large-scale methane deflagration, the model was validated by comparing the simulation with a lab-scale experiment.The model was then used to perform the analysis of large-scale LNG vapour deflagration in which flame propagation under different equivalence ratios and overpressure effects were investigated.

Model Theory Basis
When ignition occurs in a stagnant flammable natural gas cloud, the flame starts to propagate away from the ignition point through a very low-speed laminar flow.As the combustion products expand and flame propagates, the flame speed is increased dramatically in a very short duration.The overpressure in vapour cloud combustion is directly coupled with the speed at which the flame front runs through the cloud.The faster the flame propagates through the flammable vapour cloud, the higher the overpressure will be, which could enhance the blast effects outside the cloud.It implies that the mode of flame propagation is very important, upon which there exist two mechanisms, i.e., detonation and deflagration.In detonative combustion, the flame front is propagated by a shock wave which compresses the mixture beyond its autoignition temperature.At the same time, the detonation wave with supersonic speed is still maintained by the heat released from the combustion and thus generates super high overpressures (1.5∼2.0MPa).For a detonation to propagate, experimental research suggests that the flammable cloud must be rather mixed homogeneously, and the ignition source of a high-explosive charge is required for initiating a detonation.Therefore, the likelihood of a detonation in the open air to occur is quite low, while the deflagration is the most common combustion mode.A deflagration, in general, is the mode of flame propagation which is determined largely by heat conduction and molecular diffusion of heat and species.Heat is generated by chemical reaction in combustion and transported ahead of the reaction zone into an unburnt zone for preheating.Since molecular diffusion is a relatively slow process, laminar flame development is slow.The maximum laminar burning velocity for methane combustion is 0.448m/s.Generally, the deflagration flame speed between very low (a few meters per second) and very high (more than 100 meters per second) is due to hot combustion products rapid expansion creating fast flow in the flame front.Therefore, the overpressure is ranged from a few to thousands pascals.
In order to investigate deflagration both qualitatively and quantitatively, commercial CFD code ANSYS FLUENT 16.0 was used to develop the deflagration model in this study.It utilized the Finite Volume Method (FVM) [32] to discretise the computational domain and equations.Partial differential equations (PDE) of continuity, energy, three-dimensional momentum, and turbulence [2,14] were involved.Besides, submodels as premixed combustion model were taken into consideration when modelling deflagration.
Turbulence in stagnant vapour cloud deflagration can be generated mainly due to combustion products expanding into unburnt zones and changing the flow mechanism.Significant research has shown that turbulence in deflagration process could dramatically accelerate the combustion rate.In certain extreme cases, the turbulence can cause the flame propagation mode to change suddenly from deflagration into detonation, namely, a Deflagration to Detonation Transition (DDT).At the macroscopic level of combustion, turbulent mixing blends the hot combustion products and unburnt mixture which can be characterized by the large-scale motion of turbulence accounting for most of the kinetic energy, while, at the microscopic level, turbulent combustion is a molecular mixing process of chemical reaction and is of great importance in chemical conversion.Only molecular diffusion and mixing can maintain the occurrence of the chemical reaction.Molecular mixing of scalar quantities occurs essentially on the smallest turbulent scales which plays a significant role in deflagration modelling.This can be characterized by small-scale eddies which are responsible for the dissipation of turbulence scalars [33].Turbulence model for deflagration simulation in this study is selected as Large Eddy Simulation (LES) model rather than RANS.LES model characterizes turbulent flows by eddies with a wide range of length and time scales and separates the turbulent fields into large-scale resolved and small-scale unresolved contributions [21,34].Spatial filtering function used in instantaneous turbulent fields removes turbulent eddies of scales smaller than the filter size.Large eddies are resolved directly to mostly account for transporting momentum, mass, energy, and other passive scalars, while small eddies are modelled to improve accuracy and to be consequently more universal.It has been already certified that LES provided substantial advantages for modelling turbulent combustion.Compared with RANS, LES predicts the scalar mixing process and dissipation rates with considerably improved accuracy in turbulent premixed or nonpremixed combustion [21,[33][34][35][36][37].The simple LES model was first proposed by Smagorinsky [38], in which the eddy viscosity is modelled by where   is turbulent viscosity;  is fluid density;  is the average rate-of-strain tensor; and   is mixing length, expressed by where  is the von Karman constant; d is the distance to the closest wall;   is the Smagorinsky constant; and Δ is the local grid scale.
This study aimed at investigating the flame propagation and overpressure effect where the effect of turbulence is dominant.It can be very challenging to model the detailed chemical reactions in deflagration, which is due to the nonlinear relation with chemical and thermodynamic states [39,40].The deflagration flames are often characterized by propagating thin reaction flame front which is much smaller than the turbulent scales and is strongly affected by turbulence (causing flame wrinkling and complex thermochemical-turbulence interactions).The complexity can be reduced by assuming the single-step irreversible chemistry and Zeldovich instability (thermal diffusion) [24].Thus, the overall flame propagation is characterised by turbulent eddies, in which the turbulent flame speed is mainly influenced by both laminar flame speed and turbulence length.A scalar variable characterising the reaction progress from burnt reactants to unburnt mixtures, c, is defined as a normalized mass fraction of products.c equals zero in an unburnt mixture, one in the burnt products, and varies between 0∼1 in the flame front.This method has been used widely to investigate the flame propagation and overpressure effect of deflagration [14,24,28].
The progress variable is evaluated by a transport equation expressed as where  is mean reaction progress variable;   is turbulent Schmidt number; and   is reaction progress source term and is expressed as where   is the density of unburnt mixture and   is turbulent flame speed.Zimont model [31] is applied in this study to compute turbulent flame speed, given by where  is model constant (0.52 for most premixed flames);   is laminar flame speed;   is root-mean-square velocity;  is thermal diffusivity of the unburnt mixture; and   is turbulence length which is modelled in LES, given by where   is the Smagorinsky constant and û is the cell characteristic length.Laminar flame speed (  ) is the key parameter in many areas of combustion science and plays an essential role in determining several important aspects of the combustion process in spark ignition.The laminar flame velocity of methane/air has attracted the most attention and has already been investigated both experimentally and numerically in numerous studies.Gülder correlation [41,42] of laminar flame speed is widely acknowledged and thus is used in this study, given by where  is the equivalence ratio, W=0.422m/s, and other parameters related to natural gas combustion , , and  are equal to 0.15, 5.18, and 1.075, respectively.Gas combustion usually has the highest laminar burning velocity in the vicinity of stoichiometric mixing ratio.The relations of laminar flame speed and the equivalence ratio of natural gas combustion are shown in Figure 1.
The simulation was performed by using the commercial code ANSYS FLUENT 16.0.The SIMPLE algorithm was employed to account for pressure-velocity coupling.Parallel computations were performed in this study.In the computational settings, approximately 20 iterations per time step were required for residual values reducing by three orders of magnitude and below an acceptable value (i.e., 1.0 × 10-5).Simulations were executed in parallel by using 144 cores on iVEC (a government-supported high-performance computing national facility located in Western Australia).

Model Validation
Deflagration experiment in a semiconfined chamber conducted by Patel, Jarvis et al. [12] was referred to, to validate the model in this study.As illustrated in Figure 2, the experiment setup was composed by a polycarbonate chamber with a dimension of 150 mm (L) x 150 mm (W) x 500 mm (H).Three rectangular obstacles of block ratio 50% (75 mm long, 150 mm wide, and 10 mm height) were located vertically at a spacing distance of 100 mm.A high-speed laser-sheet flow visualization system (9000 Hz) was employed to image and record the flame propagation.The reactants were selected as the stoichiometric methane-air mixture which was purged through the chamber and ignited at the centre of the bottom plate.Combustion was triggered by a hemispherical spark (radius 5 mm) and lasted for 5 ms.When entering the chamber, the premixed reactants were seeded with micrometresized droplets of olive oil in order to scatter laser light and track deflagration flame.Before ignition, the top of the chamber was sealed with a thin PVC membrane which was allowed to rupture and vent the pressure in deflagration process.Overpressure was monitored in the vicinity of the ignition point by using piezoelectric pressure transducers (range 0∼1 bar, response time 0.1 ms).The computational domain was established based on the experimental setup.Nearby the obstacle area, the mesh was refined in order to observe the flame development.The overall number of hexahedral cells was 1.17 million.The numerical results, including flame propagation and pressuretime history, were compared with the experimental data.was represented by dark grey and red areas in experimental and simulation images, respectively.The deflagration flame in both experiment and simulation results was propagated and fully developed within a very short period (approximately 40 ms).In the initial 20 ms after ignition, it was found that the flame propagated very slowly at almost a laminar flame speed.Once the flame reaches the first obstacle, the deflagration process was speeded up immediately with flame jetting through the second and third obstacles in 6 ms which implied the obstacle in deflagration could enhance the turbulence and thus propagate the flame dramatically.It was observed that the flame propagation process was well predicted by the CFD model.A good agreement between simulation and experiment was achieved.
Figure 4 showed the comparison between the simulation results and the measured values of overpressures and vertical flame speed.The pressure rise started from around 28 ms when the flame passed the first obstacle and then increased dramatically to thousands of pascals.The peak overpressure was approximately 12 kPa and occurred at around 35 ms when the flame just vented through the obstacles.It could be attributed to the blockage effect by presenting the three obstacles which resulted in flame jetting and producing high back pressure at the bottom of the chamber, as well as the sealing PVC membrane used to enclose chamber that could possibly build up pressure before membrane rupture.Once the flame propagates through the chamber (at approx.38 ms), the pressure was decreased significantly and was even below zero which could be due to the flame bursting out and creating a vacuum in the other chamber end.The flame jetting effect could also be observed in the flame speed diagram (Figure 4).Whenever the flame passed an obstacle, the flame speed reached a peak value.As the obstacles create blockage effect and strengthen turbulence significantly, the flame speed was increased with time and reached the maximum value up to 65 m/s after passing the third obstacle.

Designed Scenario
Through comparing with the experimental results, the deflagration model obtained a good agreement and was proposed to evaluate the LNG combustion hazard during bunkering process.Generally, two distinct types of LNG discharge events were envisaged [43] in history.Firstly, LNG is released due to accidental collision, grounding, allision, or intentional attack damaging the vessel structure and puncturing the LNG tank, which could lead to a large discharge.However, much design and operational effort have been highlighted to assure the rarity of significant leak events.Secondly, small LNG discharge could occur because of pipework, valve, or loading arm failure during transfer operations, which allows up to a full pumping rate discharge.In this study, the LNG pipe or loading arm was assumed to be ruptured and disconnected because of significant wave or unintentional ship movement during bunkering process.This could lead to a significant LNG spill in the water area.LNG evaporation could be speeded up due to the contact with open water or water curtain (used to protect hull during bunkering).The LNG vapour could be entrained by air flow and dispersed in the downwind direction.If sparks or any other unwanted combustion sources present in the area, LNG fire is more likely to take place.The authors performed and completed the consequence dynamic analysis of LNG vapour dispersion and LNG pool fire in previous studies [2,44,45].If the vapour cloud reaches the flammable limits (i.e., 5.0∼15.0%for natural gas), both deflagration and detonation could possibly occur.However, a detonation which is in a time scale of millisecond to microsecond usually requires high degrees of confinement, strong mixing with air, and large ignition sources.Deflagration is the most common combustion mode to occur when LNG vapour is sufficiently mixed with air in molecular level and combusted rapidly (i.e., in the order of seconds) but is still hazardous and lethal due to generation of high momentum impact force and high burning temperature.This research work was mainly focused on the dynamic simulation of LNG deflagration.The assumption was made that the bunkering area (including LNG bunker and client ship) was covered by the flammable mixture with a uniform gas volume fraction after LNG vapour dispersing for a while.Ignition source or spark was assumed to present at the end of the ship (i.e., sparks from engines) and trigger the combustion.The computational domain (90 m long, 70 m wide, and 30 m high) and mesh were illustrated in Figure 5 The mesh independence and time-step independence were tested in order to obtain a proper grid resolution and time scale.In Table 1, the number of cells and average wallclock time were listed.

Verification of Mesh and Time-
Step Independence.Four different grid resolutions were tested, i.e., 0.5 million, 2.0 million, 4.0 million, and 17.0 million, upon which burning area and flame propagation were compared.As illustrated in Figure 6, the burning area at the bottom plane (XY plane) in cases of different grid resolutions was analysed.In general, the combustion was propagated very slowly in the initial 1.0 s, while the flame was developed rapidly from 1.0 s to 2.5 s.After 3.0 s, the burning area reached its maximum value.However, compared with the other cases, the coarse grid had a relatively larger burning area and could possibly overestimate the flame propagation, while the cases of the intermediate grid and fine grid gave comparable results.
The effect of grid resolution could also be compared from flame development and flame thickness, as illustrated in Figure 7.As the mesh resolution was increased, the flame surface area was refined.The flame could be qualitatively represented by progress variable, c, which ranged from zero  (unburnt reactants) to one (burnt products, unity behind the flame front).As shown in Figure 7, the flame interface became sharper as the flame thickness was decreased in relatively higher grid resolution.The simulation results could qualitatively agree with the practical phenomenon where deflagration is defined as a thin flame surface propagating with time.However, as the finer grid was relatively more time costing and computationally expensive, the intermediate grid (4.0 million) was applied in the following study.

Flame Propagation in Different Vapour
Fractions.Once the volume fraction of LNG vapour cloud reaches its flammable limit (i.e., 5%∼15%), flash fire may occur when presenting some ignition sources.In the present study, two different flash fire cases with different gas volume fractions as 10.4% (equivalence ratio 1.1) and 15% (equivalence ratio 1.7) were compared.As demonstrated in Figure 1, mixture equivalence ratio of 1.1 had the highest laminar flame speed (i.e., 0.427m/s), while the speed was 0.082m/s at the ratio of 1.7.Figures 9 and 10 illustrated the simulation results of deflagration process at different mixing ratios.It was observed that the combustion was started initially at the ignition points (assuming the ignition started at one end of LNG bunker) and then propagated to a further distance.The flame was fully developed within 3.0 s in the case of equivalence ratio 1.1, while the combustion was significantly postponed and prolonged in the mixture equivalence ratio of 1.7.The duration of combustion covering the whole domain was approximately 9.0 s which was three times longer than that in the equivalence ratio of 1.1.It was found that the combustion was postponed and more incomplete when equivalence ratio was higher.The same combustion phenomenon had been observed in the experiment of Tang, Zhang et al. [17] in which the cases of equivalence ratio 1.0∼1.4 were studied.

Pressure Dynamics in Deflagration.
Because a large portion of combustion energy (methane combustion energy 55.5MJ/kg) is released during deflagration, the overpressure is more likely to take place by developing blast wave uniformly in all directions.In contrast with the temperature contour, the maximum overpressure occurred at the flame front (thin flame sheet), as illustrated in Figure 11, with maximum overpressure approximately 300 Pa, while the pressure was relatively lower inside the combusted area.Theoretically, the vapour cloud from LNG vapour dispersion could possibly cause vapour cloud fire (deflagration) or vapour cloud explosion (detonation).However, the detonation usually involves a confined space and higher minimum ignitions energy (2.3e+08J,only 0.28mJ for deflagration) and thus generates high flame propagation speed (2∼5 times of sound speed) and overpressure (approx.1.5∼2.0MPa)[46].Because the vapour cloud was assumed to be presented in the open air, it was observed that the vapour deflagration took place with minor overpressure effect.Figure 12 showed the overpressure effect and the interaction between pressure dynamics and flame propagation by comparing the pressure field with flame surface and velocity field.It was found that the maximum velocity and pressure which is attributed to combusting/exhausted gas moving were coupled at the outer edge of flame front, with the velocity direction normal to the flame surface.

Conclusion
A dynamic model of natural-gas deflagration was developed in order to analyse the hazard consequence during ship-toship LNG bunkering.The LES turbulence model was applied to take account of small-scale eddies and to improve the accuracy.In lieu of large-scale experimental data, the CFD model of deflagration was validated using the data from a lab-scale experiment upon which pressure history, burning speed, and flame propagation were compared.The model was then applied to investigate vapour cloud fire during LNG bunkering.Two different equivalence ratios of the methane-air mixture were investigated both quantitatively and qualitatively.It was observed in the simulation that the flame developed at a very low speed initially and then propagated dramatically in a very short duration.The CFD Velocity (m/s) 0 30 Dynamic pressure (Pa) 0 300 model could reflect the facts that the most rapid flame propagation occurred when the mixture had an approximate stoichiometric mixing which produced the highest laminar flame speed, while the combustion with equivalence ratios higher than stoichiometric mixture was less robust.It was also observed in the simulation that the overpressure effect occurred at the thin flame front which was because highest velocity presented in the flame front led to a high pressure.
The findings in this study could match the understandings of the deflagration and show the validity of the CFD model in large-scale cases.It can be applied in different case scenarios and is of critical importance in the process of risk assessment.

Figure 3 Figure 3 :Figure 4 :
Figure 3: Comparison between experimental and simulated flame propagation at different times after ignition.

Figure 8
Figure 8 illustrated the graphical comparison of different time steps for calculating burning area.It was found that the time steps of 10 −3 s and 10 −4 s had comparable values, and thus the time step 0.001 s was selected in further case studies in order to save computing time.

Table 1 :
Number of cells and average wall-clock time per iteration of different grid resolutions.