Patrick Bonnet
email: patrick.bonnet@ineris.fr

Jean-Marc Lacome

AN EXPERIMENTAL SET-UP TO STUDY ACCIDENTAL INDUSTRIAL LPG RELEASES

Keywords: atmospheric dispersion, large-scale experiments, LPG, modelling

The objective of the atmospheric dispersion research project of the INERIS is to develop models of flashing releases as encountered in realistic industrial environments. Equivalent source models exist for flashing release in current long rang dispersio n models. Several factors can, however, invalidate simplified equivalent source models, especially in the very near field where obstacles can be found. To perform his project objective, the INERIS takes part in a European project called FLIE (Flashing Liquids in Industrial Environment) and also works on a project supported by the French ministry of the environment. In FLIE project, INERIS carries out large-scale trials of propane and butane releases. The French part of the project is to develop a tool able to evaluate the gas and the liquid fraction (aerosol part and liquid pool) in the near field of the flashing release. This INERIS paper aims to present the large scale experimental set-up and the main current results. The • Rain out : Six bonds equipped with continuous weight measurements are used to study the phenomena of liquid pool formation.

• Obstacle : Six surface temperature.

To develop models requires experimental validation data. Up to now, little of experiments were carried out on a great scale with such a variety of measurements. The obtained results reveal some points which are few discussed in the literature and can be useful to develop models of flashing releases.

experimental set-up is located in the INERIS site. It allows to perform propane and butane liquid releases at ambient temperature with a regulated pressure from the saturation pressure to 15 bar with an orifice (circular or rectangular shape) of an equivalent diameter from 10 mm to 25 mm. It is possible to realise free jets but also impinging jet by introducing obstacle at a maximum distance of 2 meters from the release point.

During releases, several parameters are recorded :

• Ambiant conditions : direction and speed of the wind, temperature, humidity and atmospheric pressure.

•

Release tank : pressure (the regulation of the pressure is possible), temperature at several heights in the tank and at the liquid /gas interface, weight of the tank.

• Release point : pressure and temperature.

•

In the jet : a Dual Phase Doppler Anemometer allows to measure speed and size of particle aerosol at several locations in the jet.

Introduction

Liquefied gasses are in widespread use in industry due to their chemical or physical properties. The substances involved may be flammable and explosive like LPG and also toxic and corrosive like ammonia. In parallel, numerous industrial processes use fluids that are overheated compared with ambient temperature conditions. Following a failure of the confinement around an installation, the hazards faced can be instantaneous, as with the explosion of flammable products, or differed as is the case with releases of toxic products.

As

What happens in the event of an accidental release ?

The consequences of an accidental release will depend on the one hand on the initial storage conditions, and also on the environment. This is because the presence of obstacles of sufficient size may affect the jet and change the characteristics of the toxic or flammable cloud.

At the breach level, in the case of an accidental release of liquefied gas, the product released is suddenly placed in under such temperature and pressure conditions that a part of the liquid vaporizes violently. This phenomena is generally called a "thermodynamic flash" or "flash". The sudden vaporizing of part of the liquid leads to a fragmenting of the liquid jet into fine droplets whose diameter and speed will be especially dependent on the deviation between the initial product temperature and the boiling point of the product at atmospheric pressure.

This phenomena is illustrated in figure 1.

In the absence of any obstacle, after a length of a few times the diameter from the breach, a two-phase jet develops. This jet is made up of droplets in suspension in a mixture of gas and ambient air that is drawn along. The atmospheric dispersion of this two-phase jet is notably dependent on the environment. This is illustrated in the figure 2.

In cases where no obstacles are present along the path of the jet, the different contributions to the formation of a flammable or toxic cloud are as follows:

• a quantity of product is likely quasi-instantaneously vaporized due to the flash phenomena at the breach level;

• a quantity of product will likely form a pool on the ground (commonly called rain-out) that then evaporates into the atmosphere due to the various phenomena of mass and heat transfer ;

• a quantity of product is likely drawn into the air in aerosol form.

In this case, the estimation of the main amounts of vapor and liquid consists in studying the evaporation of the formed drops and the formation of a liquid pool at ground level. This estimation is an important aspect to determine a realistic source term to provide input into an atmospheric dispersion model.

In the case of a release in the presence of obstacles, the interaction between the jet and the obstacle may notably change the calculation of the source term which will be more difficult. Actually, in addition to the three contributions mentioned previously, three phenomena could occur :

• a part of the liquid will be captured by the obstacle (run-off) ;

• a part of the droplets will evaporate in contact with the obstacle (at least during the first instants, before the obstacle cools down) ;

• a part of the drops splits, leading to the creation of smaller sized drops.

In the case of an impacting release, it is important to determine what becomes of the drops after impact. To do this, it is necessary to know the characteristics of the drops prior to impact (diameter and speed) as well as those of the wall (temperature, roughness, presence or absence of any liquid film).

The experimental set-up

The experimental set-up, developed in order to understands this phenomena, is located in the INERIS site.

The figure 4 shows the release field and the control room in the background. The figure 5 is a diagram of the experimental set-up.

The experimental set-up allows us to perform propane and butane liquid releases (the release point is presented in figure 5) at ambient temperature with a regulated pressure from the saturation pressure to 15 bar with an orifice (circular or rectangular shape : see figure 6) of an equivalent diameter from 10 mm to 25 mm.

It is possible to carry out free jets but also impinging jet by introducing obstacle at a maximum distance of 2 meters from the release point.

During releases, several parameters are recorded :

• Ambiant conditions : direction and speed of the wind, temperature, humidity and atmospheric pressure.

• Release tank : pressure (the regulation of the pressure is possible), temperature at several heights in the tank and at the liquid /gas interface, weight of the tank.

• Release point : pressure and temperature.

• Inside the jet : a Dual Phase Doppler Anemometer allows to measure speed and size of particle aerosol at several locations in the jet.

• Rain out : Six bonds equipped with continuous weight measurements are used to study the phenomena of liquid pool formation.

• Obstacle : Six surface thermocouples.

Some Results

This section deals with results obtained with propane release experiments.

Jet structure

Some results can already be presented relating mainly to jet structure :

• right after the expansion zone, the small drops are mainly found along the center whereas the largest ones are the most numerous around the edges;

• moving away from the expansion zone, an asymmetry appears progressively between the upper and lower parts of the jet; it appears that the largest drops accumulate under the jet, this is the consequence of gravity effect. ;

• It is possible to find a similitude law to determine the longitudinal speed component within the jet.

This law depends on the release configuration.

The initial results tend to show that what happens to the drops present within the jet is primarily influenced by the evaporation, by the interaction between jet turbulence and gravity. The importance of each of these influences depends apparently essentially on the size of the drops.

In other words, the measurements seem to show that there are approximately two types of behavior for drops present within the jet and which are at least partly dependent on their size :

• for the finest drops, the behavior is characterized by convection within the jet, with drops evaporating within the jet;

• for the largest drops, inertia and the gravity mean that they are less likely to follow the jet and therefore better able to exit the jet.

The calculation model developed by INERIS seems to match the results obtained in the context of this initial set of tests. The main difficulty always lies in determining the initial speed for starting the calculations. This determination will be one of the main challenges of subsequent trial runs. This is because in the next set of trials, we plan to focus the measurements (with a finer coverage) on the area closest to the leak section, in order to explore the expansion zone, as far as possible.

Mean Diameter

The figure 7 shows the number of droplets (probability function ) in the observed class versus droplet diameter. The table 1 shows the effect of taking into acount or not the greatest droplets in order to calculate the mass mean diameter.

According to these results we can see :

• it is necessary to analyse the results in terms of mass. This is an uncommon way of interpreting the results. Actually this is more common to consider that the greatest diameter is 200 microns and most of the correlation are made with this consideration.

• the under-estimation of the number of large drops leads to a speed of evaporation always too high : by using the mean diameter coming from the droplet smaller than 200 microns the diameter lost is 10%, actually a lost of 1 % seems more realistic.

• the peak corresponding to the greatest droplets might be in part the consequences of PDA problems.

The VKI researchers (our FL1E partners) think that the drops exist but there is an over-estimation due to PDA mask. If we assume that the number of wide measured droplets are over-estimated with a factor two, the difference between mean diameter keeps strong.

The problem of the measurement of the rain out

The figures 8, 9 and 10 show some results that come from impacting jet of propane. The trial experiments conditions are :

• storage pressure : 9.6 bars,

• orifice diameter : 10 mm,

• distance of the obstacle distance : 0.86 m,

• ambient conditions : temperature : 22.5 °C,HR : 65%, temperature humidity dew point : 16.2°C.

mcORPORERDiscussion and analysis

Free jet

Water of the atmosphere makes part of the rainout due to a phenomena of jet-entrainment that mix the released propane with moist air. The measured mass in the bonds ( Figure 8) result from an equilibrium between rain out of this mixing and his evaporation.

When all the phenomena (release output orifice, rain out) are in a steady-state, as seen in Figure 10, in terms of temperature, a slope can be observed on the curve. The component evaporation is certainly very weak, during the jet the air could be saturated with LPG vapours. In a first approximation the measured slope could be interpreted as the sum of phenomena that concern the mixing.

Once the jet has stopped, a compressed air is drawn to accelerate evaporation. It seems difficult to know what really happens because of the solid or liquid water state. In addition it could be assumed that LPG has been captured within the ice because this one melts by forming multiple cavities.

After total evaporation of what must be mainly LPG (at least in the propane case), the quantity (mainly water ?) of remaining liquid is measured. This data should be considered more or less as the amount, in the low range, of water that has been mixed in the jet.

Impinging jet

During the impinging jet, the previous phenomena exist, but in addition the obstacle absorb some phenomena and change the kinetic.

The obstacle must be cooled before one can observe liquid fall into the bonds. It is seen that freeze is produced on the obstacle. This ice is likely produced during the jet because of the low temperature of the obstacle.

Conclusions and perspectives

To develop models requires experimental validation data. Up to now, little of experiments were carried out on a great scale with such a variety of measurements. The obtained results reveal some points which are few discussed in the literature and can be useful to develop models of flashing releases.

Measurement of the propane rainout mass is complex. Water of the atmosphere makes part of the rainout due to a phenomena of jet-entrainment, so the rainout mass measured in the bonds correspond to LPG mass plus water mass. The special kinetic of reevaporation from the freeze water makes the estimation of the LPG mass contribution more difficult.

