CFD-Based Investigation of the Response of Mechanical Ventilation in the Case of Tunnel-Fire

The fact that many underground facilities, as large-scale spaces, differ in their geometric characteristics, continually sets new tasks for the researcher investigating the fluid phenomena. In this manner, we performed a computationally aided investigation of the new covered traffic road-communication that goes under a part of the Slovenian capital of Ljubljana. Although not yet finished – we were in position with the employed RANS transient approach (interpreting turbulence by k-ε model) using the technical drawings of the road-tunnel for the design of the computational domain, to give a CFD-forecast of the effectiveness of ventilation in the event of both, 40 and 80 MW fires.


APPLIED CFD-BASED APPROACHES THUS FAR
Due to an increasing number of accidental fires in enclosed spaces, especially in the area of modern society where our freedom is mostly expressed -traffic and tourism -there is an on-going need to undertake scientific research, aiming at better understanding these reactive flow phenomena.One of the direct benefits of such research is a suggestion for the implementation of the most effective fire suppression methods that should follow world-wide statistical data of one traffic accident occurring every 107 kilometres of covered (tunnel) road [1] and [2].
The field models for numerical research of fires in enclosed spaces (which was our applied mode as well) are established in the CFD-community over the zone-model-approach [3] and [4].Such codes (field model based) do not divide the area of interest into very few smaller control volumes only, but are based on the full solution of the fundamental physical laws of conservation.Here, the computational domain is divided into thousands of much smaller control volumes 'cells', where, after their discretisation mathematical models are "translated" into a programme code for combustion, radiation, turbulence etc.Therefore, between the mid 1980s and mid 1990s, hardware sources supported computational domains with few thousands cells [5] and [6].However, in spite of these limiting conditions highly satisfying results were accomplished in attempts of both, the validation [7] to [11] of software tools and aimed CFD-prognoses for particular explored cases of confined fluid phenomena [1] and [12].
Modern field-model codes that have been engaged in CFD-research within the last decade supported by powerful hardware can cope with the domains made out of several hundred thousands cells.All of these research works attempts that have been brought up into the CFD-community, report on good capability of the numerical approaches used in handling the reactive flows in straight enclosed traffic objects that were validated against some distinguished experiments [13].
However, in addition to the slight denivelation of a few percent, the geometry of the arbitrary tunnels was relatively a simple one [14].Therefore, in our research, we explored the "Sentvid" tunnel that goes under the city of Ljubljana and is being changed in its cross-section surface (from three-lane down to twolane road) [15].At about half of its length of 1470 m it has the bifurcation zone, where one-lane road (tunnelpart) climbs up to the local road.

TREATMENT OF THE TURBULENCE -APPLIED MATHEMATICAL MODEL IN THIS STUDY
In this study the flow phenomena are computed by the Reynolds Averaged Navier-Stokes (RANS) equations, with the turbulence k-ε model [16], representing the major characteristic of the applied CFD-investigationtool of the FLUENT; and handling the buoyancy by applying the Boussinesq approximation.This approach is not affected by the fluctuation of the initial conditions, offering a more accurate presentation of the time dependent flow -particularly the distribution of the gaseous combustion products [17].Since the investigations have shown that the maximum value for the Mach Number is of the order of 0.035, such a flow can be assumed as incompressible [18].Therefore, while crossing the reaction front, the fluid does not undergo thermal-caused expansion and the reaction makes no impact onto flow-velocity.A further assumption, to have a planar propagation front of combustion in a motionless fluid, leads to the application of the Boussinesq approximation [19] without external forces [20].Here, the flow velocity obeys the incompressible Navier-Stokes equation with a temperature-dependent force term [20].
The change of temperature is described by an advection-reaction-diffusion equation.For this incompressible gaseous reactive flow at low velocity, the governing equations of the combustion-induced flow read: Here v i denotes the average velocity component, T the mean local temperature, p the pressure, ρ the density, t the time and x i the space coordinates.The R(T) = ¼•T•(1 -T) stands for reaction rate [20] where the reciprocal value of reaction time-scale is represented by z, λ is the thermal conductivity, c p is the heat capacity at the constant pressure.Temperature T will also be used as the expression for reactionprogress-variable whose purpose is to distinguish the burned, the unburned and the partially burned state, providing an easy interpretation of flame propagation.The term -g i αΔT denotes buoyancy treated according to the Boussinesq approximation, where ΔT shows the difference between local and reference temperature.The symbol g denotes the gravity and α is the coefficient of thermal expansion.The model for the stress tensor [21], τ ij is related to the local strain rate: where we distinguish between the Newtonian stress (τ ij ) N = 2μS ij featuring molecular viscosity; and the turbulent Reynolds stress (τ ij ) T = 2μ T S ij , since the stress rate tensor S ij is defined as: and the turbulent viscosity: with k the turbulent kinetic energy and ε the dissipation rate of turbulent energy.
The applied k-ε model [22] is a two-equation eddy viscosity model [23] and [24] that uses transport equations for these two variables [25].One of these equations governs the distribution through the field of k, the local kinetic energy of the fluctuating motion.The other one yields the energy dissipation rate ε [26].
The energy term G g T

