ARTICLE IN PRESS
Journal of Loss Prevention in the Process Industries 20 (2007) 455–461 www.elsevier.com/locate/jlp
Numerical analysis of gas explosion inside two rooms connected by ducts Akinori Hashimotoa, Akiko Matsuob, a
Graduate School of Science and Technology, Keio University, 3-14-1, Hiyoshi, Kohoku-ku, Yokohama, Kanagawa 223-8522, Japan b Department of Mechanical Engineering, Keio University, 3-14-1 Hiyoshi, Kohoku-ku, Yokohama, Kanagawa 223-8522, Japan Received 1 December 2006; received in revised form 5 April 2007; accepted 5 April 2007
Abstract Simulations of gas explosion of hydrogen/air mixture inside two rooms connected by ducts are carried out. Scalar transport chemical reaction model and LES turbulence model are utilized to reduce the calculation load and to conduct real-scale analysis. The effects of ignition source locations and volume of ignited room are analyzed, and the time history of pressure and rate of pressure rise in each room are focused in this study. When the volume of the ignited room is larger than the other room, the high pressure from the other room causes a force to act on the partition to the ignited room. This study indicates that the current technique can predict specific features of gas explosions inside two rooms connected by the ducts. r 2007 Elsevier Ltd. All rights reserved. Keywords: Gas explosion; Flame propagation; Rooms connected by ducts; Differential pressure
1. Introduction Gas explosion is a frequently occurring accident and as many as 200 gas explosion accidents are reported throughout the year in Japan. Especially in the closed space such as chemical facility and atomic power plant, devastating gas explosion is more likely to happen. The damage in such case is often more severe than that of the open space since energy is not released outside the building. Furthermore, complicated building configurations and existence of windows, doors and ducts often cause the increase in damage since these dramatically change the gas flow field and cause strong turbulence. For these reasons, determining the effects of such configurations on gas explosion is required. Numerical analysis is one of the effective ways of predicting the damage of gas explosion. Miura et al. (2004) performed the numerical analysis that simulated the blast wave propagation of gas explosion accident in an
Corresponding author. Tel.: +81 45 566 1518; fax: +81 45 566 1495.
E-mail address:
[email protected] (A. Matsuo). 0950-4230/$ - see front matter r 2007 Elsevier Ltd. All rights reserved. doi:10.1016/j.jlp.2007.04.009
atomic power plant. They assumed that gas explosion is caused by rapid energy release of high explosive, the energy supply from chemical reaction of combustible gas does not occur after ignition. Firweather, Hargrave, Ibrahim, and Walker (1999) performed the numerical analysis of gas explosion in a small tube including some orifices. They utilize progress variable model as chemical reaction model, and consider the energy release of combustion gas in the tube. It is important when predicting damage of gas explosion using numerical analysis. In usual building and plants, ducts or pipelines connect many rooms. Holbrow, Lune, and Tyldesley (1999) performed the experiments of dust explosion in connected vessels, and proposed the guidance for containment and venting. The purpose of this research is to clarify the effect of the configuration in connected rooms on gas explosion considering the energy release of combustion gas. Therefore, two rooms connected by ducts are adopted as calculation model, and the volume of the ignited room and the ignition source location are set as varying parameters for this study. In the study, scalar transport chemical reaction model (Tominaga, Ito, Taniguchi, & Kobayashi, 1999) and LES turbulence model are utilized.
ARTICLE IN PRESS A. Hashimoto, A. Matsuo / Journal of Loss Prevention in the Process Industries 20 (2007) 455–461
456
2. Calculation method In this study, the structure of the flame front is not considered and it is regarded as the discontinuous boundary plane of burnt and unburned mixtures. The chemical reaction is assumed to occur at the boundary between these mixtures, and the propagation of the flame front is calculated by solving the equation of progress variable. To consider the change of gas concentration, the equation of fuel mixture fraction is solved. The turbulent effect is estimated by LES turbulence model, which mixed scale model (Lenormand & Sagaut, 2000) is employed for sub-grid scale turbulence evaluation. 2.1. Governing equations The governing equations are the equations of mass, momentum, energy, progress variable, and fuel mixture fraction for compressible fluid. Propagation of flame is modeled by progress variable equation. To account for the combustion of non-uniform mixture, equation of fuel mixture fraction is calculated. Equation of progress variable employed for LES is, qðrcÞ qðrcuj Þ q mSGS qc _c þ ¼ þo qt qxj ScSGS qxj qxj _ c ¼ rub St jrc j, o
ð1Þ
where c is the progress variable which is a normalized burnt mass fraction. In unburned mixture region, c ¼ 0, while c ¼ 1 in burnt mixture region, and 0oco1 in the region between unburned and burnt mixture. The symbol r is the density, ru is the unburned mixture density, uj is the velocity component along xi coordinate, mSGS is the viscosity in sub-grid scale, and ScSGS is the Schmidt number in sub-grid scale which takes a constant value of 0.5 in this study. Turbulent burning velocity is represented by St, which is determined by unburned mixture fraction, temperature and pressure. Equation of fuel mixture fraction employed for LES is the following: qðrxÞ qðrxuj Þ q m mSGS qx þ ¼ þ . (2) qt qxj Scx ScSGS qxj qxj Here Scx is the Schmitd number for the equation of the fuel mixture fraction. 2.2. Burning velocity Laminar burning velocity is obtained by the database of laminar burning velocity calculated using commercial software CHEMKIN. Reaction model of 10 species and 21 elementary reations is employed in the CHEMKIN calculation. The dependence of laminar burning velocity at standard atmosphere conditions on the equivalence ratio calculated by CHEMKIN and experimental result (Iijima & Takeno, 1986) are shown in Fig. 1. Both results show
Fig. 1. Dependence of laminar burning velocity on equivalence ratio.
maximum values at an equivalence ratio of 1.8 and tend to decrease as they reach to flammability limits (4.0% v/v and 75.0 v/v%). To consider the dependence of burning velocity on unburned mixture’s temperature and pressure, the model reported by Iijima and Takeno (1986) is employed in this study, and is shown below: T ub b1 P Sl ¼ Sl 0 1 þ b2 ln P0 T0 b1 ¼ 1:54 þ 0:026ðf 1Þ b1 ¼ 0:43 þ 0:003ðf 1Þ.
ð3Þ
Here Sl is the laminar burning velocity, Sl0 is the laminar burning velocity at reference conditions, P0 and T0 are the pressure and temperature at reference conditions, respectively, and b1 and b2 are the exponential coefficients. To make sure that using this laminar burning velocity model is appropriate, validation calculation was performed. The calculation configuration is a small cylindrical container whose radius and height are 0.106 and 0.111 m, respectively. The 30% hydrogen–air flammable mixture is filled inside this container. The ignition source location is at the center of this container. The comparison of time history of overpressure between experimental result (NTS, 2005) and calculation result is shown in Fig. 2. Both experimental and calculation results exhibit a similar behavior, and the pressure increases progressively until it approaches a maximum value pmax. The calculation roughly predicts the experimental result, the calculation pressure is higher than the experiment. This is due to the large heat loss at the container wall in the experimental result. The rate of pressure rise of the calculation result roughly agrees with the experimental result. The KG values of the calculation and that of the experiment are 665 and 385 bar m/s, respectively, and the calculated result is higher than the
ARTICLE IN PRESS A. Hashimoto, A. Matsuo / Journal of Loss Prevention in the Process Industries 20 (2007) 455–461
457
Fig. 3. Configuration of calculation model and the ignition source location.
Fig. 2. The comparison of overpressure between calculations result and experimental result.
experimental value because of the large heat loss at the walls in experiment. From these results, using this model is considered much reasonable and proper in this study. Various turbulent burning velocity models, which depend on unburned physical quantity and sub-grid scale turbulence are available for LES. In this study, turbulence burning velocity model reported by Flohr and Pitsch (2000) was employed. 3. Calculation model Calculation model of the two rooms connected by ducts is shown schematically in Fig. 3. It is a large-scale building which consists of two rooms connected by two ducts, and aims to model such structures like atomic power plant and chemical facility. In this figure, to the left side is a cubic room (Room A) with size H W L ¼ 5.0 5.0 5.0 m3 and to the right is another cubic room (Room B) with the same size. Each room is connected by two ducts whose cross-sectional area is H W ¼ 0.6 0.6 m2, located in the sidewalls. To reveal the gas explosion phenomenon inside these two rooms, two kinds of simulation were performed. First, to clarify the effect of the ignition source location, three locations were chosen. Point 1 is at the center of the room, point 2 is above the point 1 (the closest location to ducts among the three locations), and point 3 is located below the point 1 (the farthest place from ducts). Then, to clarify the effect of the volume of the ignited room (Room B), the length X of the Room B is changed in the range from 5.0 to 25.0 m. All calculations were performed from equivalent hydrogen–air mixture initially at rest inside this enclosure, with pressure set equal to the atmospheric pressure of P ¼ 101,325 Pa and temperature set equal to T ¼ 293.15 K.
Ignition was initiated by numerically positioning the burnt mixture’s density, temperature, and gas concentration at the location of the ignition point. All walls were assumed as no-slip adiabatic walls, and calculation grid is a structured grid with each cell length being 0.1 m. To reduce the calculation load, symmetric boundary condition was used for the center plane of the room which made the actual calculation zone to be half of the actual size. The discretization method of the convective term is the third-order space accuracy simple high-resolution upwind scheme SHUS, time integral method used is three-step Runge–Kutta method for all calculations. 4. Results 4.1. Effects of the ignition source location First, to verify the basic effect of the ignition source location and the volume of the ignited room, the simulation of gas explosion in one room is performed. The cubic room whose edge is 5 m long is the calculation target and the ignition source location shown in Fig. 3 are employed. To verify the basic effect of the volume of the ignited room, the edge of room X shown in Fig. 3, is changed to 15.0 and 25.0 m. Fig. 4 shows the time history of the pressure measured at the center of the partition. The pressure increases progressively until the pressure approaches towards the maximum value Pmax in all cases. The time of approaching the maximum pressure depends on the ignition source point, and it is fastest at point 1 and latest at point 3. Initially, the chemical energy release at point 3 is lower than those of the other cases since walls near the point 3 prevent the flame propagation. This makes the rate of pressure rise the lowest when the ignition source locations is at point 3. From Fig. 4, it is also apparent that the time of approaching the maximum pressure depends on the volume of the room. From the above discussion, the ignition source locations and the volume of the room are significant factors that affect the behavior of pressure rise in the rooms. Next, the simulations of gas explosions inside the two rooms connected by ducts are performed. The purpose of this analysis is to reveal the effect of the presence of ducts on the room’s pressure history. The series of flame
ARTICLE IN PRESS 458
A. Hashimoto, A. Matsuo / Journal of Loss Prevention in the Process Industries 20 (2007) 455–461
propagation are shown in frames a through f of Fig. 5. The flame propagates squarely in Room B and enters the ducts from the sidewall in Figs. 5a and b. The flame is accelerated by the ducts and flows into Room A, colliding with each other at the center of Room A in Figs. 5b–d. Then, the flame burns out at the upper part of Room A, while in Room B, the flame reaches the walls in Figs. 5e and f. The time history of the pressure at the center of the partition in each room is shown in Figs. 6(a) ignited at point 1, (b) ignited at point 2, and (c) ignited at point 3. The differential pressure in these figures is defined as the difference between the pressures in Rooms A and B, PAPB. Labels on frames a through f in Fig. 6(a) correspond to the same instants after ignition labeled in Fig. 5. The pressure in Room B increases progressively
Fig. 4. Dependence of pressure histories of the center of the partition in Room B with various ignition sources and room volume.
(between labels of a and b). The pressure in Room A also rapidly increases after the flame enters the room, and the pressures of both rooms become almost equal (between labels b and e). After that, the pressures in each room increases at the same rate and approach the maximum pressure (between labels e and f). The behaviors of the pressure increase are almost the same for the cases in Figs. 6(a) and (b). The most remarkable point in Fig. 6 is the pressure in Room B. It is almost the same or lower than that of the one room case when the ignition source location is at points 1 and 2 at all time, however, it is higher after about t ¼ 0.6 s in the case ignited at point 3. These imply that the ducts have an effect of diffusing the pressure in Room B for the cases in Figs. 6(a) and (b), whereas they have the effect of accelerating the pressure rise for the case in Fig. 6(c). In Fig. 6(c), the differential pressure approaches about 100 kPa, while it is almost equal to zero except the early negative phase in Figs. 6(a) and (b). To reveal the reason for this difference, the time history of the differential pressure for the cases when the ignition points are at points 1 and 3 are shown in upper parts of Figs. 7(a) and (b). For the lower parts of Fig. 7, the total mass and the burnt mixture mass in each room are shown. The total mass transportation implies that the total energy transportation is taking place, and amount of burnt mixture implies that the chemical energy transportation is taking place. Both the total and burnt mixture masses are normalized by the initial mass in each room. It is seen that the total masses hardly change for the two rooms in Fig. 7. It may be said that this is not the major cause for the difference in the values of differential pressures. On the other hand, burnt mixture mass in Room A increases and becomes higher than Room B around 0.55 s in Fig. 7(b), while in Fig. 7(a) it becomes almost the same value in Room B, and never overcome it. From this, it can be concluded that the chemical energy transportation by the burnt mixture mass
Fig. 5. Simulated dynamics of flame propagation when the ignition point is the center of the room (point1).
ARTICLE IN PRESS A. Hashimoto, A. Matsuo / Journal of Loss Prevention in the Process Industries 20 (2007) 455–461
459
Fig. 6. The time histories of overpressure in each room and differential pressure (a) point 1, (b) point 2, (c) point 3.
Fig. 7. The time histories of overpressure in each room and differential pressure (a) point 1, (b) point 3 upper: differential pressure; lower: total mass and burnt mixture mass in each room.
is the main factor for the large difference in differential pressure caused when the ignition source location is at point 3 when the volumes of both rooms are equal. 4.2. Effects of the volume of the ignited room To verify the effect of the volume of the ignited room, simulations are carried out by changing the edge length of the ignited room X as shown in Fig. 3 to 15.0 and 25.0 m. The series of flame propagation are shown in Fig. 8, frames a through f. The flame propagation in Room A is almost
the same as in Fig. 5, which is the case when Room B is a cubic room. On the other hand, the flame in Room B propagates squarely in Fig. 8a and is accelerated from the partition in Room B to the right wall in Figs. 8b–f. The flame in Room A propagates within the whole room faster than that in Room B. Figs. 9(a) and (b) respectively, show the time histories of the pressure measured at the center of the partition when X is 15.0 and 25.0 m. In Figs. 9(a) and (b), these time histories exhibit a similar behavior. In Fig. 6(a) X ¼ 5.0 m, the ducts had the effect of diffusing the pressure in Room B.
ARTICLE IN PRESS 460
A. Hashimoto, A. Matsuo / Journal of Loss Prevention in the Process Industries 20 (2007) 455–461
Fig. 8. Simulated dynamics of flame propagation when X ¼ 25.0 m (Room B is five time as large as Room A).
Fig. 9. The time histories of overpressure in each room and differential pressure (a) X ¼ 15.0 m, (b) X ¼ 25.0 m.
However, the ducts accelerate the pressure rise in the cases for Figs. 9(a) and (b). This implies that when the room volume is made larger, the effects of ducts change from diffusion of the pressure to its acceleration. As done in the previous section, when the pressure of the Room A is compared to that of the one room case, the pressure in Room A is higher in Figs. 9(a) and (b). This makes it apparent that the ducts have the effect of accelerating the pressure rise in Room B for these cases. Comparison between Fig. 6(a) (X ¼ 5.0 m; Room B is a cubic room) and Figs. 9(a) and (b) (X ¼ 15.0 and 25.0 m; Room B being an oblong room), two remarkable features are revealed. One is the periodical pressure oscillation in Room B. This is caused by the reflection of the compressive waves on the walls. When the room shape is rectangular, the generated compressive wave is more likely to be enhanced by the gas flow induced in the longitudinal
direction than the cubic room case, due to the flame propagation behavior. The other remarkable feature is the pressure peak of Room A. This is because the chemical energy supply in Room A falls below the diffusive energy transported from Rooms A to B. Upper part of Fig. 10 shows the time history of the differential pressure, and lower part of Fig. 10 shows the total and burnt mixture mass in the case of X ¼ 25.0 m. The high positive differential pressure causes the decrease in the total mass of Room A at 0.5 s, while the burnt mixture mass increases. Then the total mass in Room A becomes equal to the burnt mixture mass at about 1.5 s. This indicates that all substances are burned in Room A, and therefore the chemical energy supply in the room is equal to zero after this point. This allows the pressure to be decreased by the energy transport of total mass from Rooms A to B, creating the pressure peak in Room A.
ARTICLE IN PRESS A. Hashimoto, A. Matsuo / Journal of Loss Prevention in the Process Industries 20 (2007) 455–461
461
peak appeared in Room A and periodical pressure oscillation was also caused. It can be concluded that the ducts have a possibility of causing considerable damage and the room size has influences on the rate of the pressure rise in such explosions inside a closed space simulated in this study. Acknowledgments This study was conducted as one of the advanced basic development of JAEA (Japan Atomic Energy Agency). We would like to thank JAEA for all the supports. Reference
Fig. 10. The time histories of overpressure, total mass, and burnt mixture mass in each room (X ¼ 25.0 m) upper: differential pressure; lower : total mass and burnt mixture mass in each room.
5. Conclusions In this study, numerical simulations of gas explosions inside two rooms connected by ducts are carried out in real scale model. The rate of the pressure rise decreased when the ignition source location was at the center of the room, and in other cases increased. When the size of the ignited room (Room B) was larger than the other room (Room A), the rate of pressure rise in ignited room was lower than the other. This caused the pressure force to act on the partition to the direction of the ignited room. Here, the pressure
Firweather, M., Hargrave, G. K., Ibrahim, S. S., & Walker, G. (1999). Studies of premixed flame propagation in explosion tubes. Combustion and Flame, 116, 504–518. Flohr, P., & Pitsch, H. (2000). A turbulent flame speed closure model for LES of industrial burner flows. Center for Turbulence Research Proceedings of the Summer Program. Holbrow, P., Lune, G. A., & Tyldesley, A. (1999). Dust explosion protection in linked vessels: Guidance for containment and venting. Journal of Loss Prevention in the Process Industries, 12, 227–234. Iijima, T., & Takeno, T. (1986). Effects of temperature and pressure on burning velocity. Combustion and Flame, 65, 35–43. Lenormand, E., & Sagaut, P. (2000). Subgrid-scale models for large-eddy simulations of compressible wall bounded flows. AIAA Journal, 38(8). Miura, A., Mizukaki, T., Shiraishi, T., Matuo, A., Takayama, K., & Nojiri, I. (2004). Spread behavior of explosion in closed space. Journal of Loss Prevention in the Process Industries, 17, 81–86. NTS planning department. (2005). The corpus of hydrogen technology. (Vol. 2) NTS (InJapanese). Tominaga, T., Ito, Y., Taniguchi, N., & Kobayashi, T. (2002). Numerical simulation the flame propagation inside gas turbine combustor. Production Research, 54(1), (2) (in Japanese).