An possible perspective to that work could be to carry out laboratory experimental trials that allow us to control temperature and humidity in the ambient or confine atmosphere.

  part of its activity relating to accident hazards, INERIS must especially determine safety perimeters around industrial installations. A research project was therefore set up. The objective of the atmospheric dispersion research project of the INERIS is to develop models of flashing releases as encountered in realistic industrial environments. Equivalent source models exist for flashing release in current long rang dispersion models. Several factors can, however, invalidate simplified equivalent source models, especially in the very near field where obstacles can be found. To perform his project objective, INERIS takes part in a European project called FLIE (Flashing Liquids in Industrial Environment) and also works on a project supported by the French ministry of the environment. In FLIE project, INERIS carries out "large-scale" trials of propane and butane releases. The French part of the project is to develop a tool able to evaluate the gas and the liquid fraction (aerosol part and liquid pool) in the near field of the flashing release. The FLIE project brings together four partners (GexCon company from Norway, University of Hertfordshire from England, Von Karman Institute for fluid Dynamics from Belgium, Institut National de 1 'Environnement Industriel et des RISques from France) including INERIS for performing large scale tests. The purpose of these tests was to improve the knowledge of the dispersion of liquefied gasses in the presence of obstacles located in the immediate vicinity. This INERIS paper aims to present the large scale experimental set-up and the main current results.

Figure 7 :

 7 Figure 1 : Schematic diagram of a two-phase jet


Jean-Pierre Bigot

Abdellah Touil

Patrick Bonnet

Jean-Marc Lacome

Rain-out investigation: initial droplet size measurement

Keywords: rain-out, loss of containment, flashing liquid jets, thermal fragmentation, droplet size distribution

Droplet size distribution inside water flashing jets and corresponding rainout fraction were measured. Mass distribution showed that a few droplets are "large" (> 150 jam) and count for more than 85 % of the liquid mass in the jet because of their large individual mass. This could be due to incomplete thermal fragmentation. It could explain the rain-out falling near the orifice or pipe exit.

Introduction

The rain-out problem is part of the loss of containment scenario about liquefied gases reservoirs: how much of the released liquid phase will fall on the ground? This fraction will not directly participate to the toxic or flammable cloud.

Experimental data about this topic are seldom [START_REF] Johnson | RELEASE A model with data to predict aerosol Rainout in accidental Releases, Center for Chemical Process Safety[END_REF][START_REF] Hocquet | Loss of containment : experimental aerosol rain-out assessment[END_REF]. Models are either not sufficiently validated because of this lack of experimental data or even not in good agreement with the existing ones (Wheatley 1987[START_REF] Lanello | Aerosol research program: improved source term définition for modelling the ambient impact of accidental release of hazardous liquids[END_REF][START_REF] Epstein | A model of the dilution of a forced twophase chemical plume in a horizontal wind[END_REF][START_REF] Papadourakis | Upper and lower bounds of droplet evaporation in two-phase jets[END_REF], UDM, Johnson 1999).

A former work [START_REF] Hocquet | Loss of containment : experimental aerosol rain-out assessment[END_REF] showed that rain-out forms a few meters down-stream the accidental breach under the horizontal aerosol. However literature predicts small droplets (dd < 150 urn) will be formed through the thermal fragmentation mechanism (Brown 1962, Bushnell 1968[START_REF] Gooderum | Measurement of Mean Drop Sizes for Sprays from Superheated Water jets[END_REF][START_REF] Lienhard | The break-up of superheated liquid jets[END_REF][START_REF] Oza | An Experimental and Analytical Study of Flash-boiling Fuel Injection[END_REF][START_REF] Oza | On the Mechanism of Flashing Injection of Initially Sub-cooled Fuels[END_REF], Razzaghi 1989, Reitz 1990Park 1994, Ramsdale 1998, Witlox 2000-2001). Such small droplets would normally be corresponding author. Tel. +33 (0)4 77 42 00 29; fax : +33 (0)4 77 42 49 96 94; bigot(5).emse.fr entrained by the aerosol. The aim of this work is to understand how this non homogeneous behaviour can occur.

Experimental

We used two different experimental set-ups in order to measure drop sizes, temperature and rainout fraction for known initial conditions T res and P res . Both set-ups used the same test sections. The tested fluid was water.

Test sections: orifices and pipes

We used both stainless steel orifices and pipes (fig. 1). The 2 mm inner diameter orifice (D2) and pipe (D2L100) where mainly used. But some experiments were conducted with D5 and D8 orifices, D5L250 and D8L400 pipes.

Rainout, mass flow rate and jet temperature measurements

The experimental set-up for these measurements is mainly composed of an upstream tank, the test section, jet temperature probes and rainout capture basins (fig. 2. & 3.).

The up-stream tank is intended to maintain the tested fluid in a stable thermodynamic state all along an experiment. It consists of a jacketed stainless steel column (4.5 m high, .257 m inner diameter, .233 m volume). It allows for a few minutes steady-state experiments.

We experimented up to Po = 1.1 MPa and To = 350K (170°C). Heat is provided to the tested fluid by a 30 kW electrical heater through oil flowing in the jacket. A .12 m thick rock-wool insulation around the tank allows for a heat transfer coefficient to the surroundings as low as h = .5 W/m 2 K. Tested fluid is circulated (.7 m3/h) from the bottom of the tank to the upper gas phase. This allows obtaining + .1 K temperature uniformity for the tested fluid inside the tank.

Pressure at the test section level is maintained constant (within ±2 kPa) owing to loop controlled nitrogen introduction.

Jet temperature is measured both on its axis and just above the capture basins through Water captured in each capture basin is gravity transferred to a smaller container and then weighted. Ratio of mass in one capture basin to the total mass released from the upstream tank (obtained from mass flow rate and experiment duration measurement) is a measure of the rainout density function (± 1 % accuracy).

Droplet velocity and size measurements

The experimental set-up for these measurements is mainly composed of another upstream tank, the same test sections and a PDPA system. The test section is connected to the The measurement volume can be moved in the three directions (2.7m*.6m* 1.5m).

Thus different measurement locations can be tested during the same experiment. We generally covered the vertical plane which includes the jet axis as described in fig. 6.

Droplet size measurement accuracy

The PDPA stops when 3000 measures are collected or after 5 seconds if less than 3000 measures were collected. At 25 and 50mm downstream of the D2L100 pipe, the apparatus We estimate the Sauter mean diameter from the sample collected by the PDPA. This sample has to be large in order that the estimated value is close to the actual one. Most of the measured droplets are "small droplets" (d<j < 150 |im).

Figure 7. shows that estimated mean value for the small droplets is correct within a few micrometers when the sample has 3000 droplets. The "larger drops" play a predominant role when estimating d^ for the whole population because of their large V/S ratio. But these drops are relatively few among the whole population (~ 40 à 50 % at T res ~ 105°C, 3 to 10% when T res ~ 165°C). When the apparatus stops because he counted 3000 droplets, the sample of larger ones is still limited (100 and 300 droplets for the two situations depicted on figure zxy). It is the reason why the convergence is slower. We can estimate the uncertainty as:

± 100 |iim for a sample of 3000 droplets at high temperature (T res ~165°C) when the number of large droplets is the lower.

± 20 um for a sample of 6000 droplets at high temperature

Results and discussion

Tested reservoir conditions are summarized in figures 8. (pipes) and 9. (orifices) Mass flow rates from experiments with pipes at temperatures lightly higher than the boiling one (up to 395K) show a complete agreement (± 5%) with purely liquid flow assumption (constant linear friction coefficient: /\. c =0,016). Figure 11. represents mass velocities measured at higher temperatures (T > 395 K). The agreement with the Bernoullilike approach (or Lackmé approach: G = • 2p l ( P res -P sat (T res )} ) is quite good. It is better than using a model for initially saturated fluid [START_REF] Fauske | Flashing flows or: some practical guidelines for emergency releases[END_REF] for example) even at low subcooling (1 to 2.10 4 Pa here). This means that inflow is liquid in the main part of the pipe but there is some form of vaporisation just up-stream its end. (Brown 1962, Bushnell 1968[START_REF] Gooderum | Measurement of Mean Drop Sizes for Sprays from Superheated Water jets[END_REF][START_REF] Lienhard | The break-up of superheated liquid jets[END_REF][START_REF] Oza | An Experimental and Analytical Study of Flash-boiling Fuel Injection[END_REF][START_REF] Oza | On the Mechanism of Flashing Injection of Initially Sub-cooled Fuels[END_REF], Razzaghi 1989, Reitz 1990Park 1994, Ramsdale 1998, Witlox 2000-2001), ). The sensitivity of this fraction to the pressure seems to be low, probably lower than the measurement uncertainty.

Flow type

Rain-out measurements

Droplet diameter measurements; Number and mass distribution

qualitatively
It seems to us that such a large fraction of large drops was never mentioned in the literature up to now especially in case of thermal fragmentation. This could be due to experimental considerations : in order to obtain best measurements where most of the droplets are present, it is often better to limit the measurement range. With our apparatus configuration the signal was very weaker in the [0, 600 um] (mask C) range than in the [0, 200 um] (mask A) range. Another reason is that most of the literature is devoted to very small orifices, generally less than 500 um ( Diesel injection for example).

Change with P res et T res

The distribution which appeared as quasi-uniform in number at the lower temperatures (T res < 125°C) appears increasing with diameter when represented in mass because mass increases as diameter to the cube. The mass representation shows that the same distribution still exists at the higher temperatures for the larger droplets (we define them from the histograms as : dd >325 um). Further observing the mass histograms we are tempted to distinguish a third population: a peak between 150 \Jirn and 325 fj,m is almost absent at the lower temperatures (T res < 105°C), but grows when temperature increases. When temperature is higher than T res = 136°C its mass is more than that of « small droplets ». These comments apply to both D2 orifice and D2L100 pipe.

Conclusion

Two pilot scale experimental set-up allowed us to characterize flashing water jets. We measured mass flow rate, initial droplet size and velocity distributions, plus spatial distribution of rain-out. Test sections were both orifices and pipes, mainly 2 mm inner diameter.

Mass flow rate measurements and spatial distribution of rain-out showed that almost no internal fragmentation in the pipe occurs when reservoir sub cooling is large enough (AP / L > 1 MPa / m). In such conditions both pipes and orifice exhibit the same jet behaviour (rain-out fraction and rainout spatial distribution).

Reservoir pressure (P re s) has not a significant influence on total rain-out fraction as compared to reservoir temperature T res .

All droplet size distributions exhibit a peak in the [0 -H 150 |im] range as described in the literature. But some larger droplets are also detected. These last ones count fore more than 85 % of the liquid mass even when they count for less than 10 % of the droplets number! Quasi-uniform droplet size distribution (almost the same number of droplets in each size class) very near pipe or orifice could result from a very quick primary fragmentation mechanism.