BOUNDARY CONDITIONS AND MESHING IN COMPUTATIONAL DOMAIN -A NUMERICAL APPROACH
For transient simulations (a CFD-mode that was applied in this study) the governing equations must be discretised in both space and time [24] and [27].In choosing the numerical method we rely on the standard of the finite volumes [18], [27] and [28].The spatial discretisation of time-dependent equations employed a segregated solution method.The linearised equations result in a system of linear equations for each cell in the computational domain, containing the unknown variable at the cell centre as well as the unknown values in the surrounding neighbour cells.
This mechanism for a scalar transport equation is also used to discretise the momentum equations; in the same mode for the pressure field (if face mass fluxes were known) and the velocity field will also be obtained in the same way.In case the pressure field and face mass fluxes are not known, FLUENT [24] uses a co-located scheme, whereby pressure and velocity are both stored at cell centres.The need for interfacial values includes an application of an interpolation scheme to compute pressure and velocity from cell values.The integration over the arbitrary volume (a cell in a computational domain) can be performed yielding the discretised equation for the mass flux through an arbitrary surface of a face.
In the sequential procedure of the segregated solver, the continuity equation is used as an equation for pressure as well [24].However, the pressure does not appear explicitly for incompressible flows since the density is not directly related to pressure.Therefore, the SIMPLE-family [27] of algorithms is used for introducing pressure into the continuity equation.In this way, it supports the pressure-velocity coupling though the algorithm for the unsteady flows as for the fluid mechanics phenomenon, explored in this study.The SIMPLE uses a relationship between velocity and pressure corrections to enforce mass conservation and to obtain the pressure field.
Since in our CFD-approach, the longitudinal artificial ventilation was employed immediately, the tunnel-entrance and exits (both towards the main tunnel-line and towards the local road) characterised as open (pressure) boundaries with pressure increase due to the ventilation-caused velocity (of 4 m/s).The fuel "pool" -the simulated fire-place, has been determined by the constant max flux rate of heptane of the order 0.4545 kg/m²s [29] for simulating the fire thermal power of 40 MW over the surface of 2 m²; and for the 80 MW fire [30] the double mentioned value for the mass flux.Heptane was taken as an inflammable good as one of the most commonly used fuels in experiments and fire-tests.
The tunnel itself (the left tube), with the slope of 2.2%, has two major main-curves with a Radius of 4000 and 1500 m with both "horse-shoe" crosssection and rectangular cross-section, where the first one determines the first 1080 m of the tunnel.Starting as a three-lane tunnel, after the bifurcation (on the 720 th meter of its length) the "main stream" of the "horse-shoe-shaped" tunnel "continues" as a twolane traffic communication; the "horse-shoe-shaped" tunnel line exites and elevates to the point of about 12 m above the road-level of the main tunnel-stream, with the length of a further 400 m and a change from a one-lane to the two-lane "horse-shoe-shaped" cross-section.Therefore, the Aspect-Ratio changed in this tunnel: for three-lanes part Ap = 1.707 to the rectangular and other horse-shoe-shaped part with two lanes with Ap = 1.32.
Tunnel-entrance as open (pressure) boundary was used for initializing computational values for the velocity and pressure in the domain since the global temperature was set to the 293 K.
The tunnel housing, tunnel road and tunnel walls, were presumed to be heat transparent.This decision was based on previous research [31] and [32], using the particular thermal conductivity of a rock, where a tunnel was built through [33], the reality-oriented investigation on tunnel-construction was conducted.Particularly for the Sentvid tunnel that was built in the Karst area, the specific thermal conductivity of such limestone (λ = 2.3 W/mK) [33] and [34] was integrated in the boundary conditions.

FIRST STEP IN INVESTIGATION -THE SIMULATION OF 40 MW-FIRE WITH NATURAL VENTILATION
After the explained computational approach was validated [35], a CFD-investigation of the 40 MW fire accident in tunnel natural air-movement only was performed.This numerical experiment ran for 120 s, which is the time of real-case physical experiments [35].The findings of this investigative step assured us that after accidental fire is established (7 th second) -there is a strongly expressed propagation (113 th second) of the hot gaseous combustion products towards the actual traffic-entrance (higher geodetic position).The accompanied occurrences are thermal damaging, not only to the re-enforced concrete (concrete spalling) of the tunnel cavity, but the also destruction of the asphalt road-layer.Those results point to the extremely important task for the longitudinal (forced) ventilation that should not only induce evacuation of the smoke towards the positive tunnel direction (in the sense of traffic) but should protect the tunnel body and road-elements from thermal damages before a fire-event is localised and finally extinguished.We started this numerical experiment with an assumption that the accidental fire (both of 40 and of 80 MW) was already established and the longitudinal ventilation in the left-tube (composed out of 7 ventilator-pairs) was set to meet the criteria of the "critical velocity" -the velocity in the ventilating process of a tunnel-fire-event that is just enough to stop the propagation of the smoke in unwanted direction.
We estimated that for this covered road facility the value of the critical velocity is 4 m/s.Indeed, in the area of the fire-accident, the overall velocity gains value due to the constructive interference of buoyant forces in fluid (caused by the combustion process) and the longitudinal air movement (caused by forced ventilation).Fig. 5.In 113 th second fire is "leaving" the actual source and temperature is pointing to the inflaming of the road-asphalt; concrete spalling is inevitable; the smoke is moving towards the entrance However, the applied ventilation-velocity (due to the 16 Pa pressure-difference) in the 22 nd second of the fire-accident shows the characteristics of "backlayering" -the smoke movement along the tunnelceiling against the applied ventilation.The temperature profiles pass the concrete spalling-temperature, which might start in about the 10 th minute of a thermal load [36].A little further from the fire place is the connection of the two kinds of tunnel cross-section-parts (rectangular and "horse-shoe-shaped").The temperature iso-surfaces showing destructible thermal load in the vicinity of the tunnel-ceiling (above the fire-place).during the CFD-investigations, we recognised the latter ones [37]; here, the 4 m/s ventilation-velocity in the 88 th second of the 40 MW fire provides protection for the walls of the tunnel cavity; still, there is a minor "back-layering" to be seen on the ceiling, above the fire-place only In a further extraction of the results of our numerical experiment, we noticed that in the 88 th second during the 40 MW fire, the longitudinal ventilation was as capable (in providing the lower temperatures on the tunnel cavity) as in the 120 th second of the same fire-event.This implies that artificial longitudinal ventilation adjusted for 16 Pa of pressure difference can cope with the fire event up to 40 MW of the fire´s thermal-power -and this already after a minute-and-the-half.With this engagement, the possible RC-spalling is inhibited.However, the thermal impact (due to the developed temperature fields and due to the heat irradiance) destroys the asphalt road-layer that is, in this case, a cheaper part of this underground space facility, in terms of reconstruction.
The additional size of the temperature impact offers a view in 112 s, at the tunnel-ceiling and tunnelwalls next to the 80 MW fire-place (Fig. 9).The zones with the developed temperature of 500 K (the begin of the RC-spalling) "resist" the longitudinal ventilation, that in the 40 MW fire event was able to cope the accident.In this stage of the fire-development these occurrences present a dangerous part in the construction of the tunnel body (Fig. 9).And not even in the 126 th second, during the developing of the 80 MW fire-consequences, the applied longitudinal artificial ventilation with the velocity of 4 m/s was not able to perform a completely successful protection of the tunnel-walls and the tunnel-ceiling (Fig. 10); but the areas with the temperature of 500K that stand for the beginning of the RC thermal destruction can still be seen (Fig. 10).
The result we gained on the effectiveness of the proposed ventilation magnitude in case of 80 MW fire, encourages a further numerical experiment in order to estimate the optimal velocity of the ventilation (in case of 80 MW fire).It will offer appropriate thermal protection until accidental fire is fully under control.
However, the chosen velocity of the ventilation was able to suppress propagation of the hot gaseous combustion products in unwanted direction in both accidental cases.
By performing this study on only one object of interest of the given geometric characteristics, the CFD-based investigation on the accidental fire in the artificially ventilated tunnel "Sentvid" was conducted according to both standards and novel experimental [38] to [41] and computer aided [42] to [44] research [45] and [46].
With the results of this research-attempt, we intend to also address the civil-engineering sector [46] and [48] and enlarge data-base for the medical healthprotection [49].
The specific geometry of the object of interest -this traffic road-object built in Slovenia -was a "fine provocation" to conduct this research, expecting new answers due to the possible impact of a realityoriented enclosure (computational domain) onto large-scale fire and escorting occurrences.
The "chimney effect" of the exit up to the local road was not "strong enough" and its possible influence was inhibited due to applied longitudinal ventilation that forced gaseous combustion products towards traffic-positive direction, already after the first minute-and-a-half of numerical investigations.Therefore, the propagation of the gaseous products "followed" the tunnel line (towards a lower geodetic position) and in the cases of 40 and 80 MW-fire, the expected major influence of the tunnel-road exit was not noticed.

