- Email: [email protected]
Adaptive block-implicit methods for annular liquid jets
Adaptive block-implicit methods for annular liquid jets J. I. Ramos Department
of Mechanical
USA, and Facultad Spain
Engineering,
Carnegie Mellon
de Informcitica1E.T.S.I.
University, Pittsburgh, PA, Telecomunicacibn, Universidad de Mdlaga,
A domain-adaptive technique is used to analyze the dynamic response of inviscid, isothermal, annular liquidjets to velocityfluctuutions at the nozzle exit. The adaptive technique maps the unknown, timedependent, curvilinear geometry of the annular jet into a unit interval and yields a system of integrodifferential equations for the muss per unit length, radius, and axial and radiul velocity components of the liquidjet. The convergence length is governed by an ordinary differential equation that depends on the values of the dependent variables at the convergence point. The governing equations are solved by means of two methods of lines, which discretize the sputial coordinate but keep continuous the time, and a Newton method. A block-bidiagonal technique is also used to solve the fluid dynamics equations. It is shown that the dynamic response of liquid jets to sinusoidal velocity jluctuutions at the nozzle exit is periodic and nearly sinusoidal for small amplitudes of the velocity oscillations at the nozzle exit. For large amplitudes the fluctuations of both the convergence length and the pressure of the gases enclosed by the annulur jet are periodic but not sinusoidal. It is ulso shown thut there is a delay time between the velocity fluctuations ut the nozzle exit und those of the convergence length. This delay time is nearly independent of the amplitude of the velocity fluctuutions ut the nozzle exit and decreases as the Weber and Froude numbers and the nozzle exit angle ure decreased. Keywords:
1.
block-implicit techniques, methods of lines, annular liquid jets, block-bidiagonal
Introduction
Annular liquid jets have been the subject of intensive research in recent years, ‘-I3 because under certain conditions they can form enclosed volumes, which can be used to determine the dynamic surface tension of liquids and to burn toxic wastes. The equations that govern the fluid dynamics of annular liquid jets have been derived by the author for Newtonian and non-Newtonian fluids’4*‘5 by using Cauchy’s equations and applying interface boundary conditions. The interfaces between the liquid and the fluids surrounding and enclosed by the annular liquid jet are material surfaces where the shear stress is continuous, while the jump in normal stresses across the interfaces is balanced out by surface tension. The author derived the asymptotic equations that govern the fluid dynamics of inviscid, isothermal, annular liquid jets (Figure 1) by assuming that the jet thickness is smaller than the jet radius at the nozzle exit, averaging the governing equations across the anAddress reprint requests to Dr. Ramos at the Department of Mechanical Engineering, Carnegie Mellon University, Pittsburgh, PA 15213-3890, USA. Received 14 March ruary 1992
464
1991; revised
6 February
Appl. Math. Modelling,
1992; accepted
14 Feb-
1992, Vol. 16, September
method
nular jet, and applying the interface boundary conditions at the jet’s mean radius.‘4.‘5 The resulting set of asymptotic equations is also valid for thin, viscous annular jets except near the nozzle exit, where the velocity profile relaxes from stick boundary conditions at the nozzle walls to slip/free boundary conditions at the annular jet’s interfaces, and near the convergence point, where the radius of the annular liquid jet’s inner interface is zero, a meniscus is formed, and an axial pressure gradient may exist. The equations derived in Refs. 14 and 15 are exact for inviscid, isothermal, annular liquid membranes, i.e., annular liquid jets of zero thickness, and have analytical solutions for steady-state conditions if both the slope and the curvature of the jets are small.‘-4 For time-dependent problems, however, analytical solutions are difficult, if not impossible, to obtain, and numerical techniques must be used. In this paper the dynamic response of inviscid, isothermal, annular liquid jets to mass flow rate fluctuations at the nozzle exit is analyzed by means of three adaptive, block-implicit techniques. The domain-adaptive techniques and the numerical studies presented in this paper were motivated by some experimental results,16 which indicate that under certain conditions the convergence length of annular liquid jets may exhibit a periodic behavior. This periodic behavior may be due to periodic fluctuations in the pressure of the 0 1992 Butterworth-Heinemann
Adaptive
Figure 1.
Schematic
methods
for annular liquid jets: J. I. Ramos
terns, i.e., nonlinear dynamical systems governed by ordinary differential equations. It must be indicated that annular liquid jets absorb the gases that they enclose at a rate that depends on the pressure difference across the jet, gravitational acceleration, surface tension, annular jet’s thickness, nozzle exit angle, and diffusivity of the gases in the liquid.2.3 Therefore in order to maintain steady state, mass must be injected into the volume enclosed by the annular jet at a rate equal to the mass absorption rate by the liquid. If the mass injection rate exceeds the mass absorption rate, mass will accumulate in the volume enclosed by the liquid jet, and the pressure of the enclosed gases will increase until the annular jet becomes unstable, breaks up, and does not form an enclosed volume. On the other hand, if the mass absorption rate by the liquid is greater than the mass injection rate, the annular jet will collapse. Both the growth and the collapse of annular liquid jets are timedependent phenomena which require that the gas concentration in the liquid be determined together with the solution of the equations that govern the unsteady fluid dynamics of annular liquid jets. In this paper it is assumed that the liquid does not absorb the gases that it encloses and that mass is not injected into the volume enclosed by the annular jet. Analytical studies of annular liquid jets indicate that there is a critical value of the pressure of the gases enclosed by annular liquid jets beyond which the jet cannot form a closed volume.‘.‘.4 These analytical results were obtained by assuming steady-state conditions and annular jets of small slope and curvature. Under transient conditions it may be possible to exceed this critical pressure while still obtaining a stable, yet fluctuating, annular liquid jet.5.‘7.‘x Numerical studies of annular liquid jets subject to fluctuations in the pressure of the gases that surround the jet’ conf-irm that this is indeed possible, provided that the annular liquid jet is subjected to high-frequency unsteady pressure fluctuations. The domain-adaptive technique presented in this paper maps the unknown, time-dependent geometry of the liquid jet into a unit interval and yields a system of integrodifferential equations for the convergence length, i.e., the axial distance at which the annular liquid jet merges on the symmetry axis to become a solid jet, and the annular liquid jet’s mean radius, mass per unit length, and axial and radial velocity components. Two of the block-implicit methods presented in the paper use a method of lines technique, which discretizes the spatial derivatives, keeps continuous the time, and results in a system of nonlinear, first-order, ordinary differential equations in time that was solved by means of the Newton method. The third blockimplicit method discretizes both the spatial and the temporal derivatives and yields a block-bidiagonal matrix, which can easily be solved by forward substitution. Furthermore, one of the block-implicit methods uses a conservative, finite volume discretization of the nondimensional equations written in weak conserva-
of an annular liquid jet
gases that surround the jet, unsteady mass flow rate fluctuations at the nozzle exit, and/or fluctuations in the mass injected into the volume enclosed by the jet. The effects of fluctuations in the pressure of the gases that surround the jet have been investigated by Ramos,j who used an adaptive, block-bidiagonal technique and showed that the convergence length and the pressure coefficient are periodic functions with the same period as that of the pressure fluctuations. However, the numerical calculations presented in Ref. 4 are academic in nature, for they were based on the assumption that the pressure difference between the gases enclosed by and surrounding the liquid jet was a known function of time, and they did not account for the variations of the volume enclosed by the liquid. Furthermore, the nondimensionalization of the governing equations presented in Ref. 4 can only be used to study the dynamics of annular liquid jets in gravitational environments, because the acceleration of gravity was used to nondimensionalize the time and the radial and axial velocity components of the liquid jet. In this paper a nondimensionalization of the governing equations that is valid for annular liquid jets in gravitational and nongravitational environments is introduced, and the pressure across the annular jet is calculated by accounting for the variations in the volume enclosed by the jet. The effects of fluctuations in the mass flow rate at the nozzle exit and in the mass injected per unit time into the volume enclosed by annular liquid jets have also been studied by the author,17.1X who used a blockbidiagonal method and a method of lines technique. The numerical studies presented in Refs. 5, 17, and 18 indicate that the unsteady dynamics of inviscid, isothermal, annular liquid jets subject to fluctuations in the pressure of the gases enclosed by the annular jet, fluctuations in the mass flow rate at the nozzle exit, and fluctuations in the mass injected per unit time into the volume enclosed by the annular jet has some similarities with the nonlinear dynamics of lumped sys-
Appl.
Math.
Modelling,
1992, Vol. 16, September
465
Adaptive methods for annular liquid jets: J. I. Ramos tion law form, whereas the other two methods are based on the discretization of the nondimensional equations written in nonconservation law form. The paper has been arranged as follows. Section 2 summarizes the equations governing the time-dependent fluid dynamics of inviscid, isothermal, annular liquid jets and the nondimensionalization used throughout the paper. Section 2 also identifies the nondimensional parameters that govern the dynamics of liquid jets. In section 3 a domain-adaptive method for the numerical solution of the nondimensional equations is presented. The three block-implicit methods used to discretize the governing equations are presented in section 4, whereas the initialization of the time-dependent calculations is described in section 5. The two sections, Presentation of Results and Conclusions, put an end to the paper. 2.
Fluid dynamics equations
The dimensional equations governing the fluid dynamics of inviscid, isothermal, incompressible, annular liquid jets can be written as’4*‘5
am* i3 t*+az*(m*u*)=O
(1)
-$ (m*u*) + -$ (m*u*‘) = m*g + 2flg
$(m*P)
+
+ R* $
(P,* - Pi*) (2)
-$ (m*u*P)
(4) where the stars denote dimensional quantities: t* is time; z* is the axial coordinate measured from the nozzle exit; m* is the annular jet mass per unit length and per radian; u* is the liquid axial velocity; 3 is the cross-sectional averaged radial velocity; g is the gravitational acceleration; (T is the surface tension; R* is the annular jet’s mean radius; p: and py* denote the pressure of the gases enclosed by and surrounding the annular jet, respectively; and .l*
=Xx/[1
+
(s)‘]“’
(5)
Equations (l)-(5) were obtained by averaging the Navier-Stokes equations across the annular jet thickness and by including the kinematic and stress conditions at the jet’s inner and outer interfaces, which are material surfaces. At each gas-liquid interface the shear stress is continuous, whereas the jump in normal stresses is balanced by surface tension.‘4v’5 Equations (l)-(3) represent the conservation of mass and the conservation of axial and radial momentum,
466
Appl. Math. Modelling,
1992, Vol. 16, September
respectively, whereas (4) corresponds to the kinematic condition applied at the annular jet’s mean radius. For annular liquid jets (Figure I),
R,* = R* + b*/2
RF = R* - b*/2
(6)
where R,* and RT denote the radii of the annular liquid jet’s inner and outer interfaces, respectively, and the annular jet’s mass per unit length and per radian is governed by the following expression: m* = p*R*b*
(7)
Equations (l)-(5) are asymptotic for inviscid, isothermal, annular liquid jets of thickness b* # 0 and are exact for inviscid annular liquid membranes, i.e., for b* = 0, and their numerical solution requires the specification of initial conditions at t* = 0 and boundary conditions at the nozzle exit Z* = 0. The initial and boundary conditions may be functions of z* and t*, respectively. In this paper the following boundary conditions were employed at the nozzle exit:
m*(t*,O) = rng
u*(t*,O)= uo*+ a* sin (w*t*) (8)
Z(t*,O) = u*(t*,O)tan $
R*(t*,O) = R;
(9)
where the subscript 0 denotes conditions at the nozzle exit; m$, R,*, ut, a*, w*, and 0, are time independent; 13,denotes the angle between the symmetry axis and the velocity vector at the nozzle exit; and a* and w* are the amplitude and angular frequency, respectively, of the axial velocity fluctuations at the nozzle exit. Note that since m*(t*,O)is constant, whereas u*(t*,O) is a periodic function of time, the mass flux m*(t*,O)u*(t*,O)at the nozzle exit is a periodic function of time. Note also that for annular liquid membranes. i.e., annular jets of zero thickness, RT = R,* = R*, b* = b,* = 0, and m* is the membrane’s mass per unit length and per radian. In order to maintain steady state, mass must be injected into the volume enclosed by the jet at a rate equal to the mass absorption rate by the liquid. If mass is not injected into the volume enclosed by the annular jet, this will collapse due to mass absorption. In this paper it is assumed that the Peclet number is infinity, so that the liquid does not absorb the gases enclosed by and surrounding the annular jet,2,3 and that mass is not injected into the volume enclosed by the jet. It is further assumed that the gases surrounding and enclosed by the liquid jet are ideal and isothermal and that the gases surrounding the liquid are infinite in extent, so that p: can be assumed to be constant. For ideal gases, pf = p7Rl-T
(10)
where l? is the specific gas constant and pi* and Ti” denote the density and temperature, respectively, of the gases enclosed by the annular jet. Furthermore, for low Mach numbers, ~7 is nearly uniform, and (10) can be integrated to yield
p”V* = rnTl?p
(11)
Adaptive
where L'
L' V*
m* =
mR,@&*
= 0
p,%-R,+‘dz*
(12)
for annular liquid jets: J. I. Ramos
velocity fluctuations at the nozzle exit to a characteristic residence time. Substitution of (14) and (15) into (11) yields
0
V* and m,? denote the volume and the mass of the gases enclosed by the annular jet, respectively, and L* is the convergence length, i.e., the axial distance at which the annular jet’s inner radius vanishes, R,*(t*,L*(t*)) = 0
P~=~(m;/~R~dz.)$=&
u*
M, = rR;‘p$lRT,+
uo*
* m=-
z=L
Ro*
uo*
L = L*IR;
mX
(14)
(13
RZ of (14) and (15) into (l)-(7)
Substitution
m*
t=t*!$
R=g
b=&
mi = rn~IA4~
(13)
v*
jj-
14 = -
(28)
where
Equations (l)-( 13) can be written in nondimensional form by introducing the following parameters:
Substitution
methods
of (14) and (15) into (24) yields
C,...,,[m./~Rfdz - 1) (30)
where Cpmax= p:R;l2a
yields
am a t+$mu)=O
C,,, =
(16)
(29)
(31)
It is clear from (16)-(19) and (30) that the fluid dynamics of annular liquid jets depends on m;, which in turn depends on the mass injected into the volume enclosed by the jet. In the absence of both mass injection into the volume enclosed by the jet and mass absorption by the liquid, mi = m,(O).
(17) (18) aR :=_+uE at
az
,=,/[I
+ ($)2]“2
(19)
(20)
3.
Domain-adaptive
The annular dependent downstream formed into
numerical
technique
liquid jet geometry is curvilinear and time and has an unknown, time-dependent, boundary. This geometry can be transa unit interval by means of the mapping,5
(t,z) + (7371)
r=t
r) = ZIL
OS~~l (32)
R, = R + hi2
Ri = R - b/2
(21)
where L(t) = L*(t*)IR,* is the Jacobian of the mapping.
(22)
3.1. Governing equations in weak conservation form Substitution of (32) into (16)-(19) yields
where Fr = u~~/~R~
We = m$&‘/2aR* 0
C,,, = C, We
C, = (p” - p~)R~Zlm~u,$2
(23) (24)
Fr and We are the Froude and Weber numbers, respectively, and C,, denotes the pressure coefficient. Note that Fr and We are constant and that C,, and C,,,, may be functions of t because p” may be a function of 1. Substitution of (14) and (15) into (8) and (9) yields m(t,O) = 1
V, + (VU), = LG
(33)
where v = UL
U = (u - 7dLldr)lL
U = [m,mR,mu,mEIT F = [mu,mRu,muu,mui?]T
(34) (35)
G=
u(t,O) = u. = 1 + u sin (2nSt . t) (25)
E(t,O) = u. tan 0, where a = a*luo*
law
R(t,O) = 1
St = w*R$/2ru* 0
(36)
(26)
(27) a is the nondimensional amplitude of the axial velocity fluctuations at the nozzle exit, and St is the Strouhal number, which is the ratio of the frequency of the axial
the superscript T denotes transpose, and (16) and (19) have been used to obtain the equation for the fourth component of the vector U. Equation (33) is written in weak conservation law form, for the source term LG cannot be written in divergence form, and can be interpreted as the trans-
Appl. Math. Modelling,
1992, Vol. 16, September
467
Adaptive
methods
for annular liquid jets: J. I. Ramos
port of V with a convection velocity I?. Note that iTL represents a relative velocity. Equation (33) can be solved for V, and (34) can then be used to calculate U = V/L. However, (33) represents a system of coupled partial differential equations for V which requires the values of L and dLldT for the evaluation of U and V. Therefore an equation for L is needed. Such an equation can be obtained from (13) applied at the convergence point as follows. At the convergence point (cf. equations (13), (21), and (32)), RJT,.Z = L(t)) = 0
b(t,z = L(f)) = 2R(t,7, = L(t))
(37)
Use of (22) in (37) implies that the following algebraic equation must be satisfied at the convergence point: 2R2(t,z = L(r)) = b;
m(r,z = L(r)) R$
Differentiation of (38) with respect to r and use of (16), (19), and (32) imply that dL
ryl
T=O ,
I
01%
;
i=l Figure 2.
=:c!= I
i-l
2 Schematic
’i I+1
i
_
.
I-1
I
of the grid used in the calculations
dynamics equations, and all the dependent variables were evaluated at the grid points. Note that the nozzle exit and the convergence point correspond to i = 3/2 and i = I - l/2, respectively, where I denotes the number of grid points in the q-direction. We will use equally spaced grids throughout this section, although the finite volume formulation presented in the next paragraphs can easily be extended to nonequally spaced grids. Equation (33) can be integrated from qi_ ,I* to qi+ IIZ(Figure 2) to yield ‘I, +
a a7 ‘I,-
dL
I
112
Vdv
+
(W;+
1/2
-
(vQ)i~
I/Z
I/Z
dt=z
4R(+L)-$&(m~)
at
ZZ
77 _1 _
4R”R_&?!? a~ RX a? which is an ordinary differential
(39) equation for L.
3.2. Governing equations in nonconservation form Equation (33) can also be written as
au
~+Z@=G
law
Assuming that V and G are constant in [vi_ l12,77i+,12]r (43) can be written as av; %Ar)
+
(V4;+,/2
-
for i = 2, 3, . . . , I AT = ~;+I/z -
%A7
and I is the unit matrix and the superscript T denotes transpose. Equations (33) and (40) are subject to the following boundary conditions at the nozzle exit (cf. equations (25) and (26)): U(7,O) = [l ,l ,LI~,u~tan OolT
(42)
methods
$+
1/2
Modelling,
1992,
Vol. 16, September
=
(46)
1, where (cf. equation (34))
Uj + Ir,+ 1 2
=- 1 -huLli
(muLh+l
1
(47) - 277i+ 112 2 2L [ (mL); + hQ+ I and G+v~ has been written as a function of some of the components of the vector V. The value of VIPli2 can be determined by linear extrapolation from the values of VI_, and VI_, and can be written as V,hl/2
4.1. Conservative block-implicit method The fluid dynamics equations written in weak conservation law form (cf. equation (33)) were discretized in the partially staggered grid shown schematically in Figure 2, using finite volumes as described in the next paragraphs. The grid shown in Figure 2 was selected so that the fluxes at the nozzle exit could be implemented exactly in the finite difference form of the fluid
Math.
(45)
71~112
for i = 2, 3, . . . , I -
Appt.
1, where
+ V;iTr+,,?.- ViPlti;_l/2 = LAqGi
(41)
468
(44)
LArlG
(40)
877
Numerical
=
Since u 2 dLldr, ii 2 0, and the convective flux in (44) can be written, using upwind differences, as
where
4.
(W-m
=
(3V,_, - V,&J/2
(48)
Equation (44) involves the values of L and dLldr through (46), and the vector G can be written as
(49)
Adaptive
methods
for annular liquid jets: J. I. Ramos
nonlinear finite difference equation was solved by means of the Newton-Raphson method as follows. At the (k + 1)th iteration within the time step, (55) can be written as
where
Wh” = W” + A7Sk+’
(50) Therefore G involves first- and second-order derivatives of R with respect to 7. These derivatives can be evaluated as follows (cf. Figure 2):
(56)
where AT is the time step and the superscript n denotes the nth time level, i.e., T” = n AT. If S”+’ is quasilinearized with respect to the value of Wk, one obtains Sk” = S’! +
J(W”+’
-
WL)
where
R,+ 1,~- R;- I/Z
i = 2,3,.
. )I - 1
x
An
(58)
(51) =
L(R,+~,2
An2
-
is a Jacobian matrix. Substitution of (57) into (56) yields
2R; + Rj_,,Z) i=2,3,...
Wk+’ = (I - ATJ))‘[ W” + A7(Sk - ATJW”)] ,I - 1
The spatial derivatives of (39) can be discretized the grid schematically shown in Figure 2 as R,m 112
-
RIG
(52)
(59)
in
which was solved iteratively at each time step until the following convergence criterion was satisfied: Ii2 I- 1 x AW:. Aw, iel (60) L i=2
1
I
Note that the values of 14, R, E, and m at i = I 4 are calculated by linear extrapolation as indicated in (48). Note also that (53) depends on m, L, R, u, and 0 at i = I - 4. The values of U at i = I - 4 depend on the values of U from i = $ to i = 1 - I through (46). Therefore (46) and (53) are nonlinearly coupled. Equation (53) can be substituted into (47), and the resulting values of ii can then be substituted into (46) and (53). The resulting system of equations and (53) can then be written as CYW -= S a7
(54)
where W=(V~,V~‘,...,
V,T_,,L)T
(55)
and the vector S can easily be calculated from (36), (46), and (51)-(53). Equation (55) represents a system of nonlinearly coupled, first-order, ordinary differential equations for the vector W, because the spatial derivatives have been discretized while the independent variable r has been kept continuous, and can be solved by means of ordinary differential equations integrators. In this paper, (55) was discretized by using tirstorder backward differences in time, and the resulting
-= dr
4&
where E, is a user-specified error tolerance, AW = WA+’ - Wh, and the left-hand side of (60) represents the LZ-norm of the vector AW. Note that the quasi-linearization of S”+’ (cf. equation (57)) in the Newton-Raphson method yields a sparse, but not banded, matrix, because the kth component of the vector S, i.e., Sk, depends on V,_, , VA, Vk+ , , L, V,_,, and VIPz, as indicated in (46), (47), (53), and (55). Note also that since the column vector W has 4(Z - 2) + 1 components, the dimensions of the Jacobian matrix are (41 - 7) X (41 - 7). 4.2. Nonconservative block-implicit method A method of lines technique similar to the one presented in the previous subsection can be derived by discretizing the governing equations written in nonconservation law form (cf. section 3.2) as follows. The spatial derivatives of (40) can be discretized in an equally = 1,wherel spaced grid such that 7, = Oandv, denotes the number of grid points, by means of backward differences for the advection term, and central differences for G, and the resulting system of nonlinearly coupled, first-order, ordinary differential equations for Ui, can be written as aW -= f3r
S
(61)
where w = (UT, UT,. . . ) Wl,LLT
RX
(62)
(63)
RI - RI-, -- b,* m, - ml_, A77
(57)
Arl
Appl. Math. Modelling,
1992, Vol. 16, September
469
Adaptive
methods
for annular liquid jets: J. I. Ramos kfl
The values of u, R, E, and m at i = I can be calculated by linear extrapolation as ur = 2u,_,
- u,_,
(64)
S can easily be found from the discretization of the spatial derivatives of (40) and (63), and (63) has been included in (61). It must be noted that the derivatives aRl?~r] and a*RlQ*, which appear in the vector G (cf. equation (50)), were discretized as
These values of Uk + ’ , Lk + ’ and (dL/dr)” + ’ can then be employed to evaluate C and G, and this iterative procedure can be repeated until the following convergence criterion is satisfied:
Li=l
1 l/2
Ri+l - Ri_1 + O(AT~)
2Arl
i = 2,3,.
+
. . ,I -
1 (65)
Ri+ 1 - 2Ri + Rip 1 + O(Aq*)
AT2
i = 2,3,.
. . )z - 1
(66)
and RI was calculated from (64). Equation (61) has the same form and the same number of unknowns as (54) and was solved by means of the Newton-Raphson method until (60) was satisfied.
4.3. Nonconservative block-bidiagonal method Equation (40) can also be discretized in the equally spaced grid of the previous section by means of backward differences in time, upwind differences for the advection term, and central differences for G, and the resulting O(Ar,Aq)-accurate finite difference equation can be written as -Cl+‘uy:,l = ArG:+’
+ (I + C;+l)u;+’ + U;
i=2,3,...
,I-
Appl.
Math.
(68)
Modelling,
(70)
In all the calculations presented in this paper the first iteration was always performed with the values of U, L, and dLldr corresponding to the nth time level. It must be noted that the block-implicit methods presented in sections 4.1-4.3 solve the discretized form of the fluid dynamics equations from the nozzle exit up to the grid point closest to the convergence point and that the values of the annular liquid jet’s mass per unit length, radius, and axial and radial velocity components at the convergence point are determined by linear extrapolation (cf. equations (48) and (64)). However, the fluid dynamics equations derived in Ref. 1 are not valid at the convergence point where a meniscus is formed and a pressure gradient may exist, because axial and radial pressure gradients in the annular liquid jet were neglected in the derivation of (l)-(5). It must also be noted that the grid used in section 4.1 is different from the ones used in sections 4.2 and 4.3 and that the block-implicit method of section 4.1 is conservative, whereas the methods of sections 4.2 and 4.3 are not conservative.
5.
and the subscript i denotes the ith grid point, i.e., 77i = i An, where A7 is the spatial step size. Equation (67) represents a (banded) block-bidiagonal matrix of dimensions (41 - 8) x (41 - 8), which can be easily solved by forward substitution. However, since both H and C depend on dLldr and L, and since dLldr can be evaluated as in (63) where VI can be calculated by linear extrapolation as in (64), substitution of (63) into (41), (67), and (68) would yield a sparse system of nonlinear algebraic equations, because U, depends on the value of Ui, i = 2, 3, . . . , I - 1, through (64) and (67). Such a sparse system can be solved by means of the Newton-Raphson method presented in sections 4.1 and 4.2. However, the inversion of the Jacobian matrix of the Newton method whose dimensions are identical to those of sections 4.1 and 4.2 can be avoided if one uses the following iterative technique for both U and L. Assume that Uk, Lk, and dLkldr are known at the kth iteration. These values can be used to evaluate C and Gin (67). The solution of (67), i.e., Uk+‘, can then be used in (63) to obtain (dLldr)k+’ and
470
9 10-4
1 (67)
where C = HAr/An
(Lk+’ - LA)*
1992, Vol. 16, September
Steady state and initialization
The fluid dynamics equations and the numerical methods presented in this paper are applicable to both annular liquid jets (bg # 0) and annular liquid membranes (b$ = 0) and can be used to study steady-state and time-dependent problems. Under steady-state condition, a = 0, i.e., the boundary conditions at the nozzle exit are independent of time, and dLldr = 0. In order to determine the steady-state geometry of annular liquid jets, one can simply set all the partial derivatives with respect to r equal to zero in the equations presented in sections 2-4 and solve the resulting finite difference equations until the convergence criteria deIined in (60) and (70) are satisfied. Alternatively, one can solve the time-dependent equations subject to time-independent boundary conditions at the nozzle exit until an asymptotic steady state is reached. This time-dependent approach requires that (60) and (70) be satistied at each time step, and the asymptotic steady solution is reached whenever the following steady-state convergence criteria are satisfied: I-‘AWiT. AW, ‘I*
Iz
i=2
Ar2
5
1
E*
(71)
Adaptive methods for annular liquid jets: J. I. Ramos rA (u?+l - u;)r.
t 1;=
(u:+l
I
where AW=
_ uy) + (L”+I
A72
W”+’ - W”
_ ~721i’~ zs
exit angle, and amplitude and frequency of the fluctuations. In order to assess the influence of these parameters on the dynamic response of liquid jets, parametric studies were performed to determine the convergence length and the pressure coefficient as functions of time. The results of some of these parametric studies are illustrated in Figures 3-13, which show the (nondimensional) axial velocity component at the nozzle exit (u,, = u(T,O)), the pressure coefficient, and L’ = L/L(O) for the values of the parameters shown in Table I and b$ = 0, i.e., annular liquid membranes. The values of L(0) are equal to 12.5590 for Figures 3-5, II, and 12; 9.8652 for Figure 6; 0.3262
(73)
lZ
and c3 are user-specified error tolerances. In the calculations presented in this paper, cl = l2 = l3 = 10-4. Equation (71) was used with the methods of lines of sections 4.1 and 4.2, while (72) was employed with the block-bidiagonal technique of section 4.3. The time-dependent approach, which may be used to determine the asymptotic steady state, may require shorter computational times than the solution of the steady-state governing equations, for the discretization of the time-derivative terms increases the diagonal dominance of the matrices to be inverted when solving the fluid dynamics equations (cf. equation (67)). In this paper the steady-state fluid dynamics field was determined asymptotically from the solution of the time-dependent equations. Once steady state was achieved, 7 was reset to 0, and m,(O) = (rn;).$.,and p;(O) = (PJSS, where ss denotes steady state. These initial conditions were then used to analyze the response of annular liquid jets to fluctuations in the axial velocity at the nozzle exit in the absence of mass injection into the volume enclosed by the annular liquid jet. 6.
Presentation
(72)
c3
1
Ucl
0.6
0.2
of results
Numerical experiments were performed to determine the efficiency of the methods presented in sections 4.1-4.3 and the effects of both the time step and number of grid points on the response of annular liquid jets subject to mass flow rate fluctuations at the nozzle exit. The results of these experiments indicate that the techniques presented in sections 4.1-4.3 essentially yield the same results and that the iterative, blockbidiagonal method of section 4.3 is more efficient than the block-implicit techniques of sections 4.1 and 4.2 because these techniques require the inversion of the Jacobian matrix, whereas the method of section 4.1 involves a block-bidiagonal matrix and the vector U can easily be found by forward substitution (cf. equation (67)). The number of grid points, I, and the time step, AT, used in the calculations were varied in order to ensure grid independence. For all the calculations presented in this paper, 400 5 I 5 1600 and 0.04 5 AT 5 0.001, and the largest value of I corresponded to annular liquid jets whose steady-state convergence length was longest. The smallest value of Ar corresponds to the highest value of the Strouhal number employed in the calculations. As shown in section 2, the dynamic response of inviscid, isothermal, annular liquid jets to velocity fluctuations at the nozzle exit depends on the Froude and Weber numbers, C,,,, , initial pressure ratio, nozzle
L’
1.0
-0.2
1 0
IO
20
30
40
Figure 3. Convergence length, pressure coefficient, velocity at the nozzle exit as functions of time
50
’
and axial
L’
1.0
Uo 0.6 -
-0.2 1 0
10
20
30
40
Figure 4. Convergence length, pressure coefficient, velocity at the nozzle exit as functions of time
Appl.
Math. Modelling,
50
’
and axial
1992, Vol. 16, September
471
Adaptive
methods
for annular liquid jets: J. I. Ramos
20
10
30
40
Figure 5. Convergence length, pressure coefficient, velocity at the nozzle exit as functions of time
50
1
-0
and axial
10
20
30
50
40
Figure 8. Convergence length, pressure coefficient, velocity at the nozzle exit as functions of time
’
and axial
r 0.6
0.6
0.2
0.2
-0.2
/‘\%
1
0
10
I
20
30
40
Figure 6. Convergence length, pressure coefficient, velocity at the nozzle exit as functions of time
50
’
-I
-0.2 0
CPll
10
20
30
40
Figure 9. Convergence length, pressure coefficient, velocity at the nozzle exit as functions of time
and axial
1.0
1.0
50
’
and axia
I!
"0
“0
0.6 .
0.6
0.2 /_
-0.2
’
10
0
20
I 30
40
Figure 7. Convergence length, pressure coefficient, velocity at the nozzle exit as functions of time
472
Appl.
Math. Modelling,
50
1
and axial
1992, Vol. 16, September
Figure 10. Convergence length, pressure coefficient, velocity at the nozzle exit as functions of time
and axial
Adaptive
methods
for annular liquid jets: J. I. Ramos
for Figure
7; 19.0452 for Figure 8; 3.8236 for Figure 9; 37.3466 for Figure 10; and 10.2816 for Figure 13. Figure 3 indicates that the convergence length and
1.0
the pressure coefficient do not respond instantaneously to the velocity fluctuations at the nozzle exit. The delay times in both the convergence length and the pressure coefficient are functions of the parameters that govern the dynamics of liquid jets, as discussed in the next paragraphs. Figure 3 also indicates that the convergence length and the pressure coefficient are periodic functions that have the same frequency as the axial velocity fluctuations at the nozzle exit, that the pressure coefficient has a smaller delay time than the convergence length, and that there is a lag between the peaks of the pressure coefficient and those of the convergence length. Figures 3 and 4 show that the delay and lag times in the response of both the convergence length and the pressure coefficient are nearly independent of the amplitude of the axial velocity fluctuations at the nozzle exit, whereas the amplitude of the oscillations of both the convergence length and the pressure coefficient increases as the amplitude of the axial velocity fluctuations at the nozzle exit is increased. Figures 3 and 4 also show that the locations of the maxima and minima of the pressure coefficient are nearly independent of the amplitude of the axial velocity fluctuations at the nozzle exit. However, striking differences are observed in the convergence lengths and pressure coefficients shown in Figures 3 and 4. In particular, Figure 3 indicates that both L and C,,, are nearly sinusoidal, whereas the values of L and C,, shown in Figure 4 deviate considerably from the sinusoidal behavior of the axial velocity fluctuations at the nozzle exit. This result is a consequence of the fact that for small amplitudes of the axial velocity component at the nozzle exit the response of liquid jets is nearly linear, as illustrated in Figure 3; for large amplitudes, however, nonlinear effects are important, and the response of the liquid jet is nonlinear. Figures 3 and 5 illustrate the effects of the Strouhal number on both the convergence length and the pressure coefficient and indicate that the amplitude of the oscillations of both L and C,, decreases, whereas their frequency increases, as the Strouhal number is increased. These figures also show that the delay time in the convergence length decreases as the Strouhal number is increased. Figure 5 indicates that for St = 0.5 the pressure coefficient is nearly zero. The effect of the Froude number on the dynamic response of liquid jets is illustrated in Figures 3 and 6, which indicate that the delay time in the convergence length increases as the Froude number is increased. Figures 3 and 6 also indicate that both the convergence length and the pressure coefficient exhibit a behavior that increasingly deviates from the sinusoidal one exhibited by the axial velocity fluctuations at the nozzle exit as the Froude number is increased. This is due to the fact that the values of the convergence length at steady state decrease as the Froude number is increased. Note that since the steady-state convergence
0.6
I:lIYWL 0.2
0
10
20
30
40
50
t
Figure 11. Convergence length, pressure coefficient, and axial velocity at the nozzle exit as functions of time
II I.0
, "0
0.6
0.2
-0.;
-0.t -0
10
20
30
40
50
1
Figure 12. Convergence length, pressure coefficient, and axia velocity at the nozzle exit as functions of time
L 1.0 “0
0.6
0.2
-0.2
-0.6
% 0
10
20
30
40
Figure 13. Convergence length, pressure coefficient, velocity at the nozzle exit as functions of time
50
I
and axial
Appl.
Math.
Modelling,
1992,
Vol. 16, September
473
Adaptive
methods for annular liquid jets: J. I. Ramos Table 1.
Values
of the parameters used in the calculations P?(O)
Figure
Ff
We
00
Cpmax
P,”
a
St
3 4 5 6 7 8 9 10 11 12 13
10 10 10 10,000 10 10 10 10 10 10 10
50 50 50 50 1 100 50 50 50 50 50
0 0 0 0 0 0 -15 15 0 0 0
1 1 1 1 1 1 1 1 2 10 1
1 1 1 1 1 1 1 1 1 1 0.5
0.1 0.25 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1
0.1 0.1 0.5 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1
length is a measure of the volume enclosed by the liquid jet, the deviation from a sinusoidal behavior increases as the volume enclosed by the jet at steady state decreases. Figures 3, 7, and 8 show the effects of the Weber number on both the convergence length and the pressure coefficient, and they indicate that the delay times in both L and C,,, decrease, whereas the amplitudes of both the convergence length and the pressure coeflicient increase, as the Weber number is decreased. Figure 8 also shows that the steepness of the convergence length from its maximum to its minimum values increases as the Weber number is increased, i.e., as the convergence length at steady state and the inertia of the liquid jet are increased. Figures 3, 9, and 10 illustrate the effects of the nozzle exit angle on the dynamic response of liquid jets, and they indicate that the delay times and the amplitude of both the convergence length and the pressure coefficient increase as the nozzle exit angle is increased. This result is consistent with the steady-state convergence lengths, which increase at the nozzle exit angle is increased,‘,3.4 and indicates that the delay times decrease as the volume enclosed by the jet as steady state is decreased. Figure 10 also shows that the steepness of the convergence length from its maximum to its minimum values increases as the nozzle exit angle is increased. This result is also consistent with the results shown in Figure 8 and indicates that the steepness of the convergence length increases as the convergence length at steady state is increased. Figure 10 also shows the dynamic response of the liquid jet has not reached a periodic behavior even for times greater than 40. The effect of C,,,,, which is a measure of the importance of surface tension and pressure of the gases that surround the liquid (cf. equation (31)), on the dynamic response of liquid jets is illustrated in Figures 3, 11, and 12, which indicate that the delay time in both the convergence length and the pressure coefficient decreases as C,,,, is increased. The amplitudes of both the convergence length and the pressure coefficient decrease and increase, respectively, as C,,,, is increased. It must be noted that the value of C,,,, does not affect the convergence length at steady state for the
474
Appl.
Math.
Modelling,
1992,
Vol. 16, September
values of the parameters presented in Table 1. However, C,,,, P lays the role of an amplification factor, for it increases the pressure coefficient if there is a pressure difference between the gases enclosed by and surrounding the liquid jet (cf. equation (30)). Figures 3 and 13 illustrate the effects of the initial pressure ratio on the dynamic response of liquid jets, and they indicate that the delay times of both the convergence length and the pressure coefficient decrease as the underpressurization of the jet is increased, i.e., as the convergence length at steady state is decreased. 7.
Conclusions
A domain-adaptive technique and three block-implicit methods have been used to analyze the dynamic response of inviscid, isothermal, annular liquid jets subject to velocity fluctuations at the nozzle exit as a function of the Froude and Weber numbers, initial pressure ratio, nozzle exit angle, surface tension, and amplitude and frequency of the imposed velocity fluctuations. The domain-adaptive technique maps the timedependent, curvilinear geometry of the jet into a unit interval and yields a system of integrodifferential equations for the annular liquid jet’s mass per unit length, mean radius, and axial and radial velocity components. The integrodifferential character of the equations is due to both the equation that governs the location of the time-dependent convergence point and the condition that the mass of the gases enclosed by the annular liquid jet is constant. Two of the block-implicit methods presented in this paper use a method of lines technique, which discretizes the spatial derivatives and keeps continuous the time variable, and yield a system of nonlinearly coupled, first-order ordinary differential equations for the convergence length and the mean radius, mass per unit length, and axial and radial velocity components of the annular liquid jet. This system of ordinary differential equations has been solved by means of the Newton-Raphson method, whose Jacobian matrix is sparse but not banded. The third block-implicit method is an iterative technique, which discretizes both the spatial and the temporal derivatives and yields a block-bidiagonal matrix for the mean radius, mass per unit length, and axial
Adaptive
and radial velocity components of the annular liquid jet. This block-bidiagonal matrix is solved iteratively by forward substitution together with the ordinary differential equation for the convergence length. It has been shown that the block-bidiagonal method is more efficient than and yields the same results as the blockimplicit techniques presented in this paper. It has also been shown that both the convergence length and the pressure coefftcient are periodic functions of time when the liquid jet is subject to sinusoidal velocity fluctuations at the nozzle exit and that there is a delay time in the dynamic response of the jet. The delay time is nearly independent of the amplitude of the velocity fluctuations at the nozzle exit, and it decreases as the frequency of these fluctuations is increased. The delay time in the convergence length increases as the Froude and Weber numbers and nozzle exit angle are increased and decreases as the underpressurization of the liquid jet is increased. It is also shown that there is a lag between the response of the convergence length and that of the pressure coefficient and that both the delay time and the lag time depend on the inertia of the annular liquid jet. The convergence length and the pressure coefficient are periodic functions of time and have the same frequency as that of the axial velocity fluctuations at the nozzle exit, because the liquid jets analyzed in this paper are inviscid and isothermal. When the amplitude of the velocity fluctuations at the nozzle exit is small, the convergence length and the pressure coefficient are nearly sinusoidal, and the liquid jet behaves almost linearly. For large amplitudes, however, nonlinear effects are important, and the response of the jet is periodic but not sinusoidal. It has also been found that the steepness of the convergence length from its maximum to its minimum values increases as the Froude and Weber numbers and nozzle exit angle are increased. An increase in the Froude number corresponds to a decrease in both the force of gravity and the steady-state convergence length, whereas the convergence length at steady state increases as the Weber number and the nozzle exit angle are increased. The results shown in this paper indicate that the dynamic response of liquid jets to axial velocity fluctuations at the nozzle exit is very much dependent on the magnitude of the gravitational force and on the volume enclosed by the annular liquid jet at steady state.
methods for annular liquid jets: J. I. Ramos
Acknowledgments The calculations presented in this paper were performed at the Facultad de InformaticaiEscuela Tecnica Superior de Ingenieros de Telecommunication of the Universidad de Malaga. The author is grateful to Mr. Juan Falgueras for his help with both the graphs and the typing of this paper. References
5
12 13 14
15 16
17
18
Appl.
Ramos, J. I. Liquid curtains, I, Fluid mechanics. Chem. Eng. Sci. 1988.43, 3171-3184 Ramos, J. I. and Pitchumani, R. Liquid curtains, II, Gas absorption. Chem. Eng. Sci. 1990, 45, 1595-1604 Ramos, J. I. Annular liquid jets in zero gravity. Appl. Math. Modelling 1990, 14, 630-640 Ramos, J. I. Analytic, asymptotic and numerical studies of liquid curtains and comparisons with experimental data. Appl. Math. Modelling 1990, 14, 170-183 Ramos, J. I. Dynamic response of liquid curtains to time-dependent pressure fluctuations. App[. Muth. ModeDing 1991, 15, 126-135 Binnie, A. M. and Squire, H. B. Liquid-jets of annular cross section. Engineer (London) 1941, 171, 236-238 Baird. M. H. I. and Davidson, J. F. Annular iets. I. Fluid dvnamics. Chem. Eng. Sci. 1962, 17, 467-472 ” Baird, M. H. I. and Davidson. J. F. Annular jets. II, Gas absorption. Chem. Eng. Sci. 1962, 17, 473-480 Roidt, R. M. and Shapiro, Z. M. Liquid curtain reactor. Report No. 85M981, Westinghouse R&D Center, Pittsburgh, PA, 1985 Lee, C. P. and Wang, T. G. A theoretical model for the annular liquid jet instability. Phys. Fkds 1986. 29, 2076-2085 Hoffman, M. A., Tabahashi, R. K., and Monson, R. D. Annular liquid jet experiments. ASME J. Fluids Eng. 1980, 102, 344-349 Lee, C. P. and Wang. T. G. A theoretical model for the annular liquid jet instability-Revisited. Phys. Fhrids A 1989. 1, 967-974 Meseguer. J. The dynamics of small annular jets. ASME J. FINid. Eng. 1988, 110, 123-126 Ramos, J. I. Liquid membranes: Formulation and steady state analysis. Report No. CO/89/4, Dept. of Mech. Eng., Carnegie Mellon Univ., Pittsburgh, PA. 1989 Ramos, J. I. Annular liquid jets: Formulation and steady state analysis. Z. Anger. Moth. Mech. 1991, in press Chigier, N., Ramos. J. I., and Kihm, K. Experimental and theoretical studies of vertical annular liquid jets. Proceedings c?fthe Sixth Symposium on Energy Engineering Sciences: Flow und Trunsport in Continua, U.S. Department of Energy, Argonne National Laboratory, Argonne, IL, 4-6 May 1988, pp. 18-25 Ramos, J. I. Adaptive finite difference methods for liquid membranes, I. Mass flow rate fluctuations at the nozzle exit. Report NO. CO/91/1, Dept. of Mech. Eng., Carnegie Mellon Univ., Pittsburgh, PA, January 1991 Ramos, J. I. Adaptive finite difference methods for liquid membranes. II, Mass injection. Report No. CO/91/2, Dept. of Mech. Eng., Carnegie Mellon Univ., Pittsburgh, PA, January 1991
Math. Modelling,
1992, Vol. 16, September
475
of Mechanical
USA, and Facultad Spain
Engineering,
Carnegie Mellon
de Informcitica1E.T.S.I.
University, Pittsburgh, PA, Telecomunicacibn, Universidad de Mdlaga,
A domain-adaptive technique is used to analyze the dynamic response of inviscid, isothermal, annular liquidjets to velocityfluctuutions at the nozzle exit. The adaptive technique maps the unknown, timedependent, curvilinear geometry of the annular jet into a unit interval and yields a system of integrodifferential equations for the muss per unit length, radius, and axial and radiul velocity components of the liquidjet. The convergence length is governed by an ordinary differential equation that depends on the values of the dependent variables at the convergence point. The governing equations are solved by means of two methods of lines, which discretize the sputial coordinate but keep continuous the time, and a Newton method. A block-bidiagonal technique is also used to solve the fluid dynamics equations. It is shown that the dynamic response of liquid jets to sinusoidal velocity jluctuutions at the nozzle exit is periodic and nearly sinusoidal for small amplitudes of the velocity oscillations at the nozzle exit. For large amplitudes the fluctuations of both the convergence length and the pressure of the gases enclosed by the annulur jet are periodic but not sinusoidal. It is ulso shown thut there is a delay time between the velocity fluctuations ut the nozzle exit und those of the convergence length. This delay time is nearly independent of the amplitude of the velocity fluctuutions ut the nozzle exit and decreases as the Weber and Froude numbers and the nozzle exit angle ure decreased. Keywords:
1.
block-implicit techniques, methods of lines, annular liquid jets, block-bidiagonal
Introduction
Annular liquid jets have been the subject of intensive research in recent years, ‘-I3 because under certain conditions they can form enclosed volumes, which can be used to determine the dynamic surface tension of liquids and to burn toxic wastes. The equations that govern the fluid dynamics of annular liquid jets have been derived by the author for Newtonian and non-Newtonian fluids’4*‘5 by using Cauchy’s equations and applying interface boundary conditions. The interfaces between the liquid and the fluids surrounding and enclosed by the annular liquid jet are material surfaces where the shear stress is continuous, while the jump in normal stresses across the interfaces is balanced out by surface tension. The author derived the asymptotic equations that govern the fluid dynamics of inviscid, isothermal, annular liquid jets (Figure 1) by assuming that the jet thickness is smaller than the jet radius at the nozzle exit, averaging the governing equations across the anAddress reprint requests to Dr. Ramos at the Department of Mechanical Engineering, Carnegie Mellon University, Pittsburgh, PA 15213-3890, USA. Received 14 March ruary 1992
464
1991; revised
6 February
Appl. Math. Modelling,
1992; accepted
14 Feb-
1992, Vol. 16, September
method
nular jet, and applying the interface boundary conditions at the jet’s mean radius.‘4.‘5 The resulting set of asymptotic equations is also valid for thin, viscous annular jets except near the nozzle exit, where the velocity profile relaxes from stick boundary conditions at the nozzle walls to slip/free boundary conditions at the annular jet’s interfaces, and near the convergence point, where the radius of the annular liquid jet’s inner interface is zero, a meniscus is formed, and an axial pressure gradient may exist. The equations derived in Refs. 14 and 15 are exact for inviscid, isothermal, annular liquid membranes, i.e., annular liquid jets of zero thickness, and have analytical solutions for steady-state conditions if both the slope and the curvature of the jets are small.‘-4 For time-dependent problems, however, analytical solutions are difficult, if not impossible, to obtain, and numerical techniques must be used. In this paper the dynamic response of inviscid, isothermal, annular liquid jets to mass flow rate fluctuations at the nozzle exit is analyzed by means of three adaptive, block-implicit techniques. The domain-adaptive techniques and the numerical studies presented in this paper were motivated by some experimental results,16 which indicate that under certain conditions the convergence length of annular liquid jets may exhibit a periodic behavior. This periodic behavior may be due to periodic fluctuations in the pressure of the 0 1992 Butterworth-Heinemann
Adaptive
Figure 1.
Schematic
methods
for annular liquid jets: J. I. Ramos
terns, i.e., nonlinear dynamical systems governed by ordinary differential equations. It must be indicated that annular liquid jets absorb the gases that they enclose at a rate that depends on the pressure difference across the jet, gravitational acceleration, surface tension, annular jet’s thickness, nozzle exit angle, and diffusivity of the gases in the liquid.2.3 Therefore in order to maintain steady state, mass must be injected into the volume enclosed by the annular jet at a rate equal to the mass absorption rate by the liquid. If the mass injection rate exceeds the mass absorption rate, mass will accumulate in the volume enclosed by the liquid jet, and the pressure of the enclosed gases will increase until the annular jet becomes unstable, breaks up, and does not form an enclosed volume. On the other hand, if the mass absorption rate by the liquid is greater than the mass injection rate, the annular jet will collapse. Both the growth and the collapse of annular liquid jets are timedependent phenomena which require that the gas concentration in the liquid be determined together with the solution of the equations that govern the unsteady fluid dynamics of annular liquid jets. In this paper it is assumed that the liquid does not absorb the gases that it encloses and that mass is not injected into the volume enclosed by the annular jet. Analytical studies of annular liquid jets indicate that there is a critical value of the pressure of the gases enclosed by annular liquid jets beyond which the jet cannot form a closed volume.‘.‘.4 These analytical results were obtained by assuming steady-state conditions and annular jets of small slope and curvature. Under transient conditions it may be possible to exceed this critical pressure while still obtaining a stable, yet fluctuating, annular liquid jet.5.‘7.‘x Numerical studies of annular liquid jets subject to fluctuations in the pressure of the gases that surround the jet’ conf-irm that this is indeed possible, provided that the annular liquid jet is subjected to high-frequency unsteady pressure fluctuations. The domain-adaptive technique presented in this paper maps the unknown, time-dependent geometry of the liquid jet into a unit interval and yields a system of integrodifferential equations for the convergence length, i.e., the axial distance at which the annular liquid jet merges on the symmetry axis to become a solid jet, and the annular liquid jet’s mean radius, mass per unit length, and axial and radial velocity components. Two of the block-implicit methods presented in the paper use a method of lines technique, which discretizes the spatial derivatives, keeps continuous the time, and results in a system of nonlinear, first-order, ordinary differential equations in time that was solved by means of the Newton method. The third blockimplicit method discretizes both the spatial and the temporal derivatives and yields a block-bidiagonal matrix, which can easily be solved by forward substitution. Furthermore, one of the block-implicit methods uses a conservative, finite volume discretization of the nondimensional equations written in weak conserva-
of an annular liquid jet
gases that surround the jet, unsteady mass flow rate fluctuations at the nozzle exit, and/or fluctuations in the mass injected into the volume enclosed by the jet. The effects of fluctuations in the pressure of the gases that surround the jet have been investigated by Ramos,j who used an adaptive, block-bidiagonal technique and showed that the convergence length and the pressure coefficient are periodic functions with the same period as that of the pressure fluctuations. However, the numerical calculations presented in Ref. 4 are academic in nature, for they were based on the assumption that the pressure difference between the gases enclosed by and surrounding the liquid jet was a known function of time, and they did not account for the variations of the volume enclosed by the liquid. Furthermore, the nondimensionalization of the governing equations presented in Ref. 4 can only be used to study the dynamics of annular liquid jets in gravitational environments, because the acceleration of gravity was used to nondimensionalize the time and the radial and axial velocity components of the liquid jet. In this paper a nondimensionalization of the governing equations that is valid for annular liquid jets in gravitational and nongravitational environments is introduced, and the pressure across the annular jet is calculated by accounting for the variations in the volume enclosed by the jet. The effects of fluctuations in the mass flow rate at the nozzle exit and in the mass injected per unit time into the volume enclosed by annular liquid jets have also been studied by the author,17.1X who used a blockbidiagonal method and a method of lines technique. The numerical studies presented in Refs. 5, 17, and 18 indicate that the unsteady dynamics of inviscid, isothermal, annular liquid jets subject to fluctuations in the pressure of the gases enclosed by the annular jet, fluctuations in the mass flow rate at the nozzle exit, and fluctuations in the mass injected per unit time into the volume enclosed by the annular jet has some similarities with the nonlinear dynamics of lumped sys-
Appl.
Math.
Modelling,
1992, Vol. 16, September
465
Adaptive methods for annular liquid jets: J. I. Ramos tion law form, whereas the other two methods are based on the discretization of the nondimensional equations written in nonconservation law form. The paper has been arranged as follows. Section 2 summarizes the equations governing the time-dependent fluid dynamics of inviscid, isothermal, annular liquid jets and the nondimensionalization used throughout the paper. Section 2 also identifies the nondimensional parameters that govern the dynamics of liquid jets. In section 3 a domain-adaptive method for the numerical solution of the nondimensional equations is presented. The three block-implicit methods used to discretize the governing equations are presented in section 4, whereas the initialization of the time-dependent calculations is described in section 5. The two sections, Presentation of Results and Conclusions, put an end to the paper. 2.
Fluid dynamics equations
The dimensional equations governing the fluid dynamics of inviscid, isothermal, incompressible, annular liquid jets can be written as’4*‘5
am* i3 t*+az*(m*u*)=O
(1)
-$ (m*u*) + -$ (m*u*‘) = m*g + 2flg
$(m*P)
+
+ R* $
(P,* - Pi*) (2)
-$ (m*u*P)
(4) where the stars denote dimensional quantities: t* is time; z* is the axial coordinate measured from the nozzle exit; m* is the annular jet mass per unit length and per radian; u* is the liquid axial velocity; 3 is the cross-sectional averaged radial velocity; g is the gravitational acceleration; (T is the surface tension; R* is the annular jet’s mean radius; p: and py* denote the pressure of the gases enclosed by and surrounding the annular jet, respectively; and .l*
=Xx/[1
+
(s)‘]“’
(5)
Equations (l)-(5) were obtained by averaging the Navier-Stokes equations across the annular jet thickness and by including the kinematic and stress conditions at the jet’s inner and outer interfaces, which are material surfaces. At each gas-liquid interface the shear stress is continuous, whereas the jump in normal stresses is balanced by surface tension.‘4v’5 Equations (l)-(3) represent the conservation of mass and the conservation of axial and radial momentum,
466
Appl. Math. Modelling,
1992, Vol. 16, September
respectively, whereas (4) corresponds to the kinematic condition applied at the annular jet’s mean radius. For annular liquid jets (Figure I),
R,* = R* + b*/2
RF = R* - b*/2
(6)
where R,* and RT denote the radii of the annular liquid jet’s inner and outer interfaces, respectively, and the annular jet’s mass per unit length and per radian is governed by the following expression: m* = p*R*b*
(7)
Equations (l)-(5) are asymptotic for inviscid, isothermal, annular liquid jets of thickness b* # 0 and are exact for inviscid annular liquid membranes, i.e., for b* = 0, and their numerical solution requires the specification of initial conditions at t* = 0 and boundary conditions at the nozzle exit Z* = 0. The initial and boundary conditions may be functions of z* and t*, respectively. In this paper the following boundary conditions were employed at the nozzle exit:
m*(t*,O) = rng
u*(t*,O)= uo*+ a* sin (w*t*) (8)
Z(t*,O) = u*(t*,O)tan $
R*(t*,O) = R;
(9)
where the subscript 0 denotes conditions at the nozzle exit; m$, R,*, ut, a*, w*, and 0, are time independent; 13,denotes the angle between the symmetry axis and the velocity vector at the nozzle exit; and a* and w* are the amplitude and angular frequency, respectively, of the axial velocity fluctuations at the nozzle exit. Note that since m*(t*,O)is constant, whereas u*(t*,O) is a periodic function of time, the mass flux m*(t*,O)u*(t*,O)at the nozzle exit is a periodic function of time. Note also that for annular liquid membranes. i.e., annular jets of zero thickness, RT = R,* = R*, b* = b,* = 0, and m* is the membrane’s mass per unit length and per radian. In order to maintain steady state, mass must be injected into the volume enclosed by the jet at a rate equal to the mass absorption rate by the liquid. If mass is not injected into the volume enclosed by the annular jet, this will collapse due to mass absorption. In this paper it is assumed that the Peclet number is infinity, so that the liquid does not absorb the gases enclosed by and surrounding the annular jet,2,3 and that mass is not injected into the volume enclosed by the jet. It is further assumed that the gases surrounding and enclosed by the liquid jet are ideal and isothermal and that the gases surrounding the liquid are infinite in extent, so that p: can be assumed to be constant. For ideal gases, pf = p7Rl-T
(10)
where l? is the specific gas constant and pi* and Ti” denote the density and temperature, respectively, of the gases enclosed by the annular jet. Furthermore, for low Mach numbers, ~7 is nearly uniform, and (10) can be integrated to yield
p”V* = rnTl?p
(11)
Adaptive
where L'
L' V*
m* =
mR,@&*
= 0
p,%-R,+‘dz*
(12)
for annular liquid jets: J. I. Ramos
velocity fluctuations at the nozzle exit to a characteristic residence time. Substitution of (14) and (15) into (11) yields
0
V* and m,? denote the volume and the mass of the gases enclosed by the annular jet, respectively, and L* is the convergence length, i.e., the axial distance at which the annular jet’s inner radius vanishes, R,*(t*,L*(t*)) = 0
P~=~(m;/~R~dz.)$=&
u*
M, = rR;‘p$lRT,+
uo*
* m=-
z=L
Ro*
uo*
L = L*IR;
mX
(14)
(13
RZ of (14) and (15) into (l)-(7)
Substitution
m*
t=t*!$
R=g
b=&
mi = rn~IA4~
(13)
v*
jj-
14 = -
(28)
where
Equations (l)-( 13) can be written in nondimensional form by introducing the following parameters:
Substitution
methods
of (14) and (15) into (24) yields
C,...,,[m./~Rfdz - 1) (30)
where Cpmax= p:R;l2a
yields
am a t+$mu)=O
C,,, =
(16)
(29)
(31)
It is clear from (16)-(19) and (30) that the fluid dynamics of annular liquid jets depends on m;, which in turn depends on the mass injected into the volume enclosed by the jet. In the absence of both mass injection into the volume enclosed by the jet and mass absorption by the liquid, mi = m,(O).
(17) (18) aR :=_+uE at
az
,=,/[I
+ ($)2]“2
(19)
(20)
3.
Domain-adaptive
The annular dependent downstream formed into
numerical
technique
liquid jet geometry is curvilinear and time and has an unknown, time-dependent, boundary. This geometry can be transa unit interval by means of the mapping,5
(t,z) + (7371)
r=t
r) = ZIL
OS~~l (32)
R, = R + hi2
Ri = R - b/2
(21)
where L(t) = L*(t*)IR,* is the Jacobian of the mapping.
(22)
3.1. Governing equations in weak conservation form Substitution of (32) into (16)-(19) yields
where Fr = u~~/~R~
We = m$&‘/2aR* 0
C,,, = C, We
C, = (p” - p~)R~Zlm~u,$2
(23) (24)
Fr and We are the Froude and Weber numbers, respectively, and C,, denotes the pressure coefficient. Note that Fr and We are constant and that C,, and C,,,, may be functions of t because p” may be a function of 1. Substitution of (14) and (15) into (8) and (9) yields m(t,O) = 1
V, + (VU), = LG
(33)
where v = UL
U = (u - 7dLldr)lL
U = [m,mR,mu,mEIT F = [mu,mRu,muu,mui?]T
(34) (35)
G=
u(t,O) = u. = 1 + u sin (2nSt . t) (25)
E(t,O) = u. tan 0, where a = a*luo*
law
R(t,O) = 1
St = w*R$/2ru* 0
(36)
(26)
(27) a is the nondimensional amplitude of the axial velocity fluctuations at the nozzle exit, and St is the Strouhal number, which is the ratio of the frequency of the axial
the superscript T denotes transpose, and (16) and (19) have been used to obtain the equation for the fourth component of the vector U. Equation (33) is written in weak conservation law form, for the source term LG cannot be written in divergence form, and can be interpreted as the trans-
Appl. Math. Modelling,
1992, Vol. 16, September
467
Adaptive
methods
for annular liquid jets: J. I. Ramos
port of V with a convection velocity I?. Note that iTL represents a relative velocity. Equation (33) can be solved for V, and (34) can then be used to calculate U = V/L. However, (33) represents a system of coupled partial differential equations for V which requires the values of L and dLldT for the evaluation of U and V. Therefore an equation for L is needed. Such an equation can be obtained from (13) applied at the convergence point as follows. At the convergence point (cf. equations (13), (21), and (32)), RJT,.Z = L(t)) = 0
b(t,z = L(f)) = 2R(t,7, = L(t))
(37)
Use of (22) in (37) implies that the following algebraic equation must be satisfied at the convergence point: 2R2(t,z = L(r)) = b;
m(r,z = L(r)) R$
Differentiation of (38) with respect to r and use of (16), (19), and (32) imply that dL
ryl
T=O ,
I
01%
;
i=l Figure 2.
=:c!= I
i-l
2 Schematic
’i I+1
i
_
.
I-1
I
of the grid used in the calculations
dynamics equations, and all the dependent variables were evaluated at the grid points. Note that the nozzle exit and the convergence point correspond to i = 3/2 and i = I - l/2, respectively, where I denotes the number of grid points in the q-direction. We will use equally spaced grids throughout this section, although the finite volume formulation presented in the next paragraphs can easily be extended to nonequally spaced grids. Equation (33) can be integrated from qi_ ,I* to qi+ IIZ(Figure 2) to yield ‘I, +
a a7 ‘I,-
dL
I
112
Vdv
+
(W;+
1/2
-
(vQ)i~
I/Z
I/Z
dt=z
4R(+L)-$&(m~)
at
ZZ
77 _1 _
4R”R_&?!? a~ RX a? which is an ordinary differential
(39) equation for L.
3.2. Governing equations in nonconservation form Equation (33) can also be written as
au
~+Z@=G
law
Assuming that V and G are constant in [vi_ l12,77i+,12]r (43) can be written as av; %Ar)
+
(V4;+,/2
-
for i = 2, 3, . . . , I AT = ~;+I/z -
%A7
and I is the unit matrix and the superscript T denotes transpose. Equations (33) and (40) are subject to the following boundary conditions at the nozzle exit (cf. equations (25) and (26)): U(7,O) = [l ,l ,LI~,u~tan OolT
(42)
methods
$+
1/2
Modelling,
1992,
Vol. 16, September
=
(46)
1, where (cf. equation (34))
Uj + Ir,+ 1 2
=- 1 -huLli
(muLh+l
1
(47) - 277i+ 112 2 2L [ (mL); + hQ+ I and G+v~ has been written as a function of some of the components of the vector V. The value of VIPli2 can be determined by linear extrapolation from the values of VI_, and VI_, and can be written as V,hl/2
4.1. Conservative block-implicit method The fluid dynamics equations written in weak conservation law form (cf. equation (33)) were discretized in the partially staggered grid shown schematically in Figure 2, using finite volumes as described in the next paragraphs. The grid shown in Figure 2 was selected so that the fluxes at the nozzle exit could be implemented exactly in the finite difference form of the fluid
Math.
(45)
71~112
for i = 2, 3, . . . , I -
Appt.
1, where
+ V;iTr+,,?.- ViPlti;_l/2 = LAqGi
(41)
468
(44)
LArlG
(40)
877
Numerical
=
Since u 2 dLldr, ii 2 0, and the convective flux in (44) can be written, using upwind differences, as
where
4.
(W-m
=
(3V,_, - V,&J/2
(48)
Equation (44) involves the values of L and dLldr through (46), and the vector G can be written as
(49)
Adaptive
methods
for annular liquid jets: J. I. Ramos
nonlinear finite difference equation was solved by means of the Newton-Raphson method as follows. At the (k + 1)th iteration within the time step, (55) can be written as
where
Wh” = W” + A7Sk+’
(50) Therefore G involves first- and second-order derivatives of R with respect to 7. These derivatives can be evaluated as follows (cf. Figure 2):
(56)
where AT is the time step and the superscript n denotes the nth time level, i.e., T” = n AT. If S”+’ is quasilinearized with respect to the value of Wk, one obtains Sk” = S’! +
J(W”+’
-
WL)
where
R,+ 1,~- R;- I/Z
i = 2,3,.
. )I - 1
x
An
(58)
(51) =
L(R,+~,2
An2
-
is a Jacobian matrix. Substitution of (57) into (56) yields
2R; + Rj_,,Z) i=2,3,...
Wk+’ = (I - ATJ))‘[ W” + A7(Sk - ATJW”)] ,I - 1
The spatial derivatives of (39) can be discretized the grid schematically shown in Figure 2 as R,m 112
-
RIG
(52)
(59)
in
which was solved iteratively at each time step until the following convergence criterion was satisfied: Ii2 I- 1 x AW:. Aw, iel (60) L i=2
1
I
Note that the values of 14, R, E, and m at i = I 4 are calculated by linear extrapolation as indicated in (48). Note also that (53) depends on m, L, R, u, and 0 at i = I - 4. The values of U at i = I - 4 depend on the values of U from i = $ to i = 1 - I through (46). Therefore (46) and (53) are nonlinearly coupled. Equation (53) can be substituted into (47), and the resulting values of ii can then be substituted into (46) and (53). The resulting system of equations and (53) can then be written as CYW -= S a7
(54)
where W=(V~,V~‘,...,
V,T_,,L)T
(55)
and the vector S can easily be calculated from (36), (46), and (51)-(53). Equation (55) represents a system of nonlinearly coupled, first-order, ordinary differential equations for the vector W, because the spatial derivatives have been discretized while the independent variable r has been kept continuous, and can be solved by means of ordinary differential equations integrators. In this paper, (55) was discretized by using tirstorder backward differences in time, and the resulting
-= dr
4&
where E, is a user-specified error tolerance, AW = WA+’ - Wh, and the left-hand side of (60) represents the LZ-norm of the vector AW. Note that the quasi-linearization of S”+’ (cf. equation (57)) in the Newton-Raphson method yields a sparse, but not banded, matrix, because the kth component of the vector S, i.e., Sk, depends on V,_, , VA, Vk+ , , L, V,_,, and VIPz, as indicated in (46), (47), (53), and (55). Note also that since the column vector W has 4(Z - 2) + 1 components, the dimensions of the Jacobian matrix are (41 - 7) X (41 - 7). 4.2. Nonconservative block-implicit method A method of lines technique similar to the one presented in the previous subsection can be derived by discretizing the governing equations written in nonconservation law form (cf. section 3.2) as follows. The spatial derivatives of (40) can be discretized in an equally = 1,wherel spaced grid such that 7, = Oandv, denotes the number of grid points, by means of backward differences for the advection term, and central differences for G, and the resulting system of nonlinearly coupled, first-order, ordinary differential equations for Ui, can be written as aW -= f3r
S
(61)
where w = (UT, UT,. . . ) Wl,LLT
RX
(62)
(63)
RI - RI-, -- b,* m, - ml_, A77
(57)
Arl
Appl. Math. Modelling,
1992, Vol. 16, September
469
Adaptive
methods
for annular liquid jets: J. I. Ramos kfl
The values of u, R, E, and m at i = I can be calculated by linear extrapolation as ur = 2u,_,
- u,_,
(64)
S can easily be found from the discretization of the spatial derivatives of (40) and (63), and (63) has been included in (61). It must be noted that the derivatives aRl?~r] and a*RlQ*, which appear in the vector G (cf. equation (50)), were discretized as
These values of Uk + ’ , Lk + ’ and (dL/dr)” + ’ can then be employed to evaluate C and G, and this iterative procedure can be repeated until the following convergence criterion is satisfied:
Li=l
1 l/2
Ri+l - Ri_1 + O(AT~)
2Arl
i = 2,3,.
+
. . ,I -
1 (65)
Ri+ 1 - 2Ri + Rip 1 + O(Aq*)
AT2
i = 2,3,.
. . )z - 1
(66)
and RI was calculated from (64). Equation (61) has the same form and the same number of unknowns as (54) and was solved by means of the Newton-Raphson method until (60) was satisfied.
4.3. Nonconservative block-bidiagonal method Equation (40) can also be discretized in the equally spaced grid of the previous section by means of backward differences in time, upwind differences for the advection term, and central differences for G, and the resulting O(Ar,Aq)-accurate finite difference equation can be written as -Cl+‘uy:,l = ArG:+’
+ (I + C;+l)u;+’ + U;
i=2,3,...
,I-
Appl.
Math.
(68)
Modelling,
(70)
In all the calculations presented in this paper the first iteration was always performed with the values of U, L, and dLldr corresponding to the nth time level. It must be noted that the block-implicit methods presented in sections 4.1-4.3 solve the discretized form of the fluid dynamics equations from the nozzle exit up to the grid point closest to the convergence point and that the values of the annular liquid jet’s mass per unit length, radius, and axial and radial velocity components at the convergence point are determined by linear extrapolation (cf. equations (48) and (64)). However, the fluid dynamics equations derived in Ref. 1 are not valid at the convergence point where a meniscus is formed and a pressure gradient may exist, because axial and radial pressure gradients in the annular liquid jet were neglected in the derivation of (l)-(5). It must also be noted that the grid used in section 4.1 is different from the ones used in sections 4.2 and 4.3 and that the block-implicit method of section 4.1 is conservative, whereas the methods of sections 4.2 and 4.3 are not conservative.
5.
and the subscript i denotes the ith grid point, i.e., 77i = i An, where A7 is the spatial step size. Equation (67) represents a (banded) block-bidiagonal matrix of dimensions (41 - 8) x (41 - 8), which can be easily solved by forward substitution. However, since both H and C depend on dLldr and L, and since dLldr can be evaluated as in (63) where VI can be calculated by linear extrapolation as in (64), substitution of (63) into (41), (67), and (68) would yield a sparse system of nonlinear algebraic equations, because U, depends on the value of Ui, i = 2, 3, . . . , I - 1, through (64) and (67). Such a sparse system can be solved by means of the Newton-Raphson method presented in sections 4.1 and 4.2. However, the inversion of the Jacobian matrix of the Newton method whose dimensions are identical to those of sections 4.1 and 4.2 can be avoided if one uses the following iterative technique for both U and L. Assume that Uk, Lk, and dLkldr are known at the kth iteration. These values can be used to evaluate C and Gin (67). The solution of (67), i.e., Uk+‘, can then be used in (63) to obtain (dLldr)k+’ and
470
9 10-4
1 (67)
where C = HAr/An
(Lk+’ - LA)*
1992, Vol. 16, September
Steady state and initialization
The fluid dynamics equations and the numerical methods presented in this paper are applicable to both annular liquid jets (bg # 0) and annular liquid membranes (b$ = 0) and can be used to study steady-state and time-dependent problems. Under steady-state condition, a = 0, i.e., the boundary conditions at the nozzle exit are independent of time, and dLldr = 0. In order to determine the steady-state geometry of annular liquid jets, one can simply set all the partial derivatives with respect to r equal to zero in the equations presented in sections 2-4 and solve the resulting finite difference equations until the convergence criteria deIined in (60) and (70) are satisfied. Alternatively, one can solve the time-dependent equations subject to time-independent boundary conditions at the nozzle exit until an asymptotic steady state is reached. This time-dependent approach requires that (60) and (70) be satistied at each time step, and the asymptotic steady solution is reached whenever the following steady-state convergence criteria are satisfied: I-‘AWiT. AW, ‘I*
Iz
i=2
Ar2
5
1
E*
(71)
Adaptive methods for annular liquid jets: J. I. Ramos rA (u?+l - u;)r.
t 1;=
(u:+l
I
where AW=
_ uy) + (L”+I
A72
W”+’ - W”
_ ~721i’~ zs
exit angle, and amplitude and frequency of the fluctuations. In order to assess the influence of these parameters on the dynamic response of liquid jets, parametric studies were performed to determine the convergence length and the pressure coefficient as functions of time. The results of some of these parametric studies are illustrated in Figures 3-13, which show the (nondimensional) axial velocity component at the nozzle exit (u,, = u(T,O)), the pressure coefficient, and L’ = L/L(O) for the values of the parameters shown in Table I and b$ = 0, i.e., annular liquid membranes. The values of L(0) are equal to 12.5590 for Figures 3-5, II, and 12; 9.8652 for Figure 6; 0.3262
(73)
lZ
and c3 are user-specified error tolerances. In the calculations presented in this paper, cl = l2 = l3 = 10-4. Equation (71) was used with the methods of lines of sections 4.1 and 4.2, while (72) was employed with the block-bidiagonal technique of section 4.3. The time-dependent approach, which may be used to determine the asymptotic steady state, may require shorter computational times than the solution of the steady-state governing equations, for the discretization of the time-derivative terms increases the diagonal dominance of the matrices to be inverted when solving the fluid dynamics equations (cf. equation (67)). In this paper the steady-state fluid dynamics field was determined asymptotically from the solution of the time-dependent equations. Once steady state was achieved, 7 was reset to 0, and m,(O) = (rn;).$.,and p;(O) = (PJSS, where ss denotes steady state. These initial conditions were then used to analyze the response of annular liquid jets to fluctuations in the axial velocity at the nozzle exit in the absence of mass injection into the volume enclosed by the annular liquid jet. 6.
Presentation
(72)
c3
1
Ucl
0.6
0.2
of results
Numerical experiments were performed to determine the efficiency of the methods presented in sections 4.1-4.3 and the effects of both the time step and number of grid points on the response of annular liquid jets subject to mass flow rate fluctuations at the nozzle exit. The results of these experiments indicate that the techniques presented in sections 4.1-4.3 essentially yield the same results and that the iterative, blockbidiagonal method of section 4.3 is more efficient than the block-implicit techniques of sections 4.1 and 4.2 because these techniques require the inversion of the Jacobian matrix, whereas the method of section 4.1 involves a block-bidiagonal matrix and the vector U can easily be found by forward substitution (cf. equation (67)). The number of grid points, I, and the time step, AT, used in the calculations were varied in order to ensure grid independence. For all the calculations presented in this paper, 400 5 I 5 1600 and 0.04 5 AT 5 0.001, and the largest value of I corresponded to annular liquid jets whose steady-state convergence length was longest. The smallest value of Ar corresponds to the highest value of the Strouhal number employed in the calculations. As shown in section 2, the dynamic response of inviscid, isothermal, annular liquid jets to velocity fluctuations at the nozzle exit depends on the Froude and Weber numbers, C,,,, , initial pressure ratio, nozzle
L’
1.0
-0.2
1 0
IO
20
30
40
Figure 3. Convergence length, pressure coefficient, velocity at the nozzle exit as functions of time
50
’
and axial
L’
1.0
Uo 0.6 -
-0.2 1 0
10
20
30
40
Figure 4. Convergence length, pressure coefficient, velocity at the nozzle exit as functions of time
Appl.
Math. Modelling,
50
’
and axial
1992, Vol. 16, September
471
Adaptive
methods
for annular liquid jets: J. I. Ramos
20
10
30
40
Figure 5. Convergence length, pressure coefficient, velocity at the nozzle exit as functions of time
50
1
-0
and axial
10
20
30
50
40
Figure 8. Convergence length, pressure coefficient, velocity at the nozzle exit as functions of time
’
and axial
r 0.6
0.6
0.2
0.2
-0.2
/‘\%
1
0
10
I
20
30
40
Figure 6. Convergence length, pressure coefficient, velocity at the nozzle exit as functions of time
50
’
-I
-0.2 0
CPll
10
20
30
40
Figure 9. Convergence length, pressure coefficient, velocity at the nozzle exit as functions of time
and axial
1.0
1.0
50
’
and axia
I!
"0
“0
0.6 .
0.6
0.2 /_
-0.2
’
10
0
20
I 30
40
Figure 7. Convergence length, pressure coefficient, velocity at the nozzle exit as functions of time
472
Appl.
Math. Modelling,
50
1
and axial
1992, Vol. 16, September
Figure 10. Convergence length, pressure coefficient, velocity at the nozzle exit as functions of time
and axial
Adaptive
methods
for annular liquid jets: J. I. Ramos
for Figure
7; 19.0452 for Figure 8; 3.8236 for Figure 9; 37.3466 for Figure 10; and 10.2816 for Figure 13. Figure 3 indicates that the convergence length and
1.0
the pressure coefficient do not respond instantaneously to the velocity fluctuations at the nozzle exit. The delay times in both the convergence length and the pressure coefficient are functions of the parameters that govern the dynamics of liquid jets, as discussed in the next paragraphs. Figure 3 also indicates that the convergence length and the pressure coefficient are periodic functions that have the same frequency as the axial velocity fluctuations at the nozzle exit, that the pressure coefficient has a smaller delay time than the convergence length, and that there is a lag between the peaks of the pressure coefficient and those of the convergence length. Figures 3 and 4 show that the delay and lag times in the response of both the convergence length and the pressure coefficient are nearly independent of the amplitude of the axial velocity fluctuations at the nozzle exit, whereas the amplitude of the oscillations of both the convergence length and the pressure coefficient increases as the amplitude of the axial velocity fluctuations at the nozzle exit is increased. Figures 3 and 4 also show that the locations of the maxima and minima of the pressure coefficient are nearly independent of the amplitude of the axial velocity fluctuations at the nozzle exit. However, striking differences are observed in the convergence lengths and pressure coefficients shown in Figures 3 and 4. In particular, Figure 3 indicates that both L and C,,, are nearly sinusoidal, whereas the values of L and C,, shown in Figure 4 deviate considerably from the sinusoidal behavior of the axial velocity fluctuations at the nozzle exit. This result is a consequence of the fact that for small amplitudes of the axial velocity component at the nozzle exit the response of liquid jets is nearly linear, as illustrated in Figure 3; for large amplitudes, however, nonlinear effects are important, and the response of the liquid jet is nonlinear. Figures 3 and 5 illustrate the effects of the Strouhal number on both the convergence length and the pressure coefficient and indicate that the amplitude of the oscillations of both L and C,, decreases, whereas their frequency increases, as the Strouhal number is increased. These figures also show that the delay time in the convergence length decreases as the Strouhal number is increased. Figure 5 indicates that for St = 0.5 the pressure coefficient is nearly zero. The effect of the Froude number on the dynamic response of liquid jets is illustrated in Figures 3 and 6, which indicate that the delay time in the convergence length increases as the Froude number is increased. Figures 3 and 6 also indicate that both the convergence length and the pressure coefficient exhibit a behavior that increasingly deviates from the sinusoidal one exhibited by the axial velocity fluctuations at the nozzle exit as the Froude number is increased. This is due to the fact that the values of the convergence length at steady state decrease as the Froude number is increased. Note that since the steady-state convergence
0.6
I:lIYWL 0.2
0
10
20
30
40
50
t
Figure 11. Convergence length, pressure coefficient, and axial velocity at the nozzle exit as functions of time
II I.0
, "0
0.6
0.2
-0.;
-0.t -0
10
20
30
40
50
1
Figure 12. Convergence length, pressure coefficient, and axia velocity at the nozzle exit as functions of time
L 1.0 “0
0.6
0.2
-0.2
-0.6
% 0
10
20
30
40
Figure 13. Convergence length, pressure coefficient, velocity at the nozzle exit as functions of time
50
I
and axial
Appl.
Math.
Modelling,
1992,
Vol. 16, September
473
Adaptive
methods for annular liquid jets: J. I. Ramos Table 1.
Values
of the parameters used in the calculations P?(O)
Figure
Ff
We
00
Cpmax
P,”
a
St
3 4 5 6 7 8 9 10 11 12 13
10 10 10 10,000 10 10 10 10 10 10 10
50 50 50 50 1 100 50 50 50 50 50
0 0 0 0 0 0 -15 15 0 0 0
1 1 1 1 1 1 1 1 2 10 1
1 1 1 1 1 1 1 1 1 1 0.5
0.1 0.25 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1
0.1 0.1 0.5 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1
length is a measure of the volume enclosed by the liquid jet, the deviation from a sinusoidal behavior increases as the volume enclosed by the jet at steady state decreases. Figures 3, 7, and 8 show the effects of the Weber number on both the convergence length and the pressure coefficient, and they indicate that the delay times in both L and C,,, decrease, whereas the amplitudes of both the convergence length and the pressure coeflicient increase, as the Weber number is decreased. Figure 8 also shows that the steepness of the convergence length from its maximum to its minimum values increases as the Weber number is increased, i.e., as the convergence length at steady state and the inertia of the liquid jet are increased. Figures 3, 9, and 10 illustrate the effects of the nozzle exit angle on the dynamic response of liquid jets, and they indicate that the delay times and the amplitude of both the convergence length and the pressure coefficient increase as the nozzle exit angle is increased. This result is consistent with the steady-state convergence lengths, which increase at the nozzle exit angle is increased,‘,3.4 and indicates that the delay times decrease as the volume enclosed by the jet as steady state is decreased. Figure 10 also shows that the steepness of the convergence length from its maximum to its minimum values increases as the nozzle exit angle is increased. This result is also consistent with the results shown in Figure 8 and indicates that the steepness of the convergence length increases as the convergence length at steady state is increased. Figure 10 also shows the dynamic response of the liquid jet has not reached a periodic behavior even for times greater than 40. The effect of C,,,,, which is a measure of the importance of surface tension and pressure of the gases that surround the liquid (cf. equation (31)), on the dynamic response of liquid jets is illustrated in Figures 3, 11, and 12, which indicate that the delay time in both the convergence length and the pressure coefficient decreases as C,,,, is increased. The amplitudes of both the convergence length and the pressure coefficient decrease and increase, respectively, as C,,,, is increased. It must be noted that the value of C,,,, does not affect the convergence length at steady state for the
474
Appl.
Math.
Modelling,
1992,
Vol. 16, September
values of the parameters presented in Table 1. However, C,,,, P lays the role of an amplification factor, for it increases the pressure coefficient if there is a pressure difference between the gases enclosed by and surrounding the liquid jet (cf. equation (30)). Figures 3 and 13 illustrate the effects of the initial pressure ratio on the dynamic response of liquid jets, and they indicate that the delay times of both the convergence length and the pressure coefficient decrease as the underpressurization of the jet is increased, i.e., as the convergence length at steady state is decreased. 7.
Conclusions
A domain-adaptive technique and three block-implicit methods have been used to analyze the dynamic response of inviscid, isothermal, annular liquid jets subject to velocity fluctuations at the nozzle exit as a function of the Froude and Weber numbers, initial pressure ratio, nozzle exit angle, surface tension, and amplitude and frequency of the imposed velocity fluctuations. The domain-adaptive technique maps the timedependent, curvilinear geometry of the jet into a unit interval and yields a system of integrodifferential equations for the annular liquid jet’s mass per unit length, mean radius, and axial and radial velocity components. The integrodifferential character of the equations is due to both the equation that governs the location of the time-dependent convergence point and the condition that the mass of the gases enclosed by the annular liquid jet is constant. Two of the block-implicit methods presented in this paper use a method of lines technique, which discretizes the spatial derivatives and keeps continuous the time variable, and yield a system of nonlinearly coupled, first-order ordinary differential equations for the convergence length and the mean radius, mass per unit length, and axial and radial velocity components of the annular liquid jet. This system of ordinary differential equations has been solved by means of the Newton-Raphson method, whose Jacobian matrix is sparse but not banded. The third block-implicit method is an iterative technique, which discretizes both the spatial and the temporal derivatives and yields a block-bidiagonal matrix for the mean radius, mass per unit length, and axial
Adaptive
and radial velocity components of the annular liquid jet. This block-bidiagonal matrix is solved iteratively by forward substitution together with the ordinary differential equation for the convergence length. It has been shown that the block-bidiagonal method is more efficient than and yields the same results as the blockimplicit techniques presented in this paper. It has also been shown that both the convergence length and the pressure coefftcient are periodic functions of time when the liquid jet is subject to sinusoidal velocity fluctuations at the nozzle exit and that there is a delay time in the dynamic response of the jet. The delay time is nearly independent of the amplitude of the velocity fluctuations at the nozzle exit, and it decreases as the frequency of these fluctuations is increased. The delay time in the convergence length increases as the Froude and Weber numbers and nozzle exit angle are increased and decreases as the underpressurization of the liquid jet is increased. It is also shown that there is a lag between the response of the convergence length and that of the pressure coefficient and that both the delay time and the lag time depend on the inertia of the annular liquid jet. The convergence length and the pressure coefficient are periodic functions of time and have the same frequency as that of the axial velocity fluctuations at the nozzle exit, because the liquid jets analyzed in this paper are inviscid and isothermal. When the amplitude of the velocity fluctuations at the nozzle exit is small, the convergence length and the pressure coefficient are nearly sinusoidal, and the liquid jet behaves almost linearly. For large amplitudes, however, nonlinear effects are important, and the response of the jet is periodic but not sinusoidal. It has also been found that the steepness of the convergence length from its maximum to its minimum values increases as the Froude and Weber numbers and nozzle exit angle are increased. An increase in the Froude number corresponds to a decrease in both the force of gravity and the steady-state convergence length, whereas the convergence length at steady state increases as the Weber number and the nozzle exit angle are increased. The results shown in this paper indicate that the dynamic response of liquid jets to axial velocity fluctuations at the nozzle exit is very much dependent on the magnitude of the gravitational force and on the volume enclosed by the annular liquid jet at steady state.
methods for annular liquid jets: J. I. Ramos
Acknowledgments The calculations presented in this paper were performed at the Facultad de InformaticaiEscuela Tecnica Superior de Ingenieros de Telecommunication of the Universidad de Malaga. The author is grateful to Mr. Juan Falgueras for his help with both the graphs and the typing of this paper. References
5
12 13 14
15 16
17
18
Appl.
Ramos, J. I. Liquid curtains, I, Fluid mechanics. Chem. Eng. Sci. 1988.43, 3171-3184 Ramos, J. I. and Pitchumani, R. Liquid curtains, II, Gas absorption. Chem. Eng. Sci. 1990, 45, 1595-1604 Ramos, J. I. Annular liquid jets in zero gravity. Appl. Math. Modelling 1990, 14, 630-640 Ramos, J. I. Analytic, asymptotic and numerical studies of liquid curtains and comparisons with experimental data. Appl. Math. Modelling 1990, 14, 170-183 Ramos, J. I. Dynamic response of liquid curtains to time-dependent pressure fluctuations. App[. Muth. ModeDing 1991, 15, 126-135 Binnie, A. M. and Squire, H. B. Liquid-jets of annular cross section. Engineer (London) 1941, 171, 236-238 Baird. M. H. I. and Davidson, J. F. Annular iets. I. Fluid dvnamics. Chem. Eng. Sci. 1962, 17, 467-472 ” Baird, M. H. I. and Davidson. J. F. Annular jets. II, Gas absorption. Chem. Eng. Sci. 1962, 17, 473-480 Roidt, R. M. and Shapiro, Z. M. Liquid curtain reactor. Report No. 85M981, Westinghouse R&D Center, Pittsburgh, PA, 1985 Lee, C. P. and Wang, T. G. A theoretical model for the annular liquid jet instability. Phys. Fkds 1986. 29, 2076-2085 Hoffman, M. A., Tabahashi, R. K., and Monson, R. D. Annular liquid jet experiments. ASME J. Fluids Eng. 1980, 102, 344-349 Lee, C. P. and Wang. T. G. A theoretical model for the annular liquid jet instability-Revisited. Phys. Fhrids A 1989. 1, 967-974 Meseguer. J. The dynamics of small annular jets. ASME J. FINid. Eng. 1988, 110, 123-126 Ramos, J. I. Liquid membranes: Formulation and steady state analysis. Report No. CO/89/4, Dept. of Mech. Eng., Carnegie Mellon Univ., Pittsburgh, PA. 1989 Ramos, J. I. Annular liquid jets: Formulation and steady state analysis. Z. Anger. Moth. Mech. 1991, in press Chigier, N., Ramos. J. I., and Kihm, K. Experimental and theoretical studies of vertical annular liquid jets. Proceedings c?fthe Sixth Symposium on Energy Engineering Sciences: Flow und Trunsport in Continua, U.S. Department of Energy, Argonne National Laboratory, Argonne, IL, 4-6 May 1988, pp. 18-25 Ramos, J. I. Adaptive finite difference methods for liquid membranes, I. Mass flow rate fluctuations at the nozzle exit. Report NO. CO/91/1, Dept. of Mech. Eng., Carnegie Mellon Univ., Pittsburgh, PA, January 1991 Ramos, J. I. Adaptive finite difference methods for liquid membranes. II, Mass injection. Report No. CO/91/2, Dept. of Mech. Eng., Carnegie Mellon Univ., Pittsburgh, PA, January 1991
Math. Modelling,
1992, Vol. 16, September
475