The "small droplets" (d d < 150 (am) form mainly in the few first centimetres with the maximum number reached between 200 and 400 mm downstream the orifice or pipe. Mass histograms show that another peak forms between 150 and 325 um. These two peaks mean that an external secondary fragmentation mechanism exists.

The former observations lead to the following scheme: a primary mechanical fragmentation mechanism exists whatever the conditions (orifice or pipe, high or low temperature). It curiously leads to a quasi-uniform distribution of the droplet diameters, some of these droplets undergo very rapidly (mostly less than 10 mm) a secondary mechanical fragmentation and probably a thermal fragmentation. This leads to a population of "small droplets" (dd < 150 urn). This mechanism is stronger at the higher temperatures.

When temperature increases another fraction of the droplets is concerned with thermal fragmentation. This leads to the "medium droplets" (150 urn < dd < 325 um).

An important residual quantity of droplets do not undergo secondary fragmentation.

They constitute the "large droplets" population » ( 325 um < dd < 625 um).

The large mass fraction of large droplets (da > 150 um).could be an explanation for our former observation that rain-out is heterogeneous in nature [START_REF] Hocquet | Loss of containment : experimental aerosol rain-out assessment[END_REF] : it is composed of droplets which fall under the almost horizontal two phase jet. Our next aim will be to prove that such size distributions can be a basis in order to improve rain-out prediction models. Future work will also have to take into account near saturation reservoir conditions and fluids other than water. T res =163,5°C/P rés = 8,2bar T res =166,7°C/P rés =10,4bar

Pipe or orifice

PDPA measurement points

  20 equally spaced (1 m) Pt temperature probes. Time constant for these probes induces a ± .75K measurement uncertainty. Mass flow rate through the test section is measured by comparing the change of the pressure head above the test section during and experiment to the change when a known amount of liquid is added. 15 capture basins ( 1m * 2 m or 1 m * 1 m) allow rain-out water capture (fig. 2. & 3.).

  up-stream tank through a 3 m long, .013 m inner diameter, void heat insulated flexible hose (fig. 4. & 5.). This up-stream tank allows the same initial thermodynamic conditions (1.1 MPa, 350K) to be obtained as the other one. It is sphere shaped. The useful volume is approximately 1. m j . Internal electrical heating is provided, as well as N2 inlet. Both are manually monitored leading to a typical chang of ± 10 4 Pa and ± .5 K just up-stream the test section during an experiment. A DANTEC Dual-PDPA 1 allows local measurements (typical 1 mm 3 volume) of both two components of individual droplet velocity and diameter. Statistics about droplet diameter lead to size distribution.

  only detects a few droplets (260 et 193 respectively) during the 5 s which are allowed. The difficulty is probably linked to the optical density of the jet in this region : beams generally cannot cross the jet without being refracted. At intermediate distances (generally comprising 200 mm downstream the orifice or pipe) 3000 droplets are measured in less than 5 s. The number of measures during 5 s decreases again farther than 4 or 600 mm. The main reason is probably the decrease in the density of the jet when going farther and farther from the source. Sometimes droplets impact the lenses of the apparatus... We thus decided to use a distance of 200 mm on the jet axis as our reference location. Another justification for this location is that the change in the size distribution is very important in the first few centimetres and clearly slower after 200 to 400 mm.

Figure 10 .

 10 Figure 10. shows mass velocity (W/A) data for orifices (A is the nominal section of

Figure 12 .

 12 Figure 12. et 13. show the spatial distribution for the rain-out downstream both orifice

Figures 17 .

 17 Figures 17. (D2L100 pipe) and 18. (D2 orifice) show that all distributions have a peak

  near a log-normal or a Rosin-Rammler distribution. Our measurements show that another population exists. At low temperature (T res = 80°C, 104° C and 124°C), one observes a quasi-uniform distribution (almost the same number of droplets in each size class) up to 600 um (PLDA limit for the used configuration). At higher temperatures (136 to 167°C) a few droplets with diameter larger than 150 urn still exist, most of them between 200 et 300 um. Droplet mass is more important than droplet number as far as rainout phenomenon is concerned. This led us to represent size distribution using the mass probability density function. Figures 19. and 20. show a completely different trend. The "small droplets" (d g < 150 pirn), which represent the larger number ( T > 90 % for T res > 136°C and still T > 45 % for T res ^ 80°C) only represent 0,5 to 17' % of the mass! And the large ones count for the complement, i.e. 83 to 99,5 % of the mass! The mass fraction for the « large droplets » in our temperature range is always the larger one. It represents almost all the mass when T res < 125°C. It then slowly decreases to approximately 85% at T res ~ 165°C (fig. 21.

  Fig. 22. and 23. represent the probability density function data at different locations

  Brown, R., & York, J.L. (1962). Sprays formed by flashing liquid jets. American Institute of Chemical Engineers Journal, Vol. 8, n°. 2, p. 149-153. Bushnell, D.M., & Gooderum, P.B. (1968). Atomisation of superheated Water jets at low ambient pressures, A.I.A.A, Vol. 5, n°. 2, 231-232. Park, B.S., & Lee, S.Y.(1994). An experimental investigation of the flash atomisation mechanism, Atomisation and Sprays, Vol. 4, n°. 2, 159-179. Ramsdale, S.A. (1998). Droplet formation and rain-out from two phase flashing jets. Rapport AEAT-2124 Issue 1. Razzaghi, M. (1989). Droplet size estimation of two phase flashing jets. Nucl. Eng. Des, Vol. 114,115-124 Reitz, R.D. (1990). A Photographic Study of Flash-Boiling Atomisation. Aerosol Science and Technology, Vol. 12, 561-569. Wheatley, C.J. (1987). Discharge of liquid ammonia to moist atmospheres, survey of experimental data and model for estimating initial conditions for dispersion calculations, SRD, R410. Witlox, H.W.M., & Holt, A. (2000). Unified Dispersion Model. Technical Reference Manual, Version 6.0, DNV. London. Witlox, W.M., & Bowen. P.J. (2001). Flashing liquid jets and two-phase dispersion. HSE Books a review.

Fig. 1 .

 1 Fig. 1.: scheme of the test sections (orifices and pipes)

Fig. 4 .

 4 Fig. 4.: scheme of the set-up for droplet size measurement

Fig. 6 :

 6 Fig. 6 : typical PDPA measurement locations

Figure 10 .

 10 Figure 10. : Mass velocity versus relative reservoir pressure (orifices)

Fig 13 .:

 13 Fig 13. : Spatial distribution of tain-out (D2L100 pipe).

Fig. 14 .

 14 Fig. 14. : influence of pressure gradient on internal 2-phase flow structure

Fig. 16 .

 16 Fig. 16. : rain-out fraction versus reservoir temperature (orifices).

Fig. 17 .Fig

 17 Fig. 17. : Number size distribution at different reservoir conditions (200 mm down-stream from D2L100 pipe)

Figure 5 :

 5 Figure 5 : Photograph of test section and PDPA (no jet).

Fig

  Fig. 6 : typical PDPA measurement locations

Fig

  Fig 12. : Spatial distribution of rain-out (D2 orifice).

FigFig

  Fig. 16. : rain-out fraction versus reservoir temperature (orifices).

Fig. 22 :

 22 Fig. 22: number size distribution at different locations

Fig. 23 :

 23 Fig. 23: number size distribution at different locations

Fig. 24 :

 24 Fig. 24: mass size distribution at different locations

2 300 mW argon lasers in two orthogonal planes, 60X Fiber Flow series, 57X80 and 58N81 detectors, BSA P70 electronic module, mask C.

' aConduiteD2-L100106-C1D,5bar

T=136,2°C/P=8bar T=165°C/P=8bar


Russo, P. M , Joulain, P. * H Y ' Wang

A Carrau

RISK ASSESSMENT OF A FIRE INVOLVING COMBUSTIBLE MATERIALS IN A WAREHOUSE

This paper presents the results of a numerical modeling into a rack fire in warehouse, and an identification of the visible flame shape when a fire is generalized over a whole warehouse after collapse of the roof. The model described here is a numerical implementation of a thermal and fluid-dynamic theory to predict the rate of flame spread. Evolution of the thermodynamic system is governed by the conservation equations for total mass, momentum (Navier-Stokes), energy and species mass fractions. Large-Eddy-Simulation method is used in which the large-scale eddies that govern the mixing of the gases are directly simulated. The sub-grid scale motion is calculated using the Smagorinsky model. The combustion process is assumed to be controlled by diffusion, permitting a mixture-fraction-based modeling approach The radiant flux vector in the energy equation is calculated by integrating the radiation intensity over all directions. The radiation intensity is found by solving a radiative transfer equation without scattering. A one-dimensional heat conduction equation for the condensed phase temperature is applied. The mass loss rate of volatiles from the solid fuel surface is assumed to obey the simple Arrhenius equation. At first, the work was focused to a fully coupled solution for a single rack fire configuration to both the fluid dynamics and radiative transport equations. The prediction for a single rack fire configuration is in good agreement with the experimental data. Then, the numerical model is extended to the determination of the flame height of a warehouse fire for comparison with the available correlations. It is found that the warehouse fire looks like a cone of flame with a rectangular base.

INTRODUCTION

Warehouse fires are potentially serious. Warehouse fires are usually limited in size but behave as intense emitters of heat, smoke and other combustion products. In the warehousing field, fire safety is fundamental for the protection of human beings and the environment. In this paper, a methodology for simulating a spread fire along the solid surface in rack configuration is presented. This aspect has always occupied a large proportion of the research effort [START_REF] Ingason | Experimental and Theoretical Study of Rack Storage Fires[END_REF][START_REF] Ingason | Heat Release Rate of Rack Storage Fires[END_REF][START_REF] Ingason | An Experimental Study of Rack Storage Fires[END_REF][START_REF] Heskestad | Flame Heights of Fuel Arrays with Combustion in Depth[END_REF] in fire science. Assessing the hazard of combustible materials requires an understanding of the effect of such geometries on flame height and temperature/ velocity fields. However, in the task for obtaining the experimental data or analytical solution for a large-scale warehouse fire, it is not possible to describe the flame-surface heat flux. Therefore, as a first step towards describing fire propagation in rack storage configuration, this paper presents an application of Large-Eddy-Simulation (LES) for solving the thermal and fluid-dynamic equations of 3D elliptic, reacting flow. For this, the Fire Dynamics Simulator (FDS) code, developed by NIST [START_REF] Mcgrattan | Fire Dynamics Simulator[END_REF], is used. The numerical model is validated against the experimental data obtained by Ingason [START_REF] Ingason | Experimental and Theoretical Study of Rack Storage Fires[END_REF][START_REF] Ingason | Heat Release Rate of Rack Storage Fires[END_REF][START_REF] Ingason | An Experimental Study of Rack Storage Fires[END_REF] for a single rack fire, and a good agreement is observed. However, it is now relatively expensive to perform a CFD calculation in the design phase for safety consideration of a warehouse. Therefore, engineering methods have become important in the development of tools to assist in the design procedure. However, in the case of warehouses, fire safety engineering techniques arc less developed [START_REF][END_REF], there is no method available which reveals whether the design is sufficiently safe and economic. As a fire is generalized over a whole warehouse after collapse of the roof, the flame height becomes the most important parameter. Because near the fire, over a distance comparable to the flame height, the radiant energy flux can be sufficiently high to threaten both the structual integrity of neighbouring buildings, and the physical safety of firefighters and personnel. According to McCaffrey [START_REF] Mccaffrey | Purely Buoyant Diffusion Flames: Some Experimental Results[END_REF][START_REF] Mccaffrey | Flame Height[END_REF], the flame can be divided into three regions: the continuous flame region, the intermittent region, and the smoke plume. In his work, various correlations have been developed for the evaluation of the velocities and temperatures in the plume zone. Many analytical studies [9-131 about the pool-like fire of cellulose products or hydrocarbons in a laboratory scale have been also conducted to establish the mean flame height correlations. Heskestad [START_REF] Heskestad | Flame Heights of Fuel Arrays with Combustion in Depth[END_REF] made some modifications to his original correlation for adapting a rack fire configuration, hi particular, Ingason [START_REF] Ingason | Experimental and Theoretical Study of Rack Storage Fires[END_REF][START_REF] Ingason | Heat Release Rate of Rack Storage Fires[END_REF][START_REF] Ingason | An Experimental Study of Rack Storage Fires[END_REF] has undertaken a theoretical and experimental study about the fires in rack storage. He introduced various correlations on the heat release rate, the flame height and the in-rack centreline temperature/ velocity. Large-scale fires (A)>20 m) have also been studied by Koseki [START_REF] Koseki | Large Scale Pool Fires : Results of Recent Experiments[END_REF], and the maximum of the radiation intensity emitted by the large-scale flame was found to vary between 10 and 90 kW/m 2 depending on the fuel properties. The second part of this numerical study aims at evaluating and improving engineering tools for use in warehouse fire safety. The first task of this work is to evaluate the capability of the current flame height correlations to analyse the warehouse fire hazard. It is found that some correlations for flame height are often used apart from the domain of validity for the study of a warehouse fire, in particular for the line type warehouse configuration.

THEORETICAL ANALYSIS

In this work, Large Eddy Simulation (LES) is used to solve the gas phase model, which consists of the conservation equations for mass, momentum (Navier-Stokes), and energy

^ + V.pu = 0, ( 1
) ot p -+ (u • V)u --Vp + pg + V • rij)SC}S , (2)
,V.q r .

