Pergamon
M©dlanics Research Communications, VoL 21, No. 6, pp. 635-644, 1994 Copyright © 1994 Elsevier Science Ltd Printed in the USA. All rights reserved 0093-6413/94 $6.00 + .00
0093-6413(94)00030-1
COOLING OF A GASEOUS NOZZI~E FLOW BY INJECTION OF LIQUID DROPLETS E. DANIEL, J.C. LORAUD, M. LARINI, R.SAUREL IUSTI/SETT, Equipe Ecoulements Diphasiques et R6actifs, UA-CNRS 1168 Universit6 de Provence, Centre de St-J6r6me, case 321 13397 Marseille France
(Received 6 October 1993; accepted for print 3 August 1994)
INTRODUCTION
A numerical simulation of two-phase flows inside a nozzle is presented. This kind of flow has been intensively studied in the last decade, especially in aerospace research and more particularly in propulsion. The study of such flows was done by Chang in two-dimensional 1 and three-dimensional 2 configurations. In these works, the flow field equations are written with an Eulerian-Eulerian approach and solved by using an explicit finite-difference scheme: each phase is considered as a continuum 3. An alternative approach consists in treating the dispersed phase as a discrete one: it is an Eulerian-Lagrangian model, where the gas phase is usually solved by the means of a finite-volume technique. The works of Golafshani et al. and Sabnis et al. provide good examples of the use of such methods4,5,6. The first approach is used in this paper. The flow fields we are induced to study in this paper concern two-phase flows constituted of hot gas and water droplets. The gas is produced by the combustion of a pyrotechnic solid and the water droplets are injected in order to cool the fluid. We assume that the two-phase feature is due to the presence of water droplets in the flow. The presence of solid particles resulting from the combustion of the pyrotechnic composition is neglected. The particularity of the presented work lies in the fact that the inlet flow is a one-phase and the dispersed phase is injected into this flow field: the two-phase flow is the result of the injection of the water droplets. So there are two kinds of flow: a one-phase flow upstream of the injection then a two-phase one. The works done by Chang are of a great interest for our purpose because they allowed to do a qualitative comparison between his previous results and ours, in the case of two-phase flow in a nozzle. In this paper we present the classical governing two-phase flow equations in an unsteady, two dimensional planar coordinates system. These equations take into account interphase transfers: drag force, heat and mass exchange. The system is solved by using the Mac-Cormack finite difference scheme. Results related to two-phase flow in a nozzle are qualitatively compared with previous works computing steady solution. Then the results of the computation of strongly transient process of water droplets injection are presented.
GOVERNING EQUATIONS
The two-fluid approach yields a system composed with six vectorial unsteady equations written in two-dimensionnal planar coordinates system. Each phase is described by three equations which are solved by using the same time-integration procedure. This system is obtained by a conservative balance analysis and was former presented by Olim et al 7 635
636
E. DANIEL, J.C. LORAUD, M. LARINI and R. SAUREL
Gas phase
3 3
~7 oc (p ~. F Hal+
P~.U)= - Ud'Fd-q-
3t
(1)
Disoersed phase
~- Pd C~d+ V . Pd ~d Ud = -F
0
03
~ ~ 0d ~.d Ud + V . ad (0d UdUd )--- F Ud + Fd- ~d VP (2)
l0 ~ - i - n + V . n U d - - ( I where ~ l = l - a . These equations are written in a weak conservative form: Ut + Fx + Gy + H=0
(3)
The subscripts d and g indicate respectively a droplet phase and a gas phase quantity, o~ is the void fraction and represents the ratio of the volume of the gas to the total volume, p is the density, P the gas pressure, u et v are the components of the velocity vectors, ~ the total energy per unit mass. There is no energy equation for the dispersed phase: this equation is in tact reduced to the expression of the droplet temperature which is supposed to be uniform and known. This is the saturation temperature, Tsat, which depends on the gas pressure. So the heat flux from the gas to a droplet is totally used to the vaporisation process. The last equation of the system (2) expresses the fact that any break-up or coalescence of the droplets is ignored. The gas is assumed to be inviscid and obeys the ideal gas law. The right-hand side of the systems (1) and (2) contains the exchange terms: drag force ( "F' d), heat transfer (q) and mass transfer (F). The expression of these terms can be easily found in the literature 8.9. The mass transfer is due to the droplet vaporization, so F it can be deduced from the relation: F = Cl with Lv the specific latent heat vaporisation. We can notice Lv the presence of the term c~d V P in the right-hand side of the motion equation of the dispersed phase and ot V P in that of gas phase. This last term is split in two terms • one is in the flux vectors (F and G) and the other in te~xn source vector H. This gradient is often neglectedl, 2 or is taken into account only in the gas phase motion equation 4,5. This is done without real
GASEOUS NOZZLEFLOWCOOLING
637
justification neither on the role nor on its value by comparison with the other transfer terms. It is the reason why we keep the gradient in the presented equations. NUMERICAL PROCEDURE
The Mac-Cormack finite difference scheme is applied to the partial differential equations in a computational domain (x,q). The transformed equations take the following form:
(4) with E=E/D, F'--(~t E + ~x F + ~y G ) / D =(rh E + fix F + fly G)/D
~=H/D where D is the Jacobian of the transformation. Applying the algorithm for equation (4) yields to a predictor-corrector scheme: predictor: ~n _ Gi.j&-Gij ~n ~n I ~* = ~n + At ~ n . ~n Fiq.j -Fi.j
a~
An
(5)
l
Corrector: =
+ Uj
+At Hj
(6)
The physical domain of integration is the entire nozzle. Symmetry condition about nozzle median plane cannot be used : the injection of the dispersed phase implies a non symmetrical flow. In this physical domain the lower wall of the nozzle can be represented by the function RL(x) and the upper wall by RU(x). This last relation is deduced from RL(x) by a symmetric condition about the nozzle median plane. The coordinate transformation is given by the following relation: y - RL(x) ~(x,y)=x/L and l](x,y)= RU(-x)--~(-x) (7)
BOUNDARY CONDITIONS
At the wall Computation of variables along the body surface is a difficult problem because if the gas phase is inviscid then only one condition is given, namely ~7. "~ =0. The treatment of boundary conditions is different for the gas phase and for the dispersed one. The droplets boundary
638
E. DANIEL, J.C. LORAUD, M. LARINI and R. SAUREL
conditions are calculated from two different methods, depending on whether the case is steady or unsteady.
gas phase Conservative equations of system (3) are calculated with Mac-Cormack scheme. To succeed in this procedure, spatial derivatives are approximated by upwind second-order scheme 10. After each stage the velocity components is computed as below: v = - fix/fly u
dispersedphase The same procedure as for the gas phase is used. But if the upwind second-order scheme is applied, the solution can not be reached because of instabilities at the wall. Then derivatives are evaluated with a first order one-sided formula: Inflow and outflow boundarv conditions
Gasphase At inflow, the flow is assumed to be adiabatic and isentropic. Compatibility equations along the direction ug-c, where c is the sound velocity in the gas, are used to obtain the velocity component u, while the time component v, is taken equal to zero.
d
with d t = O--i-+ (ug-c)
}
dP dt
pc ~-~, u =0
(8)
for the equation along the dil'ectian Ug - c.
At the exit, if the outflow is subsonic, the ambient pressure is specified. Two compatibility equations along the direction ug and Ug+Care solved: dP dt
c 2dp = C ( y _ l ) _ c 2 ( r'':)~'~ A dt 7
d P _ c d P = C (y-l) + c B+ C2A dt dt y
(9)
(10)
The terms A, B, C represent respectively the source terms of the equation of mass, momentum and energy of the gas phase. If the flow is supersonic, conservative variables are extrapolated.
Dispersedphase At the inlet or at the location of injection, intbrmation such as velocity, diameter, density number is specified. At the outflow, the droplets are only allowed to exit the nozzle, then conservative variables are extrapolated.
GASEOUS NO77J .E FLOW COOLING
639
RESULTS
] .Introduction The main purpose of this work is the numerical simulation of water droplets injection in a one-phase flow. But, no experimental or numerical data have been found, so the model used should be applied in a known case. The more convenient and studied work, is obviously linked to two-phase flow in nozzle, where the inlet flow is already a two-phase one (by contrast to an injection). But the comparison can be only qualitative for several reasons: -Firstly, the governing equations are not quite the same as generally used: the energy equation of the dispersed phase is not used here (because we study saturated liquid doplet) but the droplet number conservation equation is solved. -The initial one-phase flow is subsonic. -Finally, the calculations are performed with a void fraction ct very closed to unity. Previous works have shown that downstream from the geometrical throat a free droplet zone exists 1,2,4. This zone is due to inertia phenomena: the droplet trajectory cannot curve and this area is all the more important than the size of droplets is large. It is interesting to verify if that zone exists even if the droplet number per unit volume is very low. 2.Two-phase flow in a nozzle The aim of this calculation is to reach the steady state of a two-phase flow calculation. The initial data for the gas phase are the numerical solution for inviscid, subsonic, one-phase flow. At the throat vicinity the Mach number is about (}.65. Stagnation conditions of the flow field are taken as: Pres = 1(}34 kPa and Tres= 555 K and the ratio of specific heat 7=1.4. The initial values of the dispersed phase ale arbitrarily fixed at every location of the grid: Ud= ~v Ug with kv=l; the droplet number per unit volume is taken equal to 109 droplets/m3. Presented results correspond to two different sizes of the droplets: D0=51.tm and D0=501.tm. It yields to the following values of the void fraction: or=l-o,= (}.65 1(}-7 for D0=5~tm and O.d=l-tx= 0.65 10 -4 for D0=501am. The density of the dispersed phase is constant and taken as 4000 kg/s. The length of the nozzle is 11.5 cm with a 30 ° entrance angle and a 15 ° exit angle. The radius at inflow is 4 c m , 2 cm at the exit section and at the throat lcm. The grid is built with 75"31 points and the integration time step, At, is equal to 0.5 CFL. The steady state is reached when the maximum variation of each conservative variables between two time steps does not exceed 10-2% for the gas phase and 1(}-1% for the dispersed phase. A typical run takes about 1200 cycles and 58(1CPU seconds on an IBM-3090. Contours of the normalised droplet number per unit volume at the steady state are plotted in figure 1, related to a droplet size D0=5lam.
640
E. DANIEL, J.C. LORAUD, M. LARINI and R. SAUREL
Free
droplet
zone
Figure 1: Steady calculation of the nozzle flow: Normalised density number contours for D0=51am In figure 2, the same contours are plotted, related to droplet size D0=501am. We can note in the two cases that the free droplet zone exists even if the dispersed phase is very dilute.
Figure 2: Steady calculation of the nozzle flow: Normalised density number contours for D0=50tam
The larger the radius of the droplets is, the greater are the inertia effects. For this reason the two-phase jet is larger in the exit cross section when D0 = 5() lam than when DO = 5 lam. A transversal view of the variable (n/n0) is presented in the figure 3. This transverse section is located at 9.3cm from the entrance.
GASEOUS N07'7.I E FLOW COOLING
641
o.
x
~ m
o -
lnLtLoL
x - Steod~
f o.o
5tote 5tote
J
~.s
l.o
;.s
~.o
Figure 3: Cross section of n/n0 at the location x=9.3cm from the entrance 3, In iection of dronlets The injection of droplets in a nozzle is presented now. The nozzle is symmetric with a 15° entrance and exit angle. The length is 0.1 m, the radius at the inflow and at the outflow is 2 cm, and the throat radius is 1.3cm. To solve the equations a second order damping term is necessary to avoid non-linear instabilities, with a coefficient equal to 10-3. The calculation is unsteady that it means the results obtained after each time step integration are supposed to be exact. A one-phase inviscid subsonic solution provide data for the initial gas phase. The stagnation conditions are Pres= 4.5 bar and Tres=900K. Thermophysic quantities of the gas remain constant: y----4/3,Cp= 1032 kJ/K/kg, la~l. 1310 -4 kg/rrds and L=2.410-2 W/m/K. At time t---0 s there is no droplet in the nozzle. Injection occurs at a part of the lower wall in the divergent. The diameter of injected droplets is equal to 5(l~tm. The rate of mass of this phase is proportional to the gas one: Wd=0.3Wg which yields Wd=2 kg/s. If the area of the surface S through which the rate of mass of liquid is injected is known, then the normal velocity of injected phase can be computed by: Ud = Wd / S/Pd where Pd is the density of the liquid water taken equal to 1000 kg/s. The droplet number per unit volume can be extracted from the relation: Wd= [no Pd ~ I~0] U d S . It yields to n0=0.1 1012 droplets/m 3. Then the dispersed phase is completely defined. The boundary between the two-phase and one-phase domains requires the use of a numerical cut-off. It is assumed that the flow remains one-phase until n(x,y,t) is less than a chosen value nlim; nlim must be much smaller than no. Various values of nlim have been tested and it has been found that nlim must belong to the range 1(13-106 droplets/m 3. This methodology is only required when the mass transfer between the two phases is taken into
642
E. DANIEL, J.C. LORAUD, M. LARINI and R. SAUREL
account. Indeed, in this case, the radius of droplets must be computed at every time step and every location in the grid. The radius is given by the following relation
3 {1-c~) 4rtn Then when n(x,y,t) equals zero, R is undefined. Besides, the scheme used is second order in : R 3 =
space and time accurate, then some oscillations may appear, and also some over- or undershooting values, which may be leather than zero. So, R can not be correctly evaluated and it is the reason why we used a cut-off. In figure 4, the iso-velocity contours of droplets are plotted at the steady state. We can see that the two-phase jet curves and spreads due to the drag force between gas and droplets. The figure also allows the evaluation of the jet size.
I
Figure 4: Iniection of droplets: Enlargement view of droplets isovelocity at the steady state The droplets are injected in the flow in order to cool the gas phase. In figure 5 are plotted the temperature isolines nozzle flows with and without droplets. The major variation of the gas temperature occurs in the area where the flow is a two-phase one. With the temperature contours, we can easily identify the droplets jets.
NLthout DropLet
SLeody SLot.e Figure 5: Injection of droplets: Gas temperature contours at t=0s (without droplet) and at the steady state.
GASEOUS NOZZLE FLOW COOLING
643
In figure 6 is drawn the temperature evolution along the axis of the nozzle at several time steps. In the exit section the change in temperature reaches about 9(}K, and this occurs at a short distance from the beginning of the injection equal to 2.5 cm. In this figure we can also observe a diminution of about 10K at the geometrical throat. This is due to wave moving backward from injection to entrance section l0 (with the characteristic speed equal to ug-c, here negative). These perturbations go backward. Then, the inflow condition may be changed due to the interior domain perturbations and must be calculated by an unsteady method as the characteristic method used in this paper. o.
o.
il. | 0.000
o.~
o.~
o.~
Oi.stonce t'rom InLet lml
o.'~oo
Figure 6: Iniection of droplets: Evolution of gas temperature along the median plane of the nozzle at several times CONCLUSION A numerical simulation of a dispersed phase in an inviscid nozzle flow has been achieved, using an Eulerian-Eulerian technique. The results obtained could have been compared with some results issued from an Eulerian-Lagrangian approach4,5. The injection of droplets at saturated temperature can cool the flow on a short distance. We have shown that some modification can occur upstream the injection, up to the inflow and then modify the flow conditions at the entrance section. The numerical treatment of boundary conditions able to take into account such variations must be used.
REFERENCES [1] I.Shih Chang, One and two-phase nozzle flows, AIAA Journal, vol 18, pp1455-1461 (1980).
644
E. DANIEL, J.C. LORAUD, M. LARINI and R. SAUREL
[2] I.Shih Chang, three-dimensional two-phase supersonic nozzle flows,AlAA Journal. vol 21, n°5, pp671-678 (1983). [3] Mehdi Golafshani and Hai-Tien Loh, Computation of Two-Phase Viscous Flow in Solid Rocket Motors Using a Flux-Split Eulerian-Lagrangian Technique, AIAA/ASME/ASEE, 25th Joint Propulsion Conference, Monterey, Ca (1989). [4] E. Daniel, J.C. Loraud, M. Larini, Influence de l'injection de gouttes d'eau dans de la vapeur en ecoulement en tuy~re, International Journal of Heat and Mass Transfer, Vol 36, n°6, pp 1619-1632 (1993) [5] R. K. Madabhushi, J.S. Sabnis, F.J. de Jong, H.J.Gibeling, Navier-Stosckes analysis of aft dome flow field in solid rocket motors with submerged nozzle, AIAA-89-2780 [6] C.K. Hu and S.Z. Yu, Numerical Calculation of Gas-Liquid Two-Phase Flow in a Dump Combustor, 7th International Conference on Numerical Method in Laminar and Turbulent Flow, Stanford (1991). [7] M.Olim, O.Igra, G.Ben-Dor, A general attenuation law of planar shock waves propagating into dusty gases, Proceedings of the 16th International Symposium of Shock tube, pp 217-225 (1987) [8] M.C Yuen et L.W Chen, Heat transfer measurements of evaporating liquid droplets, International Journal of Heat and Mass Transfer, Vol 21, pp 537-542 (1977) [9] J.N Chung et S.I Olafsson, Two-Phase droplet flow convective and radiative heat transfer, International Journal of Heat and Mass Transfer, Vol 27, pp 901-910 (1983) [10] R.F. Warming et R.M.Beam, Upwind second order schemes and applications in aerodynamic flows, AIAA Journal, vol 14, n°9 (1976)