ACKNOWLEDGEMENT
Authors do thank to Mr. Marko Zibert, Dipl.Ing. and to Ass.Prof. Jurij Kazimir Modic, PhD. for decisive support in performing this study.Special gratitude goes to Prof. Karel Ciahotny, PhD. as well as to Prof. Vazlav Koza, PhD. from the Institute of Chemical Technology, Prague for immense software-support during this CFD-based research.

Fig. 1 .Fig. 2
Fig. 1.Left tube of "Sentvid" as computational domain the 1470 m long left tube with exit to the local road that runs above the tunnel; on the upper side there is the entrance from Ljubljana (from the south)

Fig. 3 .Fig. 4 .
Fig. 3.The 7 th second of the 40 MW fire -velocity fields point at not established propagation of the gases

Fig. 6 .
Fig. 6.The view of the 40 MW-fire-place towards the tunnel bifurcation-zone (deep in the sketch); upper left, is the planned exit to the local road

Fig. 7 .Fig. 8 :
Fig. 7.A side view of the 22 nd second of the longitudinal ventilation in case of 40 MW fire event; on the iso-surfaces that represent ventilation velocity of 4 m/s the temperature images are drown

Fig. 9 .Fig. 10 :
Fig. 9.The additional size of the temperature impact offers a view in 112 s, at the tunnel-ceiling and tunnel-walls next to the 80 MW fire-place Pr t denotes turbulent Prandtl Number (which is of the order of unity).The constants are given: C 1 = 0.126, C 2 = 1.92,C µ = 0.09, C e = 0.07.