(

) 3
OT Dt

The basic idea behind LES is that the eddies, accounting for most of the mixing, are large enough to be calculated with reasonable accuracy from the equations of fluid dynamics. The sub-grid scale motion is calculated by use of the Smagorinsky model [START_REF] Mcgrattan | Fire Dynamics Simulator[END_REF] which relates the unknown subgrid scale (SGS) Reynolds stresses, T.J.SGS > to the local large scale rate of strain r^-Tr*.*»^^,^ Wlth S v = -^ + ^ . (4) 3 2^o x , 5xJ

The eddy viscosity is based on the local length scale, A, and time scale, S y ,

f* l = C 2 i pA i \S,j with C-0.21. ( 5
)
Combustion process is assumed to be controlled by diffusion, permitting a mixture-fraction-based modelling approach [START_REF][END_REF] ot

The oxygen mass conservation equation can be transformed into an expression for the local heat rate using the conservation equation ( 6) for the mixture fraction and a state relation [START_REF] Mcgrattan | Fire Dynamics Simulator[END_REF].

The radiant flux vector in the energy equation is calculated as follows

-V.q r = *r(x)[£/(x)-4;r/ 4 (x)] with t/(x) = f I(x,s)rfQ. ( 7
)

J4n

Here, the radiation intensity is obtained by solving the radiative transfer equation without scattering

V.O7 + Kl = K-. ( 8
)

TC

In the condensed phase, a one-dimensional heat conduction equation is solved. The mass loss rate, rh s , of volatiles from the surface is assumed to obey the simple Arrhenius expression, th s --^e EIR1 (kg/m 2 s). The frequency factor, A, and the activation energy, E/R, are taken as 2.6 10 s and 12000 (K), respectively for this paper carton. The ignition temperature is given as 250 °C. The finite-difference technique is used to discretize the partial differential equations and the associated boundary conditions. All the spatial derivatives are approximated by a second-order central differences scheme and the flow variables are updated in time using an explicit second-order Runge-Kutta scheme [5J.

VALIDATION OF THE NUMERICAL MODELS

At first, the prediction is compared with the experimental data in the rack fire configuration in terms of heat release rate, flame height, temperature and velocity. The schematic diagram of the rack storage arrays examined experimentally by Ingason [START_REF] Ingason | Experimental and Theoretical Study of Rack Storage Fires[END_REF][START_REF] Ingason | Heat Release Rate of Rack Storage Fires[END_REF][START_REF] Ingason | An Experimental Study of Rack Storage Fires[END_REF], and the coordinate system adapted from his work in numerical simulation are shown in Figure 1 .

In the experiment of Ingason The parameters a,/3,a,b can be determined from experimental data. This correlation is a combination of the exponential and power law variations. The fire spread can be divided [START_REF] Ingason | Experimental and Theoretical Study of Rack Storage Fires[END_REF][START_REF] Ingason | Heat Release Rate of Rack Storage Fires[END_REF][START_REF] Ingason | An Experimental Study of Rack Storage Fires[END_REF] into three periods : the incipient period, the fast upward fire spread period and the horizontal fire spread period. In the second period, the fire growth along the vertical channels between adjacent stacks is rather faster than that in the incubation period. However, the third is the most hazardous period. A similar trend between the prediction and experiment as a function of the horizontal spacing, w, is observed. In the experiment, during the horizontal fire spread period, the slope of the curve of the heat release rate increases with the increase of the horizontal spacing, w. However, the numerical model does not correctly predict this trend, and the predicted incubation and take-off periods are shorter as compared to the experimental ones. We believe that the main reasons for this divergence are bound to the assumption taken for carton properties and for the simulation of the operation of the ignition source. Moreover, in the experiment, the paper cartons disappear progressively due to burning. This burning away condition is not taken into account in the numerical simulation, resulting in a big difference for the evolution of the heat release rate between prediction and experiment. A comparison between the predicted flame height and experimental data as a function of Q 2 ' 5 , is plotted in Figure 3 for different horizontal spacing, w. For the prediction, the mean flame height corresponds to the zone where the temperature is higher than 500°C during one second. The experimental data suggests the following correlation for the flame height:

0.343(T ,2/5 [START_REF] Heskestad | Engineering Relations for Fire Plumes[END_REF] which is valid only for L<H. For Q 2f5 /w<50 (kW 2/5 /m), the predicted flame height as a function of the heat release rate, Q 2I \ follows closely the correlation [START_REF] Heskestad | Engineering Relations for Fire Plumes[END_REF]. However, for Q 2l5 /w>50 (kW 2/5 /m), the slope of the predicted flame height is steeper than that obtained by the correlation. The flame height for the smallest horizontal spacing, w=75 mm, takes a different form compared to the three others. It is probably due to influence of the burners because the flame reaches the top of the rack while the burners are operational for w-75 mm. In general, the flame spread velocity increases when the horizontal spacing, w, decreases.

Following the Ingason analysis [1-3], an analytical relation for the in-rack excess centreline temperature, A7, is established as a function of the convective heat release rate, as:

-il/3 2/3
A7' = 28 [START_REF] Delichatsios | Air Entrainmcnt into Buoyant Jet Flames and Pool Fires[END_REF] where z 0 is the virtual source location, defined as: A comparison between the predicted and measured temperature along the rack excess centreline as a function of the convective heat release rate is plotted in Figure 4 for three elevations, z=2.27, 3.57 and 4.87 m. Globally, the evolution of the excess centreline temperature is consistent with a large body of the experimental data and the correlation (1 1). For (z-z 0 )/Q c 2/5 >0.2 m/kW 2/5 , the temperature increases following the power -5/3 of the abscissa. It should be noted that, at the base of the rack fire (z<1.0 m), the temperature is practically independent of the heat release rate. Therefore, the correlation ( 11) is found to be valid only far away from the base of the rack.

An analytical relation of the gas velocity along the rack centreline is also established according to Ingason's experimental data [START_REF] Ingason | Experimental and Theoretical Study of Rack Storage Fires[END_REF][START_REF] Ingason | Heat Release Rate of Rack Storage Fires[END_REF][START_REF] Ingason | An Experimental Study of Rack Storage Fires[END_REF] as a function of the convective heat release rate, as: M = 3.54 [START_REF] Moorhouse | Scaling Criteria for Pool Fires Derived from Large-Scale Experiments[END_REF] A comparison between the predicted and measured gas velocity along the rack centreline as a function of the convective heat release rate is presented in Figure 5 at three elevations, z=2.27, 3.57 and 4.87 m. The numerical model tends to underestimate the centreline velocity, and however, produce a similar trend as compared to the experimental one. In general, the centreline velocity increases as the horizontal spacing, w, decreases. As for the centreline temperature, the correlation ( 13) is not valid at the base of the rack fire (z<1.0 m).

FIRE GENERALIZED OVER A WHOLE WAREHOUSE

It can be seen that in a warehouse with rack storage, the flame reaches rapidly the ceiling as the fire growth rate becomes extremely important. Both quantity and type of combustible materials, once ignited, play a major role on the flame height and on the emission of toxic chemical species. However, the complexity of the problem is such that a full simulation of the scenario including a fire spread in the incipient period, and a potential post-flashover fire encompassing realistic, arbitrary conditions in a large-scale warehouse remains intractable. So far, fire safety design has been highly reliant on prescriptive rules in warehouse codes which specify the maximum allowed distance to a fire exposed personnel. In fact, as the time elapses, collapse of the roof (usually steel) may occur due to heat transfer to the structure from the fire. So, in this study, the foreseeable worst-case scenario has been chosen with the following features. I) a fire is generalized over a whole warehouse with a surface of 7500 m 2 ; 2)

presence of a large number of racks with a height of 10 m in the warehouse; 3) all the surfaces of the solid or liquid materials are burning with a mean mass loss rate of about 7 g/m 2 s. This value of mass loss rate is provided by FDS [START_REF] Mcgrattan | Fire Dynamics Simulator[END_REF] Database file for a large rack fire of plastic products. A schematic diagram of the rack disposition in the horizontal plane of a warehouse, and the coordinate system in the numerical simulation are illustrated in Figure 6. Two warehouse configurations, differing each other by the length to width ratio, IJW, are considered. The length (/>) and width (W) of the warehouse are respectively set to 150 m and 75 m for L/W=2, and 250 m and 30 m for L/W=% with a warehouse height of 13 m. The symmetry conditions are used in the computation model. According to the numerical study (cf. Figs. 78), the visible flame shape of a warehouse fire where the temperature is higher than 500 °C looks like a cone with a rectangular base. For the sides perpendicular to racks, the flame shape can be treated as a fire wall with a surface of Z,^ xZ, m 2 where L ob is the oblique side of the cone (cf. Fig. 8), given as:

2 f , (m). ( 14
)
For the sides parallel to racks, the flame shape can be taken as a triangular fire wall with a surface of 0.5 x HfX Wm 2 . The volume of the cone with a rectangular base can be obtained as:

(m 3 ).

(

) 15
In Eq.( 15), the warehouse walls are accounted for by blocking portions of the radiation through the reduction in the cone volume (configuration factor). However, for the worst-case scenario, the warehouse walls should not be taken into account (H=0). The radiation intensity from flame is obtained by

s (kW/m2) (16)
For the large-scale warehouse fire, the combustion process is controlled by the air entrainment. Based on the oxygen consumption from FDS simulation, the maximum effective heat release generated by a warehouse fire with a burning surface (racks) of 60000 m 2 for L/W-1 can reach 3000 MW for a solid material (PMMA), and 5000 MW for an hydrocarbon fuel (heptane). So, the effective heat release accounts for about 30 percent of the theoretical total heat release, Q c (that is, combustion efficiency ^ = 0.3). In general, the radiative fraction, % r , is assumed to be 0.35 for a large-scale fire. For L/W=2, the warehouse fire can be simulated by a cylindrical pool-like fire with an equivalent diameter obtained from (m).

(

) 17
The equivalent flame surface is calculated by , (m 2 ). (18)

However, for L/W=%, the cone with a rectangular base, as shown in Figures 78, should be applied for the calculation of the view factors (F). The flame surface corresponding to the cone is

S f = 2L*L + LW + WHf, (m 2 ). ( 19
)
A comparison between the predicted flame height and the correlations for L/W=2 with an equivalent diameter of 82 m (Eq.17) is presented in Figure 9. Among relationships for the flame height which have been proposed, McCaffrey's correlation [START_REF] Mccaffrey | Flame Height[END_REF] for the particular purpose of warehouse fire applications appears the most suitable to fit the prediction from the FDS simulation. This flame height correlation is derived from the temperature one as

. (m). ( 20
)
Here, K' is defined as K'= (k/c. The following values in Eq.( 20), 7>= 773 K, T 0 = 300 K, 4= 1.9, c=0.9 and 77 = 0 , are used for calculating the flame height. However, all the other flame height correlations [START_REF] Heskestad | Engineering Relations for Fire Plumes[END_REF][START_REF] Delichatsios | Air Entrainmcnt into Buoyant Jet Flames and Pool Fires[END_REF][START_REF] Thomas | The Size of Natural Fires, 9 th International Symposium on Combustion[END_REF][START_REF] Moorhouse | Scaling Criteria for Pool Fires Derived from Large-Scale Experiments[END_REF] based on a fire diameter give results far away from the predicted flame height. It should be noted that, for a warehouse fire with the presence of a large number of racks, it is difficult to determine an equivalent diameter. Particularly, for a large-scale fire with an equivalent diameter of 82 m, the flame height correlations [START_REF] Heskestad | Engineering Relations for Fire Plumes[END_REF][START_REF] Delichatsios | Air Entrainmcnt into Buoyant Jet Flames and Pool Fires[END_REF] are completely out of its intended range of validity. Based on the flame height to which the radiant energy emitted by flame is directly proportional, the radiation intensity calculated from Eq.( 16) is presented in Table I for the different products.

INERIS has developed an analytical software, FNAP [START_REF] Ineris | Feux de Nappe[END_REF], which can be used for determining the safety areas around a warehouse in fire where humans will bear radiation levels emitted from a flame to 3 kW/m 2 for unprotected people or to 5 kW/m 2 for protected fire-fighters. In order to obtain the local heat flux from this code, some parameters, such as radiation intensity ( <j> 0 ) from flame, transmission in air ( T ) and adds the relevant view factors (F), should be specified 0 = 0 0 Fr, (kW/m 2 ).

Based on the type of combustible materials, the flame height, Hf, from McCaffrey's correlation (20), the radiation intensity, (f> 0 , from Table I, and the view factors, F, from a cone configuration Eq.( 15) can be determined.

x(L) Table 1. Radiant energy intensity emitted by flame for four different products.

CONCLUSION

Globally, the predicted results for a single rack fire are found in good agreement with the experimental data [START_REF] Ingason | Experimental and Theoretical Study of Rack Storage Fires[END_REF][START_REF] Ingason | Heat Release Rate of Rack Storage Fires[END_REF][START_REF] Ingason | An Experimental Study of Rack Storage Fires[END_REF], and follow approximately the correlation trend. According to this numerical study, the flame spread is found rather fast along the solid fuel in rack storage due to chimney effects. For the rack fire, the contribution of radiation to the flame spread is dominant. As a fire is generalized over a whole warehouse after collapse of the roof, the visible flame height is about 3 times the rack (in buming) height for /yjf=8, and however, 7 times the rack height forL/W=2. As compared to the numerical results, McCaffrey correlation (20) can be used to determine the visible flame height for a warehouse fire.

However, there are still lacks in the understanding of how a rack fire is evolved spatially and temporally in an enclosure where the ventilation is rather limited. The results as presented here cannot yet be regarded as complete. It can therefore be stated, that despite the relatively simple analytical solution available today, to meet fire safety design requirement, the procedure still lack of flexibility in their use may lead to unnecessarily expensive warehouse. It would be highly desirable to relate the burning rate to the total heat flux from the flame to the fuel surface for a large-scale fire. Also, future development should relate the soot generation rate to the local oxygen and temperature fields for determining the real radiation fraction.

  kg/m 2 s) c p specific heat (J/kg K) C s Smagorinsky constant D eq equivalent diameter (m) E activation energy (J/mol) g acceleration of gravity (m/s 2 ) h enthalpy (J/kg) H total rack height (m) Hf flame height (m) I radiation intensity (W/m 2 sr) L warehouse length (m) m s mass loss rate (kg/m 2 s) p pressure (Pa) Pr, turbulent Prandtl number q c heat release rate per unit volume q r radiant flux vector (W/m 2 ) Q total heat release rate Q c convective heat release rate Sf flame surface (m 2 ) S,j local large scale rate of strain ( 1/s) T temperature (K) Too ambient temperature (K) u velocity vector (m/s) w horizontal spacing between solid fuel W warehouse width (m) z vertical co-ordinate z (> virtual source location (m) Z mixture fraction p density (kg/m 3 ) p x ambient density (kg/m 3 ) fj, t eddy viscosity (kg/m s) K absorption coefficient (1/m) A length scale (m) a,j3,a,b : correlation parameters k,c,r] : correlation parameters X E combustion efficiency Se, turbulent Schmidt number j r radiative fraction Tij Reynolds stresses

Figure 1 .

 1 Figure 1. Instantaneous snapshot of a rack fire, and the coordinate system in numerical simulation.The evolution of the rates of heat release determined in experiments [1-3J and predicted by use of the numerical model as a function of the horizontal spacing, w, is shown in Figure2. From the experimental data for Q<iO MW, Ingason[START_REF] Ingason | Experimental and Theoretical Study of Rack Storage Fires[END_REF][START_REF] Ingason | Heat Release Rate of Rack Storage Fires[END_REF][START_REF] Ingason | An Experimental Study of Rack Storage Fires[END_REF] proposes the following relation to describe the rate of heat release:

Figure 2 .

 2 Figure 2. Evolution of the convective heat release rate normalized by the height of the rack storage as function of the time for different horizontal spacing, w.

Figure 3 .

 3 Figure 3. Evolution of the non-dimensional flame height as a function of the heat release rate parameter, Q 2 ' 5 for different horizontal spacing.

  Experimental data against correlation[START_REF] Ingason | Experimental and Theoretical Study of Rack Storage Fires[END_REF][START_REF] Ingason | Heat Release Rate of Rack Storage Fires[END_REF][START_REF] Ingason | An Experimental Study of Rack Storage Fires[END_REF]

Figure 4 .Figure 5 .

 45 Figure 4. Evolution of the excess centreline temperature as a function of the convective heat release rate for different horizontal spacing, w, at z =2.27, 3.57 and 4.87m.

Figures 7 -Figure 6 .

 76 Figure 6. Schematic diagram of the rack disposition in the horizontal plane of a warehouse, and the coordinate system in numerical simulation

Figure 7 .

 7 Figure 7. Temperature contours on the median plane perpendicular to racks.

Figure 8 .

 8 Figure 8. Temperature contours in a Figure 9. Comparison of flame heights between corcross-section parallel to racks.relations and predictions for PMMA and heptane.


W Benaissa

S Elgue

M Cabassucf

N Gabas

D Carson

Transposition of an exothermic réaction to a continuous intensifiée! reactor

Keywords: process intensification, chemical reactor, safety, esterification, simulation ANHYDRIDE PROPIONIC ESTERIFICATION 2.1

The aim of this study is thé transposition of an exothermic réaction carried out in a batch reactor to a continuous intensified one. The esterification of propionic anhydride with 2butanol was selected as a case study. The décomposition reactions of thé reactants and thé products were studied by DSC whereas thé heat of solution of thé two reactants were determined by isothermal calorimetry. An extensive expérimental programme was then carried out using thé Mettler RC-1 réaction calorimeter in order to characterize thé thermal and kinetic behaviour of thé réaction under isothermal conditions. Thé thermal conversion was compared with thé molar conversion obtained by simulation based on a kinetics from thé literature. Optimal and safe operating conditions were eventually defined for a continuous heat exchanger type reactor. Dedicated software was used to assess thé feasibility of thé esterification in this reactor by calculating thé conversion and température profiles along thé process and utility lines.

INTRODUCTION

In thé event of loss of control of thé cooling System of a reactor during thé process, exothermic reactions can lead to thermal runaway, which can in turn lead to loss of confinement of toxic and/or fiammable substances. Safety thus is an important operational objective when processes are developed or optimised. Most processes in thé specialty Introduction

The esterification of propionic anhydride by 2-butanol leads to butyl propionate and propionic acid (see fig. l). This synthesis is exothermic: Ubrich et al [START_REF] Ubrich | Optimal feed profile for a second order réaction in a semi-batch reactor under safety constraints ; Expérimental study[END_REF] found a value of AHr = -62.5 ± 1 kJ.mol" 1 . The réaction is relatively simple to carry out in a homogeneous liquid phase. This réaction thus constitutes a relevant case in thé studies relating to thé assessment of thé chemical risk (3).

CH, CH2

"i C

N y \

CH, CH2

Figure 1. Réaction schème for thé esterification of propionic anhydride by 2-butanol

According to Galvan et al. [START_REF] Galvan | The use of neural network for fitting complex kinetic data[END_REF], thé réaction rate is a fonction of catalyst: without catalyst thé réaction is second-order but has a kind of autocatalytic behaviour when sulphuric acid is présent in thé réaction mixture.

2.2

Réaction rate without catalyst When thé réaction is performed with no strong acid as catalyst, thé réaction is slow and follows second-order kinetics, fïrst order in each reactant (4):

r = k C mhp • C u \ 1]
with jfc, = 5.36178.10

2.3

Réaction with sulphuric acid as catalyst In thé présence of sulphuric acid, Zaldivar et al [START_REF] Zaldivar | The use of neural networks for thé identification of kinetic functions of complex reactions[END_REF] made thé following observations: . The réaction rate seems to be proportional to thé acid concentration, . The réaction rate increases with propionic acid concentration, causing a kind of autocatalytic behaviour (see équations [eq. 4] and [eq. 5]), . After having reached a certain concentration, propionic acid no longer influences thé réaction rate.

Since thé various theoretical réaction pathways are complex, an empirical model was devised assuming thé existence of two catalysts (caû,cat2) [START_REF] Dhanuka | Kinetics of thé liquid phase esterification of carboxylic acids with alcohols in thé présence of acid catalysts : reinterpretation of published data[END_REF]. The transformation of thé initial catalyst was developed by taking into account acidity function [START_REF] Haldar | Expérimental studies on parametric sensitivity of a batch reactor[END_REF].

Thus, thé réaction rate of thé main réaction can be written as (4) :

r x = (*, + k 2 • C catl )• C a n hp • C 2b + k 3 • C cat2 • C [eq. 3]
The réaction rate due to thé formation of thé second catalyst is also taken into account:

Lastly, thé expression of thé acidity function is:

Mr = -(A • C caa + p 2 • CacP)-[ p 3 +1 1 I [eq. 5]
Réaction rate constants follow Arrhenius law:

k t =ko i • exp(-Ea j /(R • Tr)) [eq. 6]
The kinetic parameters are given in table 1. A fïrst séries of experiments was carried out on thé pure reagents and thé réaction products. The samples were heated from room température to 400 °C at a rate of 5 °C/min. Table 2 gives thé onset température and thé température at thé peak as well as thé energy released during thé décomposition. Thèse values are compared with thé référence data bases available in thé literature. Following thé same operating conditions, a complementary test was carried out using a mixture of 1.10 mg for each of thé two reagents. The DSC curve obtained is shown in Figure 2: heat flow released during experiment (expressed in W.g" 1 ) is plotted against température (expressed in °C). 9) noticed an endothermic effect when propionic anhydride is mixed with 2butanol. The C80 isothermal calorimeter was used to déterminé thé heat of solution. 2.10" 3 mol of each reactant was placed in thé calorimeter in separate compartments. After thé calorimeter was stabilised at 40 °C, thé reactants were mixed together inside thé calorimeter. The non-catalysed esterification réaction rate is very slow, thus thé heating effect observed is due to thé solution. The results show that thé enthalpy of solution is 22.1 J.g" 1 .

KINETIC MODELLING AND SIMULATION

The kinetic model is essential in predicting thé behaviour of thé continuous intensified reactor. A computer simulation program was used to validate thé kinetic model suggested in thé literature and detailed in paragraph 2.

The standard 2 litre Mettler RC-1 reactor in thé semi batch mode was modelled. The reactor was filled with thé specified quantity of 2-butanol with or without sulphuric depending on thé experiment. The simulation takes into account thé step of adding thé second reagent (propionic anhydride) and thé second step where thé réaction was allowed to go to completion after thé addition. The foliowing assumptions were made:

. The mixture température is maintained constant by thé jacket régulation, . The System is agitated perfectly: température as well as chemical and physical properties of réaction mixture are uniform, • Réaction mixture is a homogenous liquid phase, . The variation of volume due to thé mixture and dilution are neglected.

The dynamic formulation of thé model leads to a system of differential and algebraic équations. The simulation produces thé concentration profiles and thé volume vs time.

EXPERIMENTAL STUDY IN THE RC-1 CALORIMETER

Expérimental procédure

Ten experiments were carried out in thé réaction calorimeter RC1 for jacket températures ranging from 30 to 70 °C and sulphuric acid concentrations ranging from 0.3 % to 0.8 % (expressed as a percentage of thé weight of butyl alcohol). The 2-butanol was introduced into thé reactor as well as thé proper quantity of sulphuric acid. The stirring speed was fixed at 200 RPM. When thé mixture température was equilibrated, 6.86 mol of propionic anhydride was pumped into thé reactor at a constant flow rate over a 15 minutes period. During thé réaction, thé réaction température was maintained constant by controlling thé température jacket. However, thé second reagent was added in several times when thé instantaneous heat release was likely to exceed thé cooling capacity of thé calorimeter.

Expérimental results

The results were recorded by thé WinRC® (10) software. During acquisition, it records on line ail thé pertinent thermal parameters. After thé experiment, it calculâtes thé heat flux due to thé réaction as a fonction of time (Qr in W) and calculâtes thé enthalpy of réaction for a spécifie area chosen by thé user (ÀHr en kJ.mol" 1 ). The software also calculâtes thé thermal conversion (X en %) within thé intégration limits. If we assume that heat released is due only to thé main réaction, thermal conversion can be compared with molar conversion. Thus, thé expérimental thermal conversion is compared with thé molar conversion obtained by simulation. An example is given in figure 3. The curves are in good agreement. We thus deem thé kinetic model valid.

ESTERIFICATION IN THE CONTINUOUS INTENSIFIED REACTOR

6.1

Heat-Exchanger/Reactor OPR (Open Plate Reactor) The thermal and kinetic behaviour of thé esterification of propionic anhydride by 2-butanol has now been characterized. The objective is then to study thé operating conditions of this réaction in thé heat/exchanger OPR (Open Plate Reactor) developed by Alfa Laval Vicarb. The OPR is a small, multifunctional, continuous reactor which is built like a multi-plate heat exchanger divided into sections (11). As shown in figure 4, each section is made up of a reactive plate where thé réaction mixture flows, surrounded by two cooling plates containing thé utility fiuid. This new concept of "reactor/ heat exchanger" makes it possible to perform complex chemical reactions with a very accurate thermal control. The OPR appears particularly well-suited to process intensification, as it makes possible an increase of reactant concentration and a réduction of solvent consumption at thé same time. The pilot is composed of three sections and has a total capacity of 1.5 L. Previous studies hâve shown that thé OPR can be consider as plug flow for a flow-rate of reagents equal to 50 L.h" 1 . The operating system allows a maximum cooling fluid flow-rate of 5 m 3 ^"1 . §-. -• ;.

Spécifie simulations

The use of thé data processing tool makes it possible to asses thé feasibility and potential of carrying out thé esterification in this reactor. The main parameters to be modified are thé initial sulphuric acid concentration, inlet température and cooling fluid flow-rate. It is then possible to déterminé operating conditions in terms of molar conversion while maintaining proper heat transfer and an overall température lower than 100 °C. One resuit is presented in figure 6 and shows that for a reactor composed of three sections, corresponding to a 3 minute résidence time, conversion can reach 45 % for this réaction. A kinetic model taken from thé literature for thé production of butyl propionate by reacting propionic anhydride with 2-butanol has been validated by experiments done in a 2 litres, Mettler RC-1 semi-batch réaction calorimeter. The kinetic model was then used as input to a computer simulation of an Open Plate Reactor to déterminé thé feasibility of using this new compact continuous reactor to produce butyl propionate. The operating conditions were optimised using thé results of thé computer simulations and show that thé transposition is indeed feasible. Expérimental vérification in thé OPR pilot Reactor are under way.

Figure 2 .

 2 Figure 2. Microcalorimetry DSC for a mixture of propionic anhydride and 2-butanoI 3.2 Isothermal Calorimetry C80 Alos et al. (9) noticed an endothermic effect when propionic anhydride is mixed with 2butanol. The C80 isothermal calorimeter was used to déterminé thé heat of solution. 2.10" 3 mol of each reactant was placed in thé calorimeter in separate compartments. After thé calorimeter was stabilised at 40 °C, thé reactants were mixed together inside thé calorimeter. The non-catalysed esterification réaction rate is very slow, thus thé heating effect observed is due to thé solution. The results show that thé enthalpy of solution is 22.1 J.g" 1 .

Figure 3 .

 3 Figure 3. Conversion profiles results from simulation and WinRC® (Tr=50°C, w SUi =0,3%)

Figure 4 .Figure 5 .

 45 Figure 4. Simplifiée représentation of thé "heat-exchanger/reactor" OPR 6.2Data processing tool A spécifie computer simulation program has been developed in paralleled to thé development of thé reactor. Thé software takes into account thé spécifie geometry of thé intensified reactor and provides thé température profiles of thé cooling fluid as well as thé température profile and conversion profile along thé réaction line. Preliminary réaction tests carried out on thé OPR hâve been used to validate thé simulation framework. In this way, thé more interesting results hâve been obtained from thé study of thé oxidation of sodium thiosulfate by hydrogen peroxide since this réaction présents fast kinetics. By means of température records ail along thé reactor and of thé final réaction yield, like showed in figure5, thé thiosulfate oxidation réaction allows to successfully validate thé simulation framework.

Figure 6 .

 6 Figure 6. Température profiles and molar conversion along thé reactor OPR for thé propionic anhydride esterification

  7

Table 1 . Kinetic parameters équations from [eq. 1] to [eq. 5]. (5)

 1

	i	to,	t(J.mor 1 )	Pi
	1	5.36178 10 7a	80478,64	2,002.1G -1
	2	2.8074 10 10b	79159,5	3,205.10' 2
	3	3.9480 10 10a	69974,6	-21,3754
	4	1.4031 10 8 a	76617,2	12706

L.moî ! .s" ] 3 PREUMINARY CALORIMETRIC STUDIES 3.1 Microcalorimetry DSC DSC

  (Differential Scanning Calorimeter) was used to study thé thermal stability of various samples. Only a few milligrams of each were needed.

Table 2 .

 2 Tonset : température of décomposition (°C) b Tmax : température corresponding to thé maximum energy released by thé decomposition(°C) 0 AHjec : energy of décomposition (J.g") d Tonsetyb : température of décomposition, data from bibliography[START_REF] Perry | Perry's Chemical Engineers' Handbook[END_REF].Thé first peak gives thé energy released by thé main réaction and thé second area corresponds to thé energy of décomposition. According to table 2 and figure2, thé température of thé réaction mixture must not exceed 150 °C to avoid any décomposition.

	Pure Products	Results of thé DSC experiments on pure products. Mass (mg) Tonset a Tmax b (°C) AH de c c (J/g)	Tonset bib	d
			(°C)			(°C)
	2-Butanol	2.54	183.5	187.0	-42	184.6
	Propionic anhydride	2.22	-155	167.5	-47	154
	Propionic acid	1.82	188.4	209.8	-35	186.6
	Butyl Propionate	2.05	-170	174.2-209.2	-70	-

a DSC:

2-Butanol and Propionic anhydride

	200	2S0
	Température (°C)

ACKNOWLEDGEMENT

The authors would like to thank Alfa Laval in supply of thé "heat-exchanger/reactor" and Data Processing Tool software which hâve made possible this study. OPR


Martial Mosnier

Benjamin Daudonnet

Jérôme Renard

Guy Mavrothalassitis

Frédéric Mercier

ENSIB G Mavrothalassitis Lees

Metallic pressure vessels failures

à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Metallic vessels are of common use in thé industrial field: they can be operated at very différent pressure levels to store or to transport gas or pressurized liquid (such as LPG or LNG), to dry, or as steam boiler... etc. As a conséquence, metallic enclosures are also widely involved in many accidents, sometimes with catastrophic conséquences.

Although thé causes of vessel accidents hâve been widely investigated, thé exact rôle of some parameters still remains unclear. Therefore, forecasting thé effect of such an accident on thé neighborhood of thé vessel is usually achieved with thé help of handbooks, that sometimes overestimate effects, such as thé ones of BLEVE. Dynamic crack formation and fracture propagation in thé enclosure wall are among those well identified, but unfamiliar parameters.

Two main causes need to be considered concerning thé vessel rupture. The first one is a fast rise in pressure inside thé vessel due to an overfilling or an overheating of thé vessel. This overheating is itself caused by an external fire, a runaway réaction or an internai explosion. For this reason, thé pressure inside thé vessel exceeds thé dynamic pressure limit of rupture. The second cause of rupture is a réduction of thé vessel strength caused by corrosion, overheating, material defect or external impact. This second aspect will not be investigated in this paper. A sudden rise in pressure, which leads to thé ruin of thé container at thé post-ultimate time, can be due to a human mistake, an uncontrolled chemical réaction, a physical phenomenon like BLEVE, a physical and chemical phenomenon like a rise in mass, an explosion, .... Table I brings together a panel of accidents with their causes, resulting phenomena and conséquences on both material and human.

Accidents diversity hinders thé clear forecast of thé event and its effects on thé vessel environment. Enclosure wall cracking is an important parameter to study: thé conditions of its start, thé way and thé speed at which it propagates can give essential information. If some physical laws already permit to describe some cracking propagation cases, other phenomena such as brittle rupture are unfamiliar. Moreover cracking of wall pressure vessel has to be coupled with fluid discharge.

2. Case STUDY: Bordes accident (09/05/2000) [START_REF] Caumont | Direction des Risques Accidentels, Unité Identification et Hiérarchisation des Risques[END_REF]

Context éléments

This accident involved gas cylinders explosions (525 butane 13 kg bottles, 35 butane 6 kg bottles, 140 propane 13 kg bottles, 77 propane 35 kg bottles). Early in thé morning, a 777 gas cylinders carrying truck was delivering a factory. The driver noticed that one of thé towing tires was in fire. After unsuccessfully attempting to put out fire, he went back for help. The fîrst bottles, exposed to thé fiâmes heat, exploded, causing even more bottles to explode. A security perimeter was put into place. The fire was only under control some 4 hours later. Nobody has been injured. A car wash, a dépôt, offices and some houses were damaged. Some bottle débris were found as far as one kilometer of thé accident.

Analysis of fragments

Fragments and débris total mass was about 9400 kg, with 80% of thé fragments within an 80 meters area around thé truck. Ten more percents were lying between 80 and 100 meters around thé truck. The fragments inventory is not complète, since some projectiles with unknown shapes and weights hâve clearly been witnessed flying. Some of them impacted further than 300 and 480 meters, a distance that corresponded to houses. Some fragments were even found at 8 or 900 meters.

Plotting thé mass vs. thé projection distance for long range fragments (over 100 meters) displays a graph (Figure 1), that can be divided into three zones:

•

The most important ones in terms of samples involves 7 on thé 10 listed fragments from distances between 150 and 250 meters and mass comprised between 2.5 and 8 kg.

•

Two other fragments with a mass lower than 2.5 kg are to be found beyond a 250 meters zone. Other fragments were reported at further distances, but, unfortunately, without describing their shapes. Thé biggest fragments were found at a distance less than 180 meters: they resuit from thé rupture of thé enclosure into two pièces, and, as such, are représentative of a vessel ductile fracture.

Due to their shape, flat fragments hâve fallen much further. In général, they consist in smaller débris (1 to 4 kg) created by seemingly brittle bursts of vessels.

Thé last kind of fragments to be find were security designed parts like flanges and heads, with a mean mass of 2 kg. They fell at very variable distances between 180 to 350 meters, depending on thé stiffiiess of thé pressure rise by thé time they were activated.

Explosion of gas cylinders

LPG bottles that exploded during this accident are commonly used bottle, viz. when full, LPG occupies about 85 % of the space inside thé cylinder and so, liquid phase is in equilibrium with its vapor phase, that détermines thé pressure inside thé enclosure.

Depending on how a liquefied gas bottle is heated, two main processes are susceptible to produce a capacity burst:

• If thé wall in contact with thé gas phase is directly heated by fiâmes, convective heat transfers between gas and wall are insuffîcient to counterbalance thé incoming heat flux. Wall température will locally increase, causing a diminution of its mechanical résistance. The drop off, in conjunction with thé pressure rise, can be such that thé wall can not resist thé pressure inside thé enclosure and fails. This phenomenon results in a fracture, that can be thé initiator of a burst for brittle materials, and thé precursor of a breach for ductile ones. When thé resuit is a breach in thé gas part, a torch fîre can be created. Otherwise, thé burst may produce a fireball with fragments projections. As told above in Par. 2.2, several tens of fragments corresponding to that event were found after Bordes accident.

Modeling of pressure metallic enclosures failures

Bordes accident illustrâtes thé complexity of pressure vessel rupture phenomenon, and thé various conséquences in terms of fragmentations and projections.

• Before its failure, an enclosure can suffer high déformation that will absorb part of thé phenomenon energy; • Wall cracks can also start and propagate, causing an internai pressure discharge, and finally stop; « Fragments shapes are dépendent of thé bottle burst.

Thèse considérations point out thé static or dynamic aspect of thé rupture, that is linked to thé pressure évolution. Consequently, depressurization during thé enclosure rupture has to be taken in account. This clearly indicates thé importance of fluid structure coupling.

Global description of crack and depressurization modeling

Even if several différent models describe pressure vessel failure, few models présent a global approach of thé problem. Many of them insist on thé mechanical point of view, others on thé depressurization, few on thé coupling of both phenomena.

Many différent représentations of crack propagation hâve been developed. The basis of crack development models takes place on thé use of usual failure criteria. Thus thé différence between those models mostly stands on thé diversity of thèse criteria. As examples, some of them display thé maximal stress criterion in which thé crack is supposed to propagate perpendicularly to thé maximal hoop stress. Some other use thé minimal strain density where thé crack propagates in thé direction of the minimal strain density energy. Other criterion of interest in thé crack propagation models is thé maximal energy release rate, in which thé crack is suppose to develop in direction of thé highest energy release value, and so on... The Crack Opening Displacement where thé development of crack is considered linked to its width as also to be mentioned in this review.

Finite éléments can handle crack propagation with two methods as described by [START_REF] Recho | Rupture par fissuration des structures[END_REF]:

•
The first one consists in adding a relaxation node that corresponds to a crack tip; •

In thé second, a zone is creating around a note to delimit thé crack tip.

Studies [START_REF] Needleman | Dynamic crack growth in a nonlocal progressively cavitating solid[END_REF][START_REF] Head | Some methods for thé évaluation of K from finite élément data[END_REF] using finite élément method emphasize thé influence of thé characteristic length and thé mesjï résolution: thé average crack speed increases for finer mesh while it decreases for increasing value of material characteristics such as length, thé nature of mesh influences prédiction of the crack growth.

To reduce thé importance of thèse parameters, [START_REF] Nikishkov | An algorithm and a computer program for thé three-term asymptotic expanson of elastic-plastic crack tip stress and displacement fields[END_REF] proposes a différent method: it consists in determining thé asymptotic solution of stresses near thé crack tip, since thé knowledge of stress and displacement field near crack tip are useful for fracture criteria.

Another method used to déterminé stress field in plastic zone, near thé crack tip, is thé slip-line theory: [START_REF] Isaksson | Crack kinking under high pressure in an elastic-plastic material[END_REF] obtained interesting results with this method.

Ail thèse différent models show that thé mechanical aspect of failure can be well described. However, if crack is considered as thé main aspect of pressure vessel failure, it shall not be separated from depressurization: as internai pressure evolves with crack propagation, it is essential to study this phenomenon.

Depressurization is also characterized by thé fluid state: as fluid can either be gaseous or liquid, (or in a liquid-vapor equilibrium), thé discharge can either be a single or double phase flow.

Among thé différent models that describe depressurization, most of them try to calculate thé discharging rime of thé bottle. To do that, some of them take thé fluid state as a parameter while others only consider a monophasic outflow. The other input parameters can be thé kind of breach (varying or constant) as well as thé pressure and thé température inside thé bottle. The output data can consist in thé discharging time, and also thé pressure and température évolutions inside thé bottle.

The simplest model, elaborated by [START_REF] Xia | A simplified model for depressurization of gas-filled pressure vessels[END_REF], consists in calculating gas température and pressure for a constant flow. On thé opposite side, analytical models from [START_REF] Woodward | Liquid and gas discharge rates through holes in process vessels[END_REF] involve a fluid, liquid or gas, that can changeably flow at either sonic or subsonic speed.

Baum depressurization mode!

So far, différent models hâve been exposed, that only deal with crack propagation and/or depressurization. More in relation with Bordes accident are thé Baum models [START_REF] Baum | Rupture of a gas-pressurized cylindrical pressure vessel[END_REF][11]. They depict a mixed process where a variable breach opening results in thé propulsion of fragments behaving like missiles or rockets.

Basis of thèse models are to be found in thé several studies on rocket velocity and fragments conducted by Baum.

At thé beginning of a circumferential crack, Baum's models depict thé capacity as a cylinder containing either a gas or a diphasic gas liquid médium at internai pressure P (Figure 2a,2b).

Breach opening is introduced through thé first steps of rocket displacement. The breach area is defined as thé space between thé flying part and thé remaining one. An opening time is calculated as thé time necessary to fully open thé breach area, that takes place when rocket displacement has reached half thé enclosure radius. Thus, considering internai pressure much higher than external pressure, thé minimum opening time is equal to (for rupture at one end of thé vessel, Figure 2a):

where M is thé mass of rocket, P is thé internai pressure, r pressure vessel radius.

As above mentioned, this model cannot be compared with a crack propagation model, since no stress and only displacement is involved in thé expression of breach area development However, Baum's model is of great interest with regard to thé idea of coupling, it introduces: mass flow rate evolves with thé breach area.

Fracture and depressurization coupling models

A notion exposed in Baum's model is breach opening time: time is an important characteristic of pressure vessel rupture, since crack propagation time and gas discharge time are in compétition to influence pressure vessel behavior. Indeed, discharge can be finished before failure of vessel has been reached and reciprocally, crack propagation can stop before thé end of fluid discharge. So, to be realistic, a model has to include a crack propagation law in addition to a depressurization law.

Lenclud and Venart coupling model

On thé basis of experiments, Lenclud and Venart [START_REF] Lenclud | Single and two phase discharge from a pressurized vessel[END_REF] developed a crack opening model for a single phase discharge and studied its influence on crack speed depressurization. As previously mentioned, they also reported various two phase flow models. They exposed models for every stage of thé phenomenon as well, viz. crack opening displacement and depressurization. For this last step, single or two phase flow can be considered. The crack opening in the flrst step of rocketing fragments.

They consider thé crack opening as thé resuit of both elastic and plastic displacements.

Total crack opening is supposed to be thé summation of both elastic and plastic displacement. Thus thé crack area is given by:

=nab = na -+ A ( 2
)
where a is thé half crack length, b is thé half crack width, C represents thé plastic component -= C (Fig. 3b) and A thé elastic one (Fig. 3a).

W J

Figure 3a. Elastic déformation. Figure 3b. Plastic déformation.

Crack velocity is then introduced in thé latter analytical expression. Crack velocity results from high speed cinematography experiments. From thé authors analysis, crack propagation is divided into steps, that may be single or double.

As an example, thé crack length 2a can be expressed at each time step as:

0<t<t l : t x <t<t 2 :
t>t n : 2a = 21 f with vi ,V2 crack velocity, C2 experimentally determined constant, and t time.

Depressurization in this model is asserted to be best represented by a single shocked phase flow.

Erdogan coupling model

A more elaborated model coupling fracture dynamics and fluid mechanics has been established by [START_REF] Erdogan | Crack propagation and arrest in pressurized containers[END_REF]: it represents a gas pressurized cylindrical tank rupture. Erdogan noticed thé relevant points of thé failure phenomenon that are to be considered:

•

The possibility for a crack to stop or to resuit in a catastrophic failure dépends on thé nature of thé fluid discharge due to thé internai pressure évolution; •

The distinction between brittle and ductile failure.

He then insisted on thé difficulty to estimate thé boundary conditions at thé breach. Moreover this estimation is complicated by thé crack propagation. To simplify thé problem, he proposed thé foliowing assumptionsi Given thèse hypothesis, thé main aspect of his model stands on thé choice of a dynamic rupture criterion. Such criterion relies both thé dynamic nature and thé crack driving force. Erdogan proposed to express crack accélération and velocity as a fonction of thé crack driving force. Still thé crack driving force can be defined using différent expressions, that are dépendent on thé nature of crack, and consequently, brittle and ductile fractures are to be distinguished.

According to Erdogan, thé most realistic parameters defining thé crack driving force for a brittle fracture is thé stress intensity factor K-Kc, or thé strain energy release rate G-Gc.

Ideally brittle materials (with no inertia effect) see their dynamics fracture criterion as da a relationship between crack velocity -and K: dt

dt K>K C (4)
0 for K<K with d n experimentally determined constants.

However, neither instantaneous velocity changes nor inertia effects can be modeled with this expression. Thus this fîrst criterion (crack velocity) must be complemented a d 2 a second one (crack accélération --):

b x (K-K c ) for K>K c A c ) -U2\-K-~ A c ) JOr A. < A c or ^J h( c ) for K>K C ^2 [-i: c 2 ) 2
/or K<K C with b n experimentally determined constants.

The second formulation (eq. ( 6)) is based on strain energy release rate, but instead of using thé G parameter, K 2 expression is used, since it is assumed that G~K 2 .

As such, K appears to be thé main parameter. The knowledge of thé stress mtensity factor is sufFicient déterminé crack velocity and accélération.

For a pressurized cylindrical shell suffering an axial crack, Erdogan approximated thé stress intensity factor by: ^4ÏÛ

0.5A + 0.6 + 0.4e for Â<5 rQ , (8) 1.761VA-1.9 /or X>5 and:

/I = [l2(l-D 2 j , (shell parameter)

This brittle fracture model proposed by Erdogan includes thé main crack characteristics:

• A crack driving force (thé well-known stress factor intensity); « An expression of crack évolution (velocity, accélération) and thé spécial stress intensity factor calculated for cylindrical vessel.

In addition, this model is provided with experimentally determined constants.

For ductile fractures, a différent crack driving force is introduced to take into account thé ductile behavior of thé material: thé "résistance pressure" (p r ) Thèse équations are naturally solved numerically. End of calculations is either consécutive to a nil crack velocity with thé crack length inferior to a limit (a<0,8 L/2) or to a crack length that reaches this limit.

Erdogan's models provide a good agreement with reality, for some simulated trends are conform to those expected. Thus, results of thé quasi-brittle model are in accordance with a real brittle behavior: no crack stop has been numerically observed, and crack velocity was still increasing, even after crack length reached 80% of thé total cylinder length.

The cylinder length is an important parameter of thé ductile model, that seems to strongly influence crack arrest Actually, thé longer thé cylinder is, thé more crack would tend to stop before reaching thé total cylinder length. Obviously, thé value of constants also influences quantitatively thé results. As for shell parameter, thé choice is not a problem since it is defined from crack length, cylinder radius and length.

However, even if few studies consider thé coupling between fracture and gas dynamics, and if Erdogan's results are interesting, his study is not exempted from reproach. The most évident is thé lack of comparison with experiments: if his model seems to be qualitatively good (thé description of failure behavior is realistic), quantitatively there is no référence. For instance, différent results are presented for différent constants, but Erdogan does not mention which of them would give thé best correspondence with reality.

Conclusion

Post accidentai analysis of pressurized vessels ruptures hâve compiled lots of data on their conséquences, that help understanding thé main causes of their failures. However, thé coupling between thé mechanical aspect of rupture and depressurization remains largely unexplored, especially for a smooth quantification of thé effects in thé environment of thé capacity.

So 5 there is still a place in research to refine thèse aspects, with experiments providing relations between parameters such as pressure, time, brittle point of différent materials, température... This expérimental work could supply data (material behavior, crack propagation and discharge time...) for simultaneous numerical model development, following thé approach initiated by Erdogan.

«Figure 2a .Figure 2b .

 2a2b Figure 2a. Rupture of a gas-pressurized vessel.

Figure 4 .

 4 Figure 4. Initial pressure vessel conditions.

Table I :

 I List of accidents with associated data. (ARIA data base [1], compiled by thé BARPI ; and INERIS).

	Place and year	Type, capacity of the vessels	Used products		Causes of the accident		Resulting phenomena	Material damages	Human conséquences
	Alma (U.S.) (1958)	80 m 3 cylindrical vessel	LPG (butane)	Overpressure		UCVE		?	1dead injured
	Brinkley(U.S.) (1959)	Tank trailer 3 vessels (<60m 3 )	LPG (butane)	Human mistake		BLEVE		Considérable	1dead
	Feyzin (France) (1966)	1200 and 2000 m 3 sphères	LPG (propane)		Human mistake		leak, inflammation, explosion, fragmentation	7	1 dead
	Crescent City (U.S.) (1970)	Wagon	LPG		Sparks due to derailment		BLEVE		Wagon fragments at 500 m	66 injured
	Saint-Amand (France)(1973)	Tank trailer	LPG (propane)		Leak due to road accident, ignition due to friction on ground	BLEVE		Truck projection on a house, propagating thé fire	?
	Kingman (U.S.) (1973)	75 m 3 tank wagon	LPG		Human mistake		Gas inflammation BLEVE	Wagon projection at 365 m; 5 buildings in fire	13 dead (operator and 12 firemen) 95 injured
	Kurtkoy (Turkey) (1997)	aérosol cans and 101 cylindrical vessel	LPG (butane)	Spark due to electricity return after power failure		Explosion of the drums and BLEVE of the vessel	Collars and cylinder heads Windows in a 500 m radius area projection as far as 80 m. Broken	?
	Qatar (1977)	23 0001 cryogénie vessel 50 100	LPG 150	200	Bad quality welds and insufficient maintenance 250 300	350	400	VCE	450	Totally destroyed factory; damages as far as 2 km	7 dead
	Texas City (U.S.) (1978)	800 m 3 sphère	LPG (isobutane)	Overfilling and cracking on welded joint (leak)	BLEVE		Broken Windows in a 3.5 km radius area	7 dead 10 injured
	San-Ixhuatepec (Mexico) (1984)	4 sphères of 1600 and 2400 m 3 , 48 cylindrical vessels	LPG		Rupture of a pipe under 24 bars		BLEVE and domino effects	Fragments further than 1200 m	dead 7000 injured
	Asahikawa (Japan) (1988)	Recompression unit and filling device of gas cylinder	LPG		Procédure overriding resulting in leaking		Explosion (unknown cause of ignition) and fire	1 207 destroyed gas cylinders	3 dead injured
	Sydney (Australia) (1990)	Aerial vessels tank trailers, cylinders ( 160-148-88-3 1m 3 ),	LPG		Ignition source: electrical spark. propagates a fire An opened floodgate	Succession of BLEVEs	Shock wave felt at 3km; cylinders hundreds of small	9
	(South Korea) Ulsan (1990)	Vessel	LPG (butane)	static electricity discharge		Protection System fire and vessel explosion leakage that leads to fire	?	?
	Sainte Sophie (Canada) (2001)	Gas cylinders with capacities up to 180 kg	LPG (propane)		Short circuit		BLEVE		Ail thé vessels (2-300) and a tank trailer	injured

distance (m)

Figure

1

. Fragments mass vs. falling distance (Bordes accident).

  • When thé liquid phase is directly heated or when thé whole capacity is surrounded by fiâmes, thé input flux on thé wall is almost totally transferred to thé liquid phase, and contributes to thé phase change from liquid to vapor, following thé saturated vapor curve. Since for liquefied gas, thé vapor phase spécifie volume is generally about thousand times thé liquid phase spécifie volume, this vaporization results in a huge pressure rise inside thé bottle, that can exceed thé wall résistance. Moreover, if thé capacity wall is surrounded by fire, wall mechanical résistance of thé gas side can fall down at thé same time. Depending on both pressure rise slope and wall material nature, thé response of thé enclosure can consist in a burst or in thé split of thé vessel.

LPG bottles hâve security flanges that are designed to liberate thé bottle head for high pressures, and thus help reducing thé conséquences of thé rupture. Under thé effect of thé gas pressure, thé part sealing thé nose of thé tap distends, bursts and libérâtes thé gas. So thé pressure drops inside thé enclosure, causing more liquid to vaporize. Ail this gas ignites and provokes thé éjection of thé whole set: flange, tap and cap.

  constants ( O < n < 1 / 0,05 < fi < 0.15 / a > 0) cr ys represents thé yield strength of the material, a y thé flow stress and p y pressure corresponding to a fully "yielded cylinder": hère, parameters characterizing yielding are included in formulas. So thé dynamic fracture criterion can be written as:

			C \(P~ Pr) f Or P > Pr
	df2		P r )-C 2 (P-Pr) 2	M	P<Pr
	or	^T = \	C ^P-P^		f°r P>Pr	(12)
		dt2	\c3(p 2 -p2r)-c 4 (p 2 -p 2 r) 2	for	p< pr
	with c n experimentally determined constants.
	where		p y 'y	excr	a y =o ys (1

r with n, p ,a