- Email: [email protected]
A computer simulation of turbulent jet flow
Compeersand Fledds Vol.9.pp.205-221
0045-793018110601-02051502.0010
PergwonPressLid.,1981.PrintedinGreatBritain
A COMPUTER
SIMULATION OF TURBULENT A. V. J. EDWARDSand C. L. MORFEV
JET FLOW
Institute of Sound and Vibration Research, University of Southampton, Southampton, England
(Received 14 December 1979) Al~traet--Tbe evolution of a round turbulent jet is simulated numerically under the constraint of axial symmetry. The vortex sheet shed from the orifice is represented by vortex ring elements, with a velocity field cut-off to control close encounters. Large-scale vortex clusters form in the model jet, similar to those inferred from laboratory flow visualization experiments. However, comparisons of statistical properties reveal significant differences between the axisymmetric model flow and real turbulent jets. 1. I N T R O D U C T I O N
Flow visualization experiments on plane turbulent mixing layers[l] and round turbulent jets [2] have revealed the occurrence of fairly regular large-scale features, commonly described as orderly structure. The typical length scale (wavelength) is around 0.3x, where x is the downstream distance from the start of the shear layer. Such features have been observed over a range of Reynolds numbers Re~ of order 3 x 104-3 x 10613]. We demonstrate in this paper that a similar orderly structure emerges from a highly idealized model jet flow in which two of the three vorticity components are suppressed. The model flow consists of an axisymmetric starting jet, which issues from a circular opening of radius R0 in a rigid plane boundary (Fig. 1). The normal velocity over the opening jumps, at t = 0, from zero to a constant value U0, and an axisymmetric vortex sheet is shed from the lip of the opening. The jet and ambient fluids are assumed inviscid and incompressible, with identical density pp. Instability of the vortex sheet leads to a chaotic, quasi-turbulent (but still axisymmetric) flow field, with some qualitative resemblance to a real turbulent jet. The evolution of the large-scale orderly structure in the model is described in Section 3, and statistical properties of the model jet (as it approaches a steady state) are discussed in Sections 4 and 5. We begin in Section 2 by describing the calculation method in some detail; for convenience all quantities will be normalized in what follows to unit values of p0, R0 and Up. The objective of the study is not to demonstrate that certain statistical features of turbulence in jets can be quantitatively reproduced in an axisymmetric model: such agreement is not to be expected. Rather, the aim is to set up a dynamically self-consistent axisymmetric model, to observe its statistical properties and to contrast these with real jet measurements. The axisymmetric jet model, thus, becomes a tool for studying the large-scale consequences of complicated vortex motions.t Numerical flow experiments of a similar type--but directed at the small-scale aspects of turbulent motion--have previously been described by Batchelor[4].
\ Source distribution to
/ /
:
~'~.~ ~ , ~
maintain boundary condition over opening
d••m
~ n R i n g vortices ring vortex sheet
--~ Image vortices to maintain boundary condition over solid surface ~
f
Uniform flow of velocity U
!/ / /
Fig. I. Jet flow simulation using sources and vortices. tFurther results dealing with the local pressure and radiated sound (including comparisons with experiment) will be reported in a later paper. The far-field sound predictions are summarized in Ref. [5]. 205
206
A.V.J. EDWARDSand C. L. MORFEY 2. THE L U M P E D - V O R T E X FLOW R E P R E S E N T A T I O N
We follow Davies et a/.[6], Chorin and Bernard [7] and Clements[8] and represent the continuous vortex sheet by an array of discrete elemental ring vortices.'t Each vortex element moves under the influence of all the others, including the set of image vortices reflected in the orifice plane (Fig. I). The prescribed normal velocity U0 over the opening is simulated by a source disk of appropriate strength. An unwanted consequence of discretizing the vortex sheet is the appearance of unrealistic induced velocities close to each vortex element. In particular, this causes problems when two vortices come close together. The need for some kind of smoothing of the discrete-vortex velocity field was pointed out by Chorin and Bernard [7]; they proposed a cut-off, or upper limit, for the velocity induced by each vortex, and this has been adopted in the present study. The precise manner in which the cut-off is applied turns out to be important; we return to this point later. 2.1 The shedding process In general, if the fluid crosses the orifice plane with a uniform normal velocity U(t), the assumption that a vortex sheet leaves the lip parallel to the axis implies a "rate of circulation shedding Ou F= - f o u-z-dyoy
:-21 U2(t)"
(1)
Here cylindrical coordinates (x, y, ~b) are being used, and the circulation is in the ~b direction. It follows from (1) that if the circulation shed in the interval (t, t + 6t) is lumped into a ring vortex, the circulation around the ring is F = ~U2(t)6t + O(6t) 2. Likewise the centroid of the circulation shed in this interval is located at (x, y) = (U(t)St/3, 1) at the end of the time interval, apart from O(6t) 2 terms.:~ Thus, the strength and location of each vortex element as it is first shed are determined by the efflux velocity U(t). 2.2 The motion of a vortex element An individual vortex element is subjected in general to three superimposed velocity fields. These are the velocity induced by all the other vortex elements (including images); the velocity induced by the vortex itself; and the velocity field of the source disk. (a) The velocity field of other vortices Associated with each elemental vortex ring there is a velocity field, whose axial and radial components at a point (x, y) are given by u = FR(RJ - yK),
(2)
v = FRxK.
(3)
Here F is the circulation of the elemental vortex, R is the radius of the elemental vortex, and (x, y) are axial and radial coordinates with the centre of the vortex, ring as origin. The quantities J and K are defined by
(4)
1 E(k) J = lrr23 d - F)' -1 1 [ K --- 7rr23~ 2K(k)
(2-ke)E(k)} (1 k 2)
"
(5)
Here K(k), E(k) are complete elliptic integrals of the first and second kind respectively, k 2 = 1 - rm2/r22,and (rm, rE) are the distances of the nearest and farthest points on the vortex ring from the point (x, y). Equations (4) and (5) are equivalent to expressions given by Lamb[12]. tFor a review of the discrete vortex calculation method, see Refs. [9-11]. ~The vortex sheet is assumed to leave the orifice axially. A detailed derivation appears in Ref. [22].
A computer simulation of turbulent jet flow
207
(b) The self-induced velocity of a vortex ring A ring of zero cross section has an infinite self-induced velocity, and a better approximation is clearly called for. We recognize that initially each ring represents a tubular segment of vortex sheet. Calculating the self-induced velocity of the centroid of such a segment gives F /1 +
16',
(6t
(6)
This expression strictly applies only so long as the sheet remains cylindrical. It is best regarded as a correction term which tends to zero as 6t ~ 0, but which appears worth including in the interest of greater accuracy close to the orifice, where the sheet begins to break up. (c) The velocity lield of the source disk Two methods are available for meeting the exit-plane boundary condition on normal velocity. One is to use image vortices (as in Fig. 1) together with a source disk of uniform strength 2U(t) per unit area. The alternative is to dispense with image vortices (thereby saving considerable computing time) and to use a non-uniform source disk of strength 2 U ( t ) 2uv(y, t), where uv is the normal velocity induced on the exit plane by the real vortices in x > 0. In practice a combination of these methods was used. Up to t = 15, with relatively few vortex elements present, the image method was used for maximum accuracy. Beyond t = 15, an approximate version of the second method was used: rather than updating the source distribution at each step, use was made of the fact that uv varies very little with time after about t = 10. Thus, an average value, taken over the interval 10 < t < 15, was inserted in the source strength expression. The source in the orifice plane is represented by a series of concentric rings of appropriate strengths. Obviously, from the point of view of the economical use of computer time the smaller the number of rings used to represent the continuous source distribution the better. However, using too few rings gives rise to errors near the orifice, in particular at the initial shedding location. A small study led to the conclusion that 30 rings inside the orifice and ten more outside constituted a reasonable representation of the required source distribution. It is necessary to include a number of more widely spaced rings (sinks rather than sources) outside the orifice, in order to maintain the boundary condition u = 0 for y > 1 when the vortices are limited to the x > 0 half space. 2.3 The need for a velocity cut-off Discrete vortex representations of flows with continuously distributed vorticity produce unrealistic velocity fields close to the vortex elements. In particular, this causes problems when two vortices come close together. To ensure a smooth roll-up to the starting vortex, without the strong mutual repulsions otherwise associated with close encounters, we follow Chorin and Bernard[7] in applying a cut-off or upper limit to the induced velocity associated with each vortex. In the early stages of the present work, and in a similar study of a starting jet by Acton[13], the importance of conserving momentum during the cut-off process was overlooked. As a result, spurious vortex accelerations were introduced whose effect on the calculated pressure field was particularly serious. Some of the results to be presented:~ were obtained with the momentum-conserving cut-off factor derived in[14], namely
f(rO = rl2/rc 2,
(rl < rc = cut-off radius);
f(rD = 1,
(rl <- rc).
(7)
The above weighting factor is applied to both the radial and axial velocity components induced tNote that t represents the nondimensional time (normalized by RolUo). Since F - 6t as 8t tends to zero, the log term is typically a small correction. ~Speciflcally, those with 8t = 0.2.
208
A.V.J. EDWARDSand C. L. MORFEY
by either member of the vortex pair on the other. However, many of the earlier results are included for comparison, particularly where momentum conservation did not prove critical. The choice of r~ is a compromise between the conflicting requirements of a smooth initial vortex roll-up on one hand, and the need to minimize the invocation of cut-off on the other, as discussed in Ref. [14]. The results presented are based on r,. = 0.05.
2.4 The integration scheme (a) Linear multistep procedure (Adams-Bash[orth ) Suppose that the velocity u(t) is fitted in the vicinity of t = 0 by the order k polynomial k+l
~(t) = ~ aitH;
(8)
i=1
then to determine the (k + I) coefficients ai we must equate ti(t) to u(t) at (k + 1) points (times). These may be taken as t = 0 , - 1 . . . . . - k where the time scale is normalized on the interval between successive points. The resulting (k + 1) equations take the form k+¿
Zai(1-j)i-I=ui,
]=1,2 ..... k+l
(9)
i=l
where uj represents the value of u taken ( j - 1 ) intervals back from t = 0. We use the term "back-velocity" to refer to the velocity of a vortex at some previous time-step. Finally, (xl - x0) = Ax is estimated from eqn (8) as 1
Ax-'
k+l
y0
k+l
(10)
t~(t)dt = ~•: a_~_ i - ~ b~u~,say.
The coefficeints bi are listed below for orders up to k = 5.
k
Dbl
Db2
0 1 2 3 4 5
i 3 23 55 1901 4277
-1 - 16 -59 -2774 -7923
Db3
5 37 2616 9982
Db4
-9 - 1274 -7298
Db5
251 2877
Db6
-475
D 1 2 12 24 720 1440
To obtain bi, the values in each row must be divided by the corresponding number D in the final column.
(b) Estimation of the correction order required In the computation of Ax from eqn (10) one may not invoke the use of an arbitrarily high correction order without suitable justification. In the present case we have, in fact, a highdensity classical many-body problem so that Ax is of the order of the separation of the elemental vortices. One might suspect that use of more than a very few back-velocities of a particular vortex in the calculation of its subsequent position would be erroneous, since one would be imposing continuity upon velocity derivatives calculated from greatly differing past configurations. In other words, velocity derivatives of too high an order become spuriously large causing eqn (10) to diverge. It was found that in the present case the use of more than one back-velocity gave rise to spurious behaviour of the elemental vortices. Accordingly k was taken as I, and vortex
209
A computer simulation of turbulent jet flow
coordinates at each step were calculated from 1 x.+,(i) = x~(i)+ ~(3u~(i)- u._,(i)) fit. 1
y.+,(i) = y.(i)+ ~(3v.(i)- v~_,(i)) 8t
(11)
where u.(i) and v.(i) are, respectively, the axial and radial velocities of the ith vortex at the nth time increment. 3. V I S U A L
CHARACTERISTICS
OF THE
SIMULATED
FLOW
In this section we briefly describe the model starting jet in terms of the evolving vorticity distribution as viewed in a plane containing the jet axis. The development of the simulated flow (from time t = 0) may be divided into three fairly distinct phases: (i) The formation of a starting vortex. (ii) The detachment of the starting vortex. (iii) The steady-state situation. Figure 2 shows the state of the simulated flow at four different times, calculated with 8t = 0.1. The same calculation using 8t = 0.2 yields similar vortex patterns (Fig. 3) despite differences of detail; it appears that large scale features of the vorticity distribution are not sensitive to St. The first frame (t = 6.2) shows the starting vortex about to detach, leaving a trail of smaller vortex clusters which show up clearly in the second frame. Thereafter, small clusters form at fairly regular intervals near the orifice, and grow by amalgamation with neighbouring clusters as they are swept downstream. Eventually the flow in any given region downstream achieves a statistically steady state, in which the general appearance of the vortex pattern no longer changes with time. This occurs when the starting vortex has passed well beyond the region concerned. The formation in real turbulent jets of similar vorticity concentrations, and their growth in scale (both size and separation) with increasing distance downstream, have been described by Brown and Roshko[l], Damms and Kficbemann[15] et al. The fact that similar features appear TIME
6.2000
61 VORTICES
TIME 16.2000
~.'~ •
~ .,~ .~. r,:.
,
TIME 26.2000
-,, ,.¢.
.~,~.
,
..,~-..
,~.'.~
,
246 VORTICES
..
154 VORTICES
~:~: • '.:a'~-
,
TIME 36.2000
..
,
358 VORTICES
" . "~-~'-..
." ::A~',..::.,,', ".~..
Fig. 2. Evolution of vortex clusters in axisymmetric jet, calculated using 8t = 0.1. The lower half of each frame contains the instantaneous velocity profile of the radial component, v, along the lip-line (y = 1.0).
210
A. V. J. EDWARDSand C. L. MORFEY TIME 6.2000
31 VORTICES
81 VORTICES
TIME 16.2000
....... :.f!:i::.
vv
TnlE 26.2000
151 VORTICES
:~ #, ,.:,,.. ~.....:,:~..:"
...,.
A A
, ~
•
: " ;.-...;.
~:~ :" ::
~
:
W
TIME 36.2000
. . . . . . . . . . . .
,
,
;" :~ "' ....:W:'::'.
A, A ~
181 VORTICES
:..?:..,.
•
.:
" '":: ..:'":...
,X~A\/;A
:...
...i.:..i:q
...~.:::
x~i//'-~x
Fig. 3. As Fig. 2, but calculated using ~St= 0.2 (momentum conserved during cut-off).
in the model jet is the main justification for studying the model flow in more detail; we can exploit the accessibility of information from the numerical model to provide insights into the large-scale aspects of the real three-dimensional flow. 4. S I N G L E - P O I N T S T A T I S T I C S OF THE V E L O C I T Y F I E L D
In this section we examine the statistics of the time-dependent velocity at various points in the model flow. For this purpose the initial development of the flow from rest is ignored, and statistical estimates of mean velocity, variance, spectra, etc. are based on the steady-state part of the time-history. Comparisons are made between model flows with different time increments St, and also between the axisymmetric model and real laboratory jets. 4.1 The amplitude density of u(t) and its moments Figure 4 provides a visual comparison of u(t) time histories for two model realizations (St =0.1,0.2) and a turbulent jet in the laboratory[16]. The position in the flow is (x = 4, y =
0.8), which lies in the initial turbulent shear layer of the real jet. The two model time-histories, though qualitatively similar, exhibit detailed differences especially in the low-velocity regions: here trace (a) is unrealistically rounded, while trace (b) shows more of the sharp downward spikes present in the real flow (c). The difference between traces (a) and (b) is attributed partly to the fact that momentum was not conserved during the cut-off process in the earlier calculation (a)--see Section 2.3. Perhaps of greater importance, however, is the fact that convergence of the flow statistics at this location downstream (x = 4) appears to require a time increment of 8t = 0.1 or less: thus trace (b)--with 8t = 0.2--retains some influence of step size and does not truly approximate the limit 8t ~ 0. Further downstream, beyond x = 6, this problem disappears (as noted below). Mean velocity profiles ~(y) are shown in Fig. 5. Figure 5(a) tests the model flow profiles for similarity in the first 5 dia. downstream; it can be seen that similarity is achieved only beyond about 2 dia. in this case (St = 0.1). With smaller and larger time increments (St = 0.044, 0.2) the distance downstream at which similarity is achieved is roughly 1 and 3 dia., respectively. These results suggest that the finite step size (with its associated cut-off) provides an
A computer simulation of turbulent jet flow I
!
I
211 i
(a)
(b)
8
0
I0
20
30
40
T
Fig. 4. Axial velocity time histories at (x =4, y =0.8). (a) Model, calculated using 8t = 0.1. (b) Model, calculated using 8t = 0.2 (momentum conserved during cut-off). (c) Measurement in turbulent jet [16], at Rex = 0.4.10 6.
artificial viscosity, and leads to a "transition region" close to the orifice similar to that observed on real jets at finite Reynolds numbers (Uodh, ~ 104- I0~). In Fig. 5(b) a typical model profile (at x = 4) is compared with the profile measured at the same location in a turbulent jet[17]. The agreement is surprisingly good considering the axisymmetric nature of the model; this point is taken up in Section 6, where entrainment is discussed. Figures 6-8 show profiles of the standard deviation (u,~s), skewness (/~,) and kurtosis (K,) of the axial velocity u(t), at positions up to 5 dia. downstream of the orifice. In this region, real jet data collapse on the similarity coordinate r/= (y - 1)/x; accordingly 7/is used as the abscissa in all these plots. Each of the three model flow realizations (St =0.044,0.1,0.2) gives peak values of u,,,,s (axial turbulence intensity) around 25-30%,t i.e. roughly double what is measured in real jets[17]. This difference may be attributed to the greater orderliness imposed by the axisymmetry constraint in the model; corresponding differences are found in the local fluctuating pressure and the far-field radiated sound[5]. The skewness profiles in Fig. 7 show some significant differences between the 8t = 0.1 and 8t = 0.2 model flows, the latter corresponding more closely with real jets[18]. The differences appear to be due in part to the different cut-off procedures mentioned above, since they persist out to x = 10 where the effect of the time increment is no longer evident in other statistics (e.g. Fig. 6). tExcept where convergence has not been achieved, e.g. at x = 4 with 8t = 0.2.
212
A. V. J. EDWARDSand C. L. MORFEY
1.20
I
51a)
i.oo o.oo
°°.° F a~
-020 L
i
°°~.eo
~
-o.4o
-o.zo
o.oo (Y-I)/X DT-O.i
o2o
i
0.4o
o.eo
I.OO •
a m
cteo 0
[]
5(b)
•
0.6O
I=
=Z 0.40
0.20 0.00
m
i
A
-0.20 0.00 0120 o.~ o~ o'Jo ,ioo ,.2o' ,.~ i.~
•
,.~' 2'oo z.m
Y/RO
X/ROt4.0 Fig. 5. Mean axial velocity profiles. (a) Similarity plot of model results, calculated using 8t =0.1. Values of x are 2, 4, 6, 8, 10 with smallest x corresponding to largest range of abscissa. (b) Comparison of model and experimental velocity profiles at x =4. Code: [~, model (St =0.1); --, experiment[17], at Re~ =0.35' 106 (d = nozzle dia.).
Similar remarks apply to the kurtosis profiles of Fig. 8, where again the model results-despite the similarity collapse--bear only a modest resemblance to real jet data[18]. 4.2 Power spectra and covariances o/ u(t) and v(t) Figure 9 shows the spectral distribution of the u(t) variance with respect to Strouhai number (S = [dlUo), for three different model flow realizations sampled at (x = 10, y = I). If the three signals were statistically identical and the sample lengths had been equal, 90% of the spectral estimates would be expected to lie within a band 6 dB wide; this is a measure of the statistical accuracy of the results, and follows from the. use of 7 overlapping data segments for analysis. However, the usable sample lengths for the three runs (St = 0.044, 0.1,0.2) were different (T = 21,84, 112, respectively); therefore the resolution of the 8t = 0.044 spectrum is relatively coarse, which could account for the apparent bias of this spectrum towards higher Strouhal numbers. (The spacing AS of the spectral estimates is equal to 8IT.) The different sample lengths arose from a reduction in total flow duration To (for economy reasons) as 8t was decreased. The computing time for each run of the flow model varied
A computer simulation of turbulent jet flow
213
u
rmi
a
0.3
0.2 ÷
"4"11
4au
0.1
+ + ÷
++
++j
I
I
-0.4
-0.2
0
I
0.2
0.4
"I] = ( Y - l ) / X
Fig. 6. Profiles of standard deviation u,~s (axial turbulence intensity), plotted against similarity coordinate 7. Model data (x= 10): O, 8t =0.0M; II, 8t =0.1; A, 8t =0.2 (momentum conserved during cut-off). Experimental data (x = 8): +, from reference [17], at Rex = 1.4.106. 2
i
Bu
I
I
I
I
t
(°) I IDA
0
II
- I r-
-2
t -0.2
J3u
2
I
I
-0.1
I
0
0,1
I
(b)
I
1
I
0,2
I
0° in •
0
6 ,L
•
i
-1
-2
i -0.2
I
-0.I
I
I
0 TI = ( Y - l ) / X
I
0.I
I
0.2
Fig. 7. Skewness profiles for axial velocity, plotted against similarity coordinate 7/. Code: triangles, x = 4; circles, x = 6; squares, x = 10. (a) Model data: solid symbols, 8t = 0.1; open symbols, lit = 0.2 (momentum conserved during cutoff). (b) Experimental data, from reference [18], at Red = 0.26.106.
approximately as (To/St) 3, with a proportionality constant of order 20/is o n the CDC 7600 at the University of London Computer Centre. Corresponding power spectra of v(t) are shown in Fig. I0 for the same three model runs; the low-frequency roll-off is faster than for u(t), but otherwise the spectra are similar. Figures 9 and 10 support the conclusion that the model flow statistics converge as the time increment is prngressively reduced. Figure 11 compares model v(t) spectra 2dia. downstream with the measurements of
214
A. V. J. EDWARDS and i
~u
i
C. L. MORFEY i
t
(o)
•
•
•
.
.
.
.
ii--
.
.
.
•
-~__ ;_ __-" . . . . . . . .
.
[]
I I ! I
L I
0 -0.2
-0.1
I
l
~
0.1
0
t
I
I
I
J
0.2
I
(b)
!
4
0
l
-0.2
I
-0.1
I
I
0
I
0.1
0.2
~1 : (Y-]) / x
Fig. 8. As Fig.7, but for kurtosisof axial velocity. Bradshaw et a/.[17] in a turbulent jet. The real jet spectra lack the strong "hump" between S = 0.2 and 1 which is characteristic of the model jet flow, and which presumably arises from the enforced circumferential coherence of the model. Finally the normalized covariance R,v(r) is presented in Fig. 12, for two points in the model flow; the corresponding turbulent-jet measurements of Lau et a/.[16] are superimposed for comparison. Qualitative agreement is evident, although descrepancies at zero delay (z = 0) point to differences in mean momentum transfer across the shear layer in the model and real jets. 5. TWO-POINT STATISTICS OF THE VELOCITY FIELD The "convected-pattern" nature of the model flow field, which is strikingly shown on cine-film and has been mentioned in Section 3, may be expressed statistically by means of two-point covariances. An example is shown in Fig. 13, where values of u at two different points seqarated by Ax have been used to form a normalized covariance R,, (Ax, z). The sequence of displaced peaks indicates a convected pattern, and the persistence of high correlations (R,, -~ 0.5) out to Ax = 2.8 shows that the pattern maintains its identity from the fixed point at 1.5 dia. downstream out to about 3 dia. Corresponding data for a real turbulent jet[19] are shown in the lower half of Fig. 13. Convection of the velocity field is still apparent (with a pattern convection speed of around 0.6,
A computer simulation of turbulent jet flow 101
.....
~
.......
i
~
"..
215
. . . . . . . .
100
10 -1
10 -2 i
i~l
%.
\,@,. "-....
10 -3
v ~
Fg~\,,,,,^"._ y , "~, ~-'.". '~/~,~ "~
10-4
..
-^10 -5
.....
q .1
~
,
, , , J,,I
,
|
1
I
~ v•
I
'
•
•
10
S
Fig. 9. Single-sided power spectral density of axial velocity at (x = 10, y = 1), with respect to Strouhal number S. Calculated from axisymmetric flow model with different time increments. Code: . . . . . . . 8t = 0.044 (T = 21); .... ,8t = 0.1 (T = 84); - - , 8t = 0.2 (T = 112).
I01 i
i
,
i
,
, I
i
i
i
i
,
i
i
,
I
i
i
,
i
i
t
i
100
10 -I i j
"°o°.
pi j
%%~. " " . . %
10-2
" ~ X "°% ... o.. .... /,~
"°.
10 -3
10 -4
.¢:. 10-5
i
t
i
i
I I
.1
,
i
I
i
i
i
i
i I
1
i
i
i
i
i
i
i
S
i
10
Fig. 10. As Fig. 9, but for radial velocity component.
as c o m p a r e d with 0.5 in the model flow); but the correlation in the real jet d e c a y s much more rapidly, falling to 0.5 by about 2 dia. d o w n s t r e a m from the nozzle exit (as o p p o s e d to 3 dia. in the model jet). 6. INTEGRAL BALANCES AND JET FLOW ENTRAINMENT An internal c h e c k on the c o m p u t e d flow field is provided by the integral equations for mass, m o m e n t u m and energy conservation. It is also of interest to c o m p a r e the flow entrainment rate in the model flow with that measured on real turbulent jets; the entrainment rate is an indirect measure of the rate at which the mean velocity profile spreads.
A. V. J. EDWARDSand C. L. MORFEY
216 100
'
,
,
'
, ' , ' 1
.
.
.
.
.
.
.
.
I
.
.
.
.
.
.
.
.
-
10-1
10-2 Y=l.2
~,
,
10-3 Y = 0.8
z
\X-~v.
-'..
10 -4 .1
.01
1
S
10
Fig. I1. Comparison of model and experimental radial-velocitypower spectra at x = 4, for two radial positions. Code: solid dymbols, model (St = 0.1); open symbols, experiment[17], at Red= 0.35. 106. 6.1
The fluxes of mass, momentum
and energy
The quantities FI = F3
y~ dy,
fo'°
F2 =
fo'°
y(p + u 2) dy,
-
= ffOy(-~ + z1~ + ~ 1 uv ) dy
(12a-c)
are measures of the mean fluxes of mass, x - m o m e n t u m and energy across any plane x = constant, out to radius YD. As such, their rates of change with x must be balanced by fluxes across y = Yo as follows aF1 = _ Yofi(Yo),
3x
OF3c~x-
y o {-pv+
aF2
ax
= -
youv(yo),
21-U2V+I - ~V
}(Yo).
(13a-c)
These equations will be applied at the following radii in the model flow: y o = 2 (13a); y0 = 1,~ (13b); y 0 ~ (13c). Note that the r.h.s, of (13b) and (13c) tend to zero in the limit Y0"-->~. For purposes of checking the momentum and energy balances, the pressure terms in F2 and F3 will be neglected; the pressure in the model flow has a standard deviation of order 0.1 and a mean close to zero, which makes its contribution small compared with the axial velocity terms. Viscous dissipation is absent in the model, so does not appear in (13c). For the mass and momentum balances at finite radius, we introduce the notation
F(x)=
fo ~
yti dy,
G(x)=
yu2 dy; fo'--
(14)
while for testing the invariance of the total energy and momentum (Yo-->~) we define
E(x)=fo y(~--~+
dr, M(x)= f f yu2dy.
(15)
Thus, we expect from mass conservation
F'(x) =
-20(2) =
[(x),
say,
(16)
217
A computer simulationof turbulentjet flow 1.00 ~
12(o)
S;
m •
m
g
e'm
0.40
0.20
g tr"
=
OlO0
•
•
-0.20
•
B
=
•
-040 -0.60
! 1
-0.08 --100 --~.O0
I --2.~0
J --2,00
L --I.50
' J -I.00-0.50
' 0.00050
•
m
i I00
i
i
i
1.500 2.00 2.50 3.00
LOg
Rw(4,.4)
I.OO _ 12(b)
0.80
m
0.600.40f
m
= "
0.20 [ => n..
" "
=
=
_ =
0.00 -0.20 -0.40
0.60 • •
-0,80
-100
= I i -3.00 -2.50 -2.00 -I I 50 -I.00
_OLSO
I I I 0.00 0.50 1.00 1.50
21.00
= 2.50 75.00
Log
Rud4, 21
Fig. 12. Normalized covariance R,o(r), for two positions: (a) x = 4, y = 0.4. (b) x = 4, y = 2. Model data: l , ~t = 0.1. Experimentaldata: - - , from reference [16], at Rex= 0.4.106. and (to the extent that pressure terms in FI, F2 are small) 1
G'(x) - - uv(l) = g(x), E'(x)- O, M'(x)- O.
(17a-c)
Computed values o f / ( x ) , F(x), g(x), G(x), M(x) and E(x) are presented in Table 1, for the steady-state model flow. As a rough check, trapezoidal integration of fix) and g(x) between x = 2 and 8 gives (0.130,-0.067); the corresponding changes in F and G are (0.127,-0.089). The overall momentum and energy flox estimates M and E are constant within about 2%. We conclude from these tests that the discrete-vortex calculation method provides reasonably accurate overall balances, except for the turbulent shear stress across y = 1 which is roughly 25% too low. This may in part be a statistical problem resulting from the limited averaging time avialable. 6.2 Jet flow
entrainment
The values of F(x) in Table 1 provide a measure of mean entrainment by the model jet flow, CAF Vol. 9, No. 2--H
218
A.V.J. EDWARDSand C. L. MORFEY 1.0
I
I
I
I
ax=.4
0.8
(a) 1.2
0.6 0.4 0.2 0 -0.2 -0.4
-0.6
0
I
I
I
I
I
1
2
3
4
5
Time delay AT
1.0
0.8
i
t
,
i
,
Separation
(b)
6X=.4
0.6 1.2 1.6
0.4 0.2 0 I
I
2
;
I
4
Time delay fiT
Fig. 13. Two-point normalized covariance R,,(Ax, r), with fixed point at (x = 3, y = 1) and variable downstream displacement Ax. (a) Model, calculated using 8t = 0.1. (b) Measurement in turbulent jet[19], at Red =0.13.106 . Table 1. Quantities involved in mass, momentum and energy balances. Calculated from steady-state model jet flow with 8t = 0.1 (cut-off does not conserve momentum) x
102/(x)
F(x)
102g(x) G ( x ) M(x)
E(x)
2.0 3.0 4.0 5.0 6.0 7.0 8.0
4.50 3.28 1.57 1.06 0.98 3.61 0.45
0.595 0.633 0.652 0.662 0.677 0.694 0.722
- 2.28 - 0.96 - 0.58 -0.40 - 1.49 - 1.53 - 1.14
0.487 0.491 0.485 0.486 0.489 0.483 0.486
0.476 0.462 0.452 0.445 0.430 0.405 0.387
0.514 0.524 0.519 0.518 0.521 0.521 0.526
which may be c o m p a r e d with laboratory measurements[2,20,21]. total mean mass flow rate across the plane x = const.
Q(x)/Q(O) - F(x)/F(O) = 2F(x),
Specifically, if Q(x) is the
(18)
since F ( 0 ) = 1/2. The result is approximate in that F(x) defines the jet flow rate only out to y = 2 . Thus eqn (18) holds for the first 4 or 5dia. of the jet (where the flow profile is contained almost entirely within y = 2), but becomes inaccurate further downstream. Figure 14 shows that the entrainment over the first 4 dia. is of similar magnitude in the model and real turbulent jets; the rate of entrainment in the model flow is initially more rapid, and then falls below the real jet value. Physically, the entrainment of flow by a jet is connected with the diffusion of the jet velocity profile (see Appendix); the observations above are, therefore, already anticipated by Fig. 5, which shows the model shear layer spreading more rapidly than the real one to begin with, then slowing down and being overtaken around x = 4.
A computer simulationof turbulentjet flow 2.0
,
,
,
,
219
,
1.8 Q Qo 1.6
1.4
1.2
1.0 0
I
I
I
I
i
2
4
6
8
10
12
X Distance from orifice
Fig. 14. Comparisonof flowentrainment in the initial regionof modeljet and real turbulentjets. Code:O, model (St = 0.1);--, Rea = 0.1.10*[2];---,
Red = (0.06-0.3). 106120]; - - , Rdd = 0.9.10s[21l.
7. CONCLUSIONS 1. Axisymmetric turbulence in a high Reynolds number starting jet has been simulated numerically with the aid of ring vortex elements. Sufficiently long runs were obtained (up to t ~ 140) to yield statistical properties of the steady-state flow for several diameters downstream of the orifice. Further detail, including results for sinusoidally-modulated jets, may be found in Ref. [22]. 2. A striking qualitative feature of the model jet flow is the formation of almost-regular clusters of vorticity; the average length scale of the clusters increases (by amalgamation) as they are swept downstream. In this respect the idealized axisymmetric jet resembles laboratory jets. 3. The mean spread rate of the initial mixing layer (and the corresponding entrainment over the first 4 dia.) are quite similar in the model and real jets, despite the absence of viscosity in the model simulation. In other respects, however, the steady-state model flow differs significantly from real turbulent jets. Turbulence intensities (axial and radial) in the shear layer are higher by a factor of 2, and two-point velocity correlations reveal a greater persistence of eddy structure with downstream separation than is measured in the laboratory. Similar findings, albeit based on rather short data samples (To < 45) have been reported by Acton[13]. Ashurst[23] has simulated a two-dimensional mixing layer with line vortex elements and also found (u,,,,, v,,,s) values of order 20-25% but was able to reduce these by allowing each line vortex to spread exponentially. 4. Further differences between the axisymmetric model flow and real turbulent jets are apparent in Fig. 4, which shows the u(t) time history near the inner edge of the initial mixing layer. The predominance of negative spikes in the laboratory data is not matched in the model, although the version with 8t = 0.2 (a time increment too large for accurate convergence at this position) is closer to the measured trace. We conclude that if the negative spikes measured by Lau et a/.[l'6] are indeed accounted for by an "axisymmetric vortex structure", as Lau and Fisher[24] have suggested, the dynamics of such a structure must differ significantly from the strictly axisymmetric inviscid model investigated here. 5. Possible future developments of the jet model include: (a) Addition of a co-flowing outer stream. (b) Addition of a plane boundary normal to the jet axis, to simulate an impinging jet. (c) Allowing the jet and ambient fluids to have different densities. (d) Allowing the jet circulation to decay downstream, e.g. by removing vortices whose radius falls below a limiting value. However, the present work provides a caution against transferring conclusions from axisymmetric turbulent models to real jet flows.
220
A.V.J. EDWARDSad C. L. MORFEY
Acknowledgements--The present study was initiated by Prof. P. O. A. L. Davies, whose guidance during the first author's Ph.D. studies is acknowledged. Financial support was provided by the Science Research Council and the University of Southampton. REFERENCES 1. G. L. Brown and A. Roshko, On density effects and large structure in turbulent mixing layers. J. Fluid Mech. 64, 775-816 (1974). 2. S. C. Crow and F. H. Champagne, Orderly structure in jet turbulence. J. Fluid Mech. 48, 547-591 (1971). 3. A. Roshko, Structure of turbulent shear flows: a new look. AIAA J. 14, 1349--1357 (1976). 4. G. K. Batchelor, Computation of the energy spectrum in homogeneous two-dimensional turbulence. Phys. Fluids 12, II.233-II.239 (1969). 5. C. L. Morley, Sound radiation from a simulated jet flow. In Mechanics of Sound Generation in Flows (Edited by E.-A. Mfiller), pp. 12-18. Springer-Verlag, Berlin (1979). 6. P. O. A. L. Davies and J. C. Hardin, Potential flow modelling of unsteady flow. In Numerical Methods in Fluid Dynamics (Edited by C. A. Brebbia and J. J. Connor). Pentech Press, London (1974). 7. A. J. Chorin and P. S. Bernard, Discretization of a vortex sheet with an example of roll-up. J. Comput. Physics 13, 423-429 (1973). 8. R. R. Clements, An inviscid model of two-dimensional vortex shedding. J. Fluid Mech. 57, 321-336 (1973). 9. R. R. Clements and D. J. Maull, The representation of sheets of vorticity by discrete vortices. Prog. Aerospace Sci. 16, 12%146 (1975). 10. R. R. Clements, Flow representation, including separated regions, using discrete vortices. In Cutup. Fluid Dynamics, AGARD LS-86 (1977). 11. D. J. Mautl, An Introduction to the Discrete Vortex Method. Cambridge University Engineering Daptartment. CUED/A-Aero/TR8 (1979). 12. H. Lamb, Hydrodynamics, 6th Edn, pp. 202-249. Dover, New York (1945). 13. E. Acton, A Modelling of Large Jet Eddies. Ph.D. Thesis, Cambridge University (1977); see also J. Fluid Mech. 98, 1-31
(1980). 14. C. L. Morley and A. V. J. Edwards, Energy and Momentum Constraints on the use of a Cutoff in Lumped-vortex Flow Models. University of Southampton, ISVR Memorandum 587 (1978). 15. S. M. Datums and D. Kiichemann, On a vortex-sheet model for the mixing between two parallel streams. Proc. R. Soc. Lurid. A 339, 451--461 (1974). 16. J. C. Lau, M. J. Fisher and H. V. Fuchs, The intrinsic structure of turbulent jets. J. Sound Vib. 22, 379-.406 (1972). 17. P. Bradshaw, D. H. Ferriss and R. F. Johnson, Turbulence in the noise-producing region of a circular jet. J. Fluid Mech. 19, 591-624 (1964). 18. B. Bose, Experimental Study o[ the Wave-like Characteristics of Turbulence in Shear Layers and Behind Grids. Ph.D. Thesis, Southampton University (1968). 19. M. J. Fisher and P. O. A. L. Davies, Correlation measurements in a non-frozen pattern of turbulence. J. Fluid Mech. 18, 97-116 (1963). 20. B. J. Hill, Measurement of local entrainment rate in the initial region of axisymmetric turbulent air jets. J. Fluid Mech. 51,773-779 (1972). 21. L Maestrello and E. McDaid, Acoustic characteristics of a high-subsonic jet. AIAA J. 9, 1058-1066 (1971). 22. A. V. J. Edwards, The Computer Simulation of an Unsteady Flow. Ph.D. Thesis, Southampton University (1978). 23. W.T.Ashurst, Numerical simulation of turbulent mixing layers via vortex dynamics. In Turbulent Shear Flows--l (Edited by F. Durst, B. F. Launder, F. W. Schmidt and J. H. Whitelaw), pp. 402-413. Springer-Verlag, Berlin (1979). 24. J. C. Lau and M. J. Fisher, The vortex-street structure of "turbulent" jets. J. Fluid Mech. 6'/, 299-337 (1975). APPENDIX. R E L A T I O N B E T W E E N E N T R A I N M E N T AND D I F F U S I O N OF MEAN VELOCITY IN A JET For a monotonic mean velocity profile ti(y), a round jet (Fig. 15a) has its mass flow distributed with respect to li according to . a~ - i
P(~)=yu ~yy
/fo
®
(~
(o)
y
(AI)
P(;)
//~
Increasing x
II
IN.'~I
. Omox
0
mx
Fig. 15. (a) Mean axial velocity profile in jet. (b) Amplitude density of mean velocity with respect to mass flow rate.
A computer simulation of turbulent jet flow
221
In other words, P(a) da is the proportion of the mean flow rate (across a plane x = const.) for which the mean velocity lies between (a, a + da), thus
oP(a) da = I. Let Q(x) and X(x) denote the total mean flux of mass and x-momentum across
(A2) x = const., then
dX = a dQ (pressure terms and u 2 - (a) 2 neglected) = aQ. P0i) da.
(A3)
Integration gives
f~
X(x)/Q(x) = J o aP(a) dfi.
(A4)
Momentum conservation implies X(x)= const.; diffusion of the mean velocity profile (Fig. 15b) implies a decrease in the integral on the right; thus, if the profile spreads, Q(x) must increase with x, meaning that fluid is entrained by the jet.
0045-793018110601-02051502.0010
PergwonPressLid.,1981.PrintedinGreatBritain
A COMPUTER
SIMULATION OF TURBULENT A. V. J. EDWARDSand C. L. MORFEV
JET FLOW
Institute of Sound and Vibration Research, University of Southampton, Southampton, England
(Received 14 December 1979) Al~traet--Tbe evolution of a round turbulent jet is simulated numerically under the constraint of axial symmetry. The vortex sheet shed from the orifice is represented by vortex ring elements, with a velocity field cut-off to control close encounters. Large-scale vortex clusters form in the model jet, similar to those inferred from laboratory flow visualization experiments. However, comparisons of statistical properties reveal significant differences between the axisymmetric model flow and real turbulent jets. 1. I N T R O D U C T I O N
Flow visualization experiments on plane turbulent mixing layers[l] and round turbulent jets [2] have revealed the occurrence of fairly regular large-scale features, commonly described as orderly structure. The typical length scale (wavelength) is around 0.3x, where x is the downstream distance from the start of the shear layer. Such features have been observed over a range of Reynolds numbers Re~ of order 3 x 104-3 x 10613]. We demonstrate in this paper that a similar orderly structure emerges from a highly idealized model jet flow in which two of the three vorticity components are suppressed. The model flow consists of an axisymmetric starting jet, which issues from a circular opening of radius R0 in a rigid plane boundary (Fig. 1). The normal velocity over the opening jumps, at t = 0, from zero to a constant value U0, and an axisymmetric vortex sheet is shed from the lip of the opening. The jet and ambient fluids are assumed inviscid and incompressible, with identical density pp. Instability of the vortex sheet leads to a chaotic, quasi-turbulent (but still axisymmetric) flow field, with some qualitative resemblance to a real turbulent jet. The evolution of the large-scale orderly structure in the model is described in Section 3, and statistical properties of the model jet (as it approaches a steady state) are discussed in Sections 4 and 5. We begin in Section 2 by describing the calculation method in some detail; for convenience all quantities will be normalized in what follows to unit values of p0, R0 and Up. The objective of the study is not to demonstrate that certain statistical features of turbulence in jets can be quantitatively reproduced in an axisymmetric model: such agreement is not to be expected. Rather, the aim is to set up a dynamically self-consistent axisymmetric model, to observe its statistical properties and to contrast these with real jet measurements. The axisymmetric jet model, thus, becomes a tool for studying the large-scale consequences of complicated vortex motions.t Numerical flow experiments of a similar type--but directed at the small-scale aspects of turbulent motion--have previously been described by Batchelor[4].
\ Source distribution to
/ /
:
~'~.~ ~ , ~
maintain boundary condition over opening
d••m
~ n R i n g vortices ring vortex sheet
--~ Image vortices to maintain boundary condition over solid surface ~
f
Uniform flow of velocity U
!/ / /
Fig. I. Jet flow simulation using sources and vortices. tFurther results dealing with the local pressure and radiated sound (including comparisons with experiment) will be reported in a later paper. The far-field sound predictions are summarized in Ref. [5]. 205
206
A.V.J. EDWARDSand C. L. MORFEY 2. THE L U M P E D - V O R T E X FLOW R E P R E S E N T A T I O N
We follow Davies et a/.[6], Chorin and Bernard [7] and Clements[8] and represent the continuous vortex sheet by an array of discrete elemental ring vortices.'t Each vortex element moves under the influence of all the others, including the set of image vortices reflected in the orifice plane (Fig. I). The prescribed normal velocity U0 over the opening is simulated by a source disk of appropriate strength. An unwanted consequence of discretizing the vortex sheet is the appearance of unrealistic induced velocities close to each vortex element. In particular, this causes problems when two vortices come close together. The need for some kind of smoothing of the discrete-vortex velocity field was pointed out by Chorin and Bernard [7]; they proposed a cut-off, or upper limit, for the velocity induced by each vortex, and this has been adopted in the present study. The precise manner in which the cut-off is applied turns out to be important; we return to this point later. 2.1 The shedding process In general, if the fluid crosses the orifice plane with a uniform normal velocity U(t), the assumption that a vortex sheet leaves the lip parallel to the axis implies a "rate of circulation shedding Ou F= - f o u-z-dyoy
:-21 U2(t)"
(1)
Here cylindrical coordinates (x, y, ~b) are being used, and the circulation is in the ~b direction. It follows from (1) that if the circulation shed in the interval (t, t + 6t) is lumped into a ring vortex, the circulation around the ring is F = ~U2(t)6t + O(6t) 2. Likewise the centroid of the circulation shed in this interval is located at (x, y) = (U(t)St/3, 1) at the end of the time interval, apart from O(6t) 2 terms.:~ Thus, the strength and location of each vortex element as it is first shed are determined by the efflux velocity U(t). 2.2 The motion of a vortex element An individual vortex element is subjected in general to three superimposed velocity fields. These are the velocity induced by all the other vortex elements (including images); the velocity induced by the vortex itself; and the velocity field of the source disk. (a) The velocity field of other vortices Associated with each elemental vortex ring there is a velocity field, whose axial and radial components at a point (x, y) are given by u = FR(RJ - yK),
(2)
v = FRxK.
(3)
Here F is the circulation of the elemental vortex, R is the radius of the elemental vortex, and (x, y) are axial and radial coordinates with the centre of the vortex, ring as origin. The quantities J and K are defined by
(4)
1 E(k) J = lrr23 d - F)' -1 1 [ K --- 7rr23~ 2K(k)
(2-ke)E(k)} (1 k 2)
"
(5)
Here K(k), E(k) are complete elliptic integrals of the first and second kind respectively, k 2 = 1 - rm2/r22,and (rm, rE) are the distances of the nearest and farthest points on the vortex ring from the point (x, y). Equations (4) and (5) are equivalent to expressions given by Lamb[12]. tFor a review of the discrete vortex calculation method, see Refs. [9-11]. ~The vortex sheet is assumed to leave the orifice axially. A detailed derivation appears in Ref. [22].
A computer simulation of turbulent jet flow
207
(b) The self-induced velocity of a vortex ring A ring of zero cross section has an infinite self-induced velocity, and a better approximation is clearly called for. We recognize that initially each ring represents a tubular segment of vortex sheet. Calculating the self-induced velocity of the centroid of such a segment gives F /1 +
16',
(6t
(6)
This expression strictly applies only so long as the sheet remains cylindrical. It is best regarded as a correction term which tends to zero as 6t ~ 0, but which appears worth including in the interest of greater accuracy close to the orifice, where the sheet begins to break up. (c) The velocity lield of the source disk Two methods are available for meeting the exit-plane boundary condition on normal velocity. One is to use image vortices (as in Fig. 1) together with a source disk of uniform strength 2U(t) per unit area. The alternative is to dispense with image vortices (thereby saving considerable computing time) and to use a non-uniform source disk of strength 2 U ( t ) 2uv(y, t), where uv is the normal velocity induced on the exit plane by the real vortices in x > 0. In practice a combination of these methods was used. Up to t = 15, with relatively few vortex elements present, the image method was used for maximum accuracy. Beyond t = 15, an approximate version of the second method was used: rather than updating the source distribution at each step, use was made of the fact that uv varies very little with time after about t = 10. Thus, an average value, taken over the interval 10 < t < 15, was inserted in the source strength expression. The source in the orifice plane is represented by a series of concentric rings of appropriate strengths. Obviously, from the point of view of the economical use of computer time the smaller the number of rings used to represent the continuous source distribution the better. However, using too few rings gives rise to errors near the orifice, in particular at the initial shedding location. A small study led to the conclusion that 30 rings inside the orifice and ten more outside constituted a reasonable representation of the required source distribution. It is necessary to include a number of more widely spaced rings (sinks rather than sources) outside the orifice, in order to maintain the boundary condition u = 0 for y > 1 when the vortices are limited to the x > 0 half space. 2.3 The need for a velocity cut-off Discrete vortex representations of flows with continuously distributed vorticity produce unrealistic velocity fields close to the vortex elements. In particular, this causes problems when two vortices come close together. To ensure a smooth roll-up to the starting vortex, without the strong mutual repulsions otherwise associated with close encounters, we follow Chorin and Bernard[7] in applying a cut-off or upper limit to the induced velocity associated with each vortex. In the early stages of the present work, and in a similar study of a starting jet by Acton[13], the importance of conserving momentum during the cut-off process was overlooked. As a result, spurious vortex accelerations were introduced whose effect on the calculated pressure field was particularly serious. Some of the results to be presented:~ were obtained with the momentum-conserving cut-off factor derived in[14], namely
f(rO = rl2/rc 2,
(rl < rc = cut-off radius);
f(rD = 1,
(rl <- rc).
(7)
The above weighting factor is applied to both the radial and axial velocity components induced tNote that t represents the nondimensional time (normalized by RolUo). Since F - 6t as 8t tends to zero, the log term is typically a small correction. ~Speciflcally, those with 8t = 0.2.
208
A.V.J. EDWARDSand C. L. MORFEY
by either member of the vortex pair on the other. However, many of the earlier results are included for comparison, particularly where momentum conservation did not prove critical. The choice of r~ is a compromise between the conflicting requirements of a smooth initial vortex roll-up on one hand, and the need to minimize the invocation of cut-off on the other, as discussed in Ref. [14]. The results presented are based on r,. = 0.05.
2.4 The integration scheme (a) Linear multistep procedure (Adams-Bash[orth ) Suppose that the velocity u(t) is fitted in the vicinity of t = 0 by the order k polynomial k+l
~(t) = ~ aitH;
(8)
i=1
then to determine the (k + I) coefficients ai we must equate ti(t) to u(t) at (k + 1) points (times). These may be taken as t = 0 , - 1 . . . . . - k where the time scale is normalized on the interval between successive points. The resulting (k + 1) equations take the form k+¿
Zai(1-j)i-I=ui,
]=1,2 ..... k+l
(9)
i=l
where uj represents the value of u taken ( j - 1 ) intervals back from t = 0. We use the term "back-velocity" to refer to the velocity of a vortex at some previous time-step. Finally, (xl - x0) = Ax is estimated from eqn (8) as 1
Ax-'
k+l
y0
k+l
(10)
t~(t)dt = ~•: a_~_ i - ~ b~u~,say.
The coefficeints bi are listed below for orders up to k = 5.
k
Dbl
Db2
0 1 2 3 4 5
i 3 23 55 1901 4277
-1 - 16 -59 -2774 -7923
Db3
5 37 2616 9982
Db4
-9 - 1274 -7298
Db5
251 2877
Db6
-475
D 1 2 12 24 720 1440
To obtain bi, the values in each row must be divided by the corresponding number D in the final column.
(b) Estimation of the correction order required In the computation of Ax from eqn (10) one may not invoke the use of an arbitrarily high correction order without suitable justification. In the present case we have, in fact, a highdensity classical many-body problem so that Ax is of the order of the separation of the elemental vortices. One might suspect that use of more than a very few back-velocities of a particular vortex in the calculation of its subsequent position would be erroneous, since one would be imposing continuity upon velocity derivatives calculated from greatly differing past configurations. In other words, velocity derivatives of too high an order become spuriously large causing eqn (10) to diverge. It was found that in the present case the use of more than one back-velocity gave rise to spurious behaviour of the elemental vortices. Accordingly k was taken as I, and vortex
209
A computer simulation of turbulent jet flow
coordinates at each step were calculated from 1 x.+,(i) = x~(i)+ ~(3u~(i)- u._,(i)) fit. 1
y.+,(i) = y.(i)+ ~(3v.(i)- v~_,(i)) 8t
(11)
where u.(i) and v.(i) are, respectively, the axial and radial velocities of the ith vortex at the nth time increment. 3. V I S U A L
CHARACTERISTICS
OF THE
SIMULATED
FLOW
In this section we briefly describe the model starting jet in terms of the evolving vorticity distribution as viewed in a plane containing the jet axis. The development of the simulated flow (from time t = 0) may be divided into three fairly distinct phases: (i) The formation of a starting vortex. (ii) The detachment of the starting vortex. (iii) The steady-state situation. Figure 2 shows the state of the simulated flow at four different times, calculated with 8t = 0.1. The same calculation using 8t = 0.2 yields similar vortex patterns (Fig. 3) despite differences of detail; it appears that large scale features of the vorticity distribution are not sensitive to St. The first frame (t = 6.2) shows the starting vortex about to detach, leaving a trail of smaller vortex clusters which show up clearly in the second frame. Thereafter, small clusters form at fairly regular intervals near the orifice, and grow by amalgamation with neighbouring clusters as they are swept downstream. Eventually the flow in any given region downstream achieves a statistically steady state, in which the general appearance of the vortex pattern no longer changes with time. This occurs when the starting vortex has passed well beyond the region concerned. The formation in real turbulent jets of similar vorticity concentrations, and their growth in scale (both size and separation) with increasing distance downstream, have been described by Brown and Roshko[l], Damms and Kficbemann[15] et al. The fact that similar features appear TIME
6.2000
61 VORTICES
TIME 16.2000
~.'~ •
~ .,~ .~. r,:.
,
TIME 26.2000
-,, ,.¢.
.~,~.
,
..,~-..
,~.'.~
,
246 VORTICES
..
154 VORTICES
~:~: • '.:a'~-
,
TIME 36.2000
..
,
358 VORTICES
" . "~-~'-..
." ::A~',..::.,,', ".~..
Fig. 2. Evolution of vortex clusters in axisymmetric jet, calculated using 8t = 0.1. The lower half of each frame contains the instantaneous velocity profile of the radial component, v, along the lip-line (y = 1.0).
210
A. V. J. EDWARDSand C. L. MORFEY TIME 6.2000
31 VORTICES
81 VORTICES
TIME 16.2000
....... :.f!:i::.
vv
TnlE 26.2000
151 VORTICES
:~ #, ,.:,,.. ~.....:,:~..:"
...,.
A A
, ~
•
: " ;.-...;.
~:~ :" ::
~
:
W
TIME 36.2000
. . . . . . . . . . . .
,
,
;" :~ "' ....:W:'::'.
A, A ~
181 VORTICES
:..?:..,.
•
.:
" '":: ..:'":...
,X~A\/;A
:...
...i.:..i:q
...~.:::
x~i//'-~x
Fig. 3. As Fig. 2, but calculated using ~St= 0.2 (momentum conserved during cut-off).
in the model jet is the main justification for studying the model flow in more detail; we can exploit the accessibility of information from the numerical model to provide insights into the large-scale aspects of the real three-dimensional flow. 4. S I N G L E - P O I N T S T A T I S T I C S OF THE V E L O C I T Y F I E L D
In this section we examine the statistics of the time-dependent velocity at various points in the model flow. For this purpose the initial development of the flow from rest is ignored, and statistical estimates of mean velocity, variance, spectra, etc. are based on the steady-state part of the time-history. Comparisons are made between model flows with different time increments St, and also between the axisymmetric model and real laboratory jets. 4.1 The amplitude density of u(t) and its moments Figure 4 provides a visual comparison of u(t) time histories for two model realizations (St =0.1,0.2) and a turbulent jet in the laboratory[16]. The position in the flow is (x = 4, y =
0.8), which lies in the initial turbulent shear layer of the real jet. The two model time-histories, though qualitatively similar, exhibit detailed differences especially in the low-velocity regions: here trace (a) is unrealistically rounded, while trace (b) shows more of the sharp downward spikes present in the real flow (c). The difference between traces (a) and (b) is attributed partly to the fact that momentum was not conserved during the cut-off process in the earlier calculation (a)--see Section 2.3. Perhaps of greater importance, however, is the fact that convergence of the flow statistics at this location downstream (x = 4) appears to require a time increment of 8t = 0.1 or less: thus trace (b)--with 8t = 0.2--retains some influence of step size and does not truly approximate the limit 8t ~ 0. Further downstream, beyond x = 6, this problem disappears (as noted below). Mean velocity profiles ~(y) are shown in Fig. 5. Figure 5(a) tests the model flow profiles for similarity in the first 5 dia. downstream; it can be seen that similarity is achieved only beyond about 2 dia. in this case (St = 0.1). With smaller and larger time increments (St = 0.044, 0.2) the distance downstream at which similarity is achieved is roughly 1 and 3 dia., respectively. These results suggest that the finite step size (with its associated cut-off) provides an
A computer simulation of turbulent jet flow I
!
I
211 i
(a)
(b)
8
0
I0
20
30
40
T
Fig. 4. Axial velocity time histories at (x =4, y =0.8). (a) Model, calculated using 8t = 0.1. (b) Model, calculated using 8t = 0.2 (momentum conserved during cut-off). (c) Measurement in turbulent jet [16], at Rex = 0.4.10 6.
artificial viscosity, and leads to a "transition region" close to the orifice similar to that observed on real jets at finite Reynolds numbers (Uodh, ~ 104- I0~). In Fig. 5(b) a typical model profile (at x = 4) is compared with the profile measured at the same location in a turbulent jet[17]. The agreement is surprisingly good considering the axisymmetric nature of the model; this point is taken up in Section 6, where entrainment is discussed. Figures 6-8 show profiles of the standard deviation (u,~s), skewness (/~,) and kurtosis (K,) of the axial velocity u(t), at positions up to 5 dia. downstream of the orifice. In this region, real jet data collapse on the similarity coordinate r/= (y - 1)/x; accordingly 7/is used as the abscissa in all these plots. Each of the three model flow realizations (St =0.044,0.1,0.2) gives peak values of u,,,,s (axial turbulence intensity) around 25-30%,t i.e. roughly double what is measured in real jets[17]. This difference may be attributed to the greater orderliness imposed by the axisymmetry constraint in the model; corresponding differences are found in the local fluctuating pressure and the far-field radiated sound[5]. The skewness profiles in Fig. 7 show some significant differences between the 8t = 0.1 and 8t = 0.2 model flows, the latter corresponding more closely with real jets[18]. The differences appear to be due in part to the different cut-off procedures mentioned above, since they persist out to x = 10 where the effect of the time increment is no longer evident in other statistics (e.g. Fig. 6). tExcept where convergence has not been achieved, e.g. at x = 4 with 8t = 0.2.
212
A. V. J. EDWARDSand C. L. MORFEY
1.20
I
51a)
i.oo o.oo
°°.° F a~
-020 L
i
°°~.eo
~
-o.4o
-o.zo
o.oo (Y-I)/X DT-O.i
o2o
i
0.4o
o.eo
I.OO •
a m
cteo 0
[]
5(b)
•
0.6O
I=
=Z 0.40
0.20 0.00
m
i
A
-0.20 0.00 0120 o.~ o~ o'Jo ,ioo ,.2o' ,.~ i.~
•
,.~' 2'oo z.m
Y/RO
X/ROt4.0 Fig. 5. Mean axial velocity profiles. (a) Similarity plot of model results, calculated using 8t =0.1. Values of x are 2, 4, 6, 8, 10 with smallest x corresponding to largest range of abscissa. (b) Comparison of model and experimental velocity profiles at x =4. Code: [~, model (St =0.1); --, experiment[17], at Re~ =0.35' 106 (d = nozzle dia.).
Similar remarks apply to the kurtosis profiles of Fig. 8, where again the model results-despite the similarity collapse--bear only a modest resemblance to real jet data[18]. 4.2 Power spectra and covariances o/ u(t) and v(t) Figure 9 shows the spectral distribution of the u(t) variance with respect to Strouhai number (S = [dlUo), for three different model flow realizations sampled at (x = 10, y = I). If the three signals were statistically identical and the sample lengths had been equal, 90% of the spectral estimates would be expected to lie within a band 6 dB wide; this is a measure of the statistical accuracy of the results, and follows from the. use of 7 overlapping data segments for analysis. However, the usable sample lengths for the three runs (St = 0.044, 0.1,0.2) were different (T = 21,84, 112, respectively); therefore the resolution of the 8t = 0.044 spectrum is relatively coarse, which could account for the apparent bias of this spectrum towards higher Strouhal numbers. (The spacing AS of the spectral estimates is equal to 8IT.) The different sample lengths arose from a reduction in total flow duration To (for economy reasons) as 8t was decreased. The computing time for each run of the flow model varied
A computer simulation of turbulent jet flow
213
u
rmi
a
0.3
0.2 ÷
"4"11
4au
0.1
+ + ÷
++
++j
I
I
-0.4
-0.2
0
I
0.2
0.4
"I] = ( Y - l ) / X
Fig. 6. Profiles of standard deviation u,~s (axial turbulence intensity), plotted against similarity coordinate 7. Model data (x= 10): O, 8t =0.0M; II, 8t =0.1; A, 8t =0.2 (momentum conserved during cut-off). Experimental data (x = 8): +, from reference [17], at Rex = 1.4.106. 2
i
Bu
I
I
I
I
t
(°) I IDA
0
II
- I r-
-2
t -0.2
J3u
2
I
I
-0.1
I
0
0,1
I
(b)
I
1
I
0,2
I
0° in •
0
6 ,L
•
i
-1
-2
i -0.2
I
-0.I
I
I
0 TI = ( Y - l ) / X
I
0.I
I
0.2
Fig. 7. Skewness profiles for axial velocity, plotted against similarity coordinate 7/. Code: triangles, x = 4; circles, x = 6; squares, x = 10. (a) Model data: solid symbols, 8t = 0.1; open symbols, lit = 0.2 (momentum conserved during cutoff). (b) Experimental data, from reference [18], at Red = 0.26.106.
approximately as (To/St) 3, with a proportionality constant of order 20/is o n the CDC 7600 at the University of London Computer Centre. Corresponding power spectra of v(t) are shown in Fig. I0 for the same three model runs; the low-frequency roll-off is faster than for u(t), but otherwise the spectra are similar. Figures 9 and 10 support the conclusion that the model flow statistics converge as the time increment is prngressively reduced. Figure 11 compares model v(t) spectra 2dia. downstream with the measurements of
214
A. V. J. EDWARDS and i
~u
i
C. L. MORFEY i
t
(o)
•
•
•
.
.
.
.
ii--
.
.
.
•
-~__ ;_ __-" . . . . . . . .
.
[]
I I ! I
L I
0 -0.2
-0.1
I
l
~
0.1
0
t
I
I
I
J
0.2
I
(b)
!
4
0
l
-0.2
I
-0.1
I
I
0
I
0.1
0.2
~1 : (Y-]) / x
Fig. 8. As Fig.7, but for kurtosisof axial velocity. Bradshaw et a/.[17] in a turbulent jet. The real jet spectra lack the strong "hump" between S = 0.2 and 1 which is characteristic of the model jet flow, and which presumably arises from the enforced circumferential coherence of the model. Finally the normalized covariance R,v(r) is presented in Fig. 12, for two points in the model flow; the corresponding turbulent-jet measurements of Lau et a/.[16] are superimposed for comparison. Qualitative agreement is evident, although descrepancies at zero delay (z = 0) point to differences in mean momentum transfer across the shear layer in the model and real jets. 5. TWO-POINT STATISTICS OF THE VELOCITY FIELD The "convected-pattern" nature of the model flow field, which is strikingly shown on cine-film and has been mentioned in Section 3, may be expressed statistically by means of two-point covariances. An example is shown in Fig. 13, where values of u at two different points seqarated by Ax have been used to form a normalized covariance R,, (Ax, z). The sequence of displaced peaks indicates a convected pattern, and the persistence of high correlations (R,, -~ 0.5) out to Ax = 2.8 shows that the pattern maintains its identity from the fixed point at 1.5 dia. downstream out to about 3 dia. Corresponding data for a real turbulent jet[19] are shown in the lower half of Fig. 13. Convection of the velocity field is still apparent (with a pattern convection speed of around 0.6,
A computer simulation of turbulent jet flow 101
.....
~
.......
i
~
"..
215
. . . . . . . .
100
10 -1
10 -2 i
i~l
%.
\,@,. "-....
10 -3
v ~
Fg~\,,,,,^"._ y , "~, ~-'.". '~/~,~ "~
10-4
..
-^10 -5
.....
q .1
~
,
, , , J,,I
,
|
1
I
~ v•
I
'
•
•
10
S
Fig. 9. Single-sided power spectral density of axial velocity at (x = 10, y = 1), with respect to Strouhal number S. Calculated from axisymmetric flow model with different time increments. Code: . . . . . . . 8t = 0.044 (T = 21); .... ,8t = 0.1 (T = 84); - - , 8t = 0.2 (T = 112).
I01 i
i
,
i
,
, I
i
i
i
i
,
i
i
,
I
i
i
,
i
i
t
i
100
10 -I i j
"°o°.
pi j
%%~. " " . . %
10-2
" ~ X "°% ... o.. .... /,~
"°.
10 -3
10 -4
.¢:. 10-5
i
t
i
i
I I
.1
,
i
I
i
i
i
i
i I
1
i
i
i
i
i
i
i
S
i
10
Fig. 10. As Fig. 9, but for radial velocity component.
as c o m p a r e d with 0.5 in the model flow); but the correlation in the real jet d e c a y s much more rapidly, falling to 0.5 by about 2 dia. d o w n s t r e a m from the nozzle exit (as o p p o s e d to 3 dia. in the model jet). 6. INTEGRAL BALANCES AND JET FLOW ENTRAINMENT An internal c h e c k on the c o m p u t e d flow field is provided by the integral equations for mass, m o m e n t u m and energy conservation. It is also of interest to c o m p a r e the flow entrainment rate in the model flow with that measured on real turbulent jets; the entrainment rate is an indirect measure of the rate at which the mean velocity profile spreads.
A. V. J. EDWARDSand C. L. MORFEY
216 100
'
,
,
'
, ' , ' 1
.
.
.
.
.
.
.
.
I
.
.
.
.
.
.
.
.
-
10-1
10-2 Y=l.2
~,
,
10-3 Y = 0.8
z
\X-~v.
-'..
10 -4 .1
.01
1
S
10
Fig. I1. Comparison of model and experimental radial-velocitypower spectra at x = 4, for two radial positions. Code: solid dymbols, model (St = 0.1); open symbols, experiment[17], at Red= 0.35. 106. 6.1
The fluxes of mass, momentum
and energy
The quantities FI = F3
y~ dy,
fo'°
F2 =
fo'°
y(p + u 2) dy,
-
= ffOy(-~ + z1~ + ~ 1 uv ) dy
(12a-c)
are measures of the mean fluxes of mass, x - m o m e n t u m and energy across any plane x = constant, out to radius YD. As such, their rates of change with x must be balanced by fluxes across y = Yo as follows aF1 = _ Yofi(Yo),
3x
OF3c~x-
y o {-pv+
aF2
ax
= -
youv(yo),
21-U2V+I - ~V
}(Yo).
(13a-c)
These equations will be applied at the following radii in the model flow: y o = 2 (13a); y0 = 1,~ (13b); y 0 ~ (13c). Note that the r.h.s, of (13b) and (13c) tend to zero in the limit Y0"-->~. For purposes of checking the momentum and energy balances, the pressure terms in F2 and F3 will be neglected; the pressure in the model flow has a standard deviation of order 0.1 and a mean close to zero, which makes its contribution small compared with the axial velocity terms. Viscous dissipation is absent in the model, so does not appear in (13c). For the mass and momentum balances at finite radius, we introduce the notation
F(x)=
fo ~
yti dy,
G(x)=
yu2 dy; fo'--
(14)
while for testing the invariance of the total energy and momentum (Yo-->~) we define
E(x)=fo y(~--~+
dr, M(x)= f f yu2dy.
(15)
Thus, we expect from mass conservation
F'(x) =
-20(2) =
[(x),
say,
(16)
217
A computer simulationof turbulentjet flow 1.00 ~
12(o)
S;
m •
m
g
e'm
0.40
0.20
g tr"
=
OlO0
•
•
-0.20
•
B
=
•
-040 -0.60
! 1
-0.08 --100 --~.O0
I --2.~0
J --2,00
L --I.50
' J -I.00-0.50
' 0.00050
•
m
i I00
i
i
i
1.500 2.00 2.50 3.00
LOg
Rw(4,.4)
I.OO _ 12(b)
0.80
m
0.600.40f
m
= "
0.20 [ => n..
" "
=
=
_ =
0.00 -0.20 -0.40
0.60 • •
-0,80
-100
= I i -3.00 -2.50 -2.00 -I I 50 -I.00
_OLSO
I I I 0.00 0.50 1.00 1.50
21.00
= 2.50 75.00
Log
Rud4, 21
Fig. 12. Normalized covariance R,o(r), for two positions: (a) x = 4, y = 0.4. (b) x = 4, y = 2. Model data: l , ~t = 0.1. Experimentaldata: - - , from reference [16], at Rex= 0.4.106. and (to the extent that pressure terms in FI, F2 are small) 1
G'(x) - - uv(l) = g(x), E'(x)- O, M'(x)- O.
(17a-c)
Computed values o f / ( x ) , F(x), g(x), G(x), M(x) and E(x) are presented in Table 1, for the steady-state model flow. As a rough check, trapezoidal integration of fix) and g(x) between x = 2 and 8 gives (0.130,-0.067); the corresponding changes in F and G are (0.127,-0.089). The overall momentum and energy flox estimates M and E are constant within about 2%. We conclude from these tests that the discrete-vortex calculation method provides reasonably accurate overall balances, except for the turbulent shear stress across y = 1 which is roughly 25% too low. This may in part be a statistical problem resulting from the limited averaging time avialable. 6.2 Jet flow
entrainment
The values of F(x) in Table 1 provide a measure of mean entrainment by the model jet flow, CAF Vol. 9, No. 2--H
218
A.V.J. EDWARDSand C. L. MORFEY 1.0
I
I
I
I
ax=.4
0.8
(a) 1.2
0.6 0.4 0.2 0 -0.2 -0.4
-0.6
0
I
I
I
I
I
1
2
3
4
5
Time delay AT
1.0
0.8
i
t
,
i
,
Separation
(b)
6X=.4
0.6 1.2 1.6
0.4 0.2 0 I
I
2
;
I
4
Time delay fiT
Fig. 13. Two-point normalized covariance R,,(Ax, r), with fixed point at (x = 3, y = 1) and variable downstream displacement Ax. (a) Model, calculated using 8t = 0.1. (b) Measurement in turbulent jet[19], at Red =0.13.106 . Table 1. Quantities involved in mass, momentum and energy balances. Calculated from steady-state model jet flow with 8t = 0.1 (cut-off does not conserve momentum) x
102/(x)
F(x)
102g(x) G ( x ) M(x)
E(x)
2.0 3.0 4.0 5.0 6.0 7.0 8.0
4.50 3.28 1.57 1.06 0.98 3.61 0.45
0.595 0.633 0.652 0.662 0.677 0.694 0.722
- 2.28 - 0.96 - 0.58 -0.40 - 1.49 - 1.53 - 1.14
0.487 0.491 0.485 0.486 0.489 0.483 0.486
0.476 0.462 0.452 0.445 0.430 0.405 0.387
0.514 0.524 0.519 0.518 0.521 0.521 0.526
which may be c o m p a r e d with laboratory measurements[2,20,21]. total mean mass flow rate across the plane x = const.
Q(x)/Q(O) - F(x)/F(O) = 2F(x),
Specifically, if Q(x) is the
(18)
since F ( 0 ) = 1/2. The result is approximate in that F(x) defines the jet flow rate only out to y = 2 . Thus eqn (18) holds for the first 4 or 5dia. of the jet (where the flow profile is contained almost entirely within y = 2), but becomes inaccurate further downstream. Figure 14 shows that the entrainment over the first 4 dia. is of similar magnitude in the model and real turbulent jets; the rate of entrainment in the model flow is initially more rapid, and then falls below the real jet value. Physically, the entrainment of flow by a jet is connected with the diffusion of the jet velocity profile (see Appendix); the observations above are, therefore, already anticipated by Fig. 5, which shows the model shear layer spreading more rapidly than the real one to begin with, then slowing down and being overtaken around x = 4.
A computer simulationof turbulentjet flow 2.0
,
,
,
,
219
,
1.8 Q Qo 1.6
1.4
1.2
1.0 0
I
I
I
I
i
2
4
6
8
10
12
X Distance from orifice
Fig. 14. Comparisonof flowentrainment in the initial regionof modeljet and real turbulentjets. Code:O, model (St = 0.1);--, Rea = 0.1.10*[2];---,
Red = (0.06-0.3). 106120]; - - , Rdd = 0.9.10s[21l.
7. CONCLUSIONS 1. Axisymmetric turbulence in a high Reynolds number starting jet has been simulated numerically with the aid of ring vortex elements. Sufficiently long runs were obtained (up to t ~ 140) to yield statistical properties of the steady-state flow for several diameters downstream of the orifice. Further detail, including results for sinusoidally-modulated jets, may be found in Ref. [22]. 2. A striking qualitative feature of the model jet flow is the formation of almost-regular clusters of vorticity; the average length scale of the clusters increases (by amalgamation) as they are swept downstream. In this respect the idealized axisymmetric jet resembles laboratory jets. 3. The mean spread rate of the initial mixing layer (and the corresponding entrainment over the first 4 dia.) are quite similar in the model and real jets, despite the absence of viscosity in the model simulation. In other respects, however, the steady-state model flow differs significantly from real turbulent jets. Turbulence intensities (axial and radial) in the shear layer are higher by a factor of 2, and two-point velocity correlations reveal a greater persistence of eddy structure with downstream separation than is measured in the laboratory. Similar findings, albeit based on rather short data samples (To < 45) have been reported by Acton[13]. Ashurst[23] has simulated a two-dimensional mixing layer with line vortex elements and also found (u,,,,, v,,,s) values of order 20-25% but was able to reduce these by allowing each line vortex to spread exponentially. 4. Further differences between the axisymmetric model flow and real turbulent jets are apparent in Fig. 4, which shows the u(t) time history near the inner edge of the initial mixing layer. The predominance of negative spikes in the laboratory data is not matched in the model, although the version with 8t = 0.2 (a time increment too large for accurate convergence at this position) is closer to the measured trace. We conclude that if the negative spikes measured by Lau et a/.[l'6] are indeed accounted for by an "axisymmetric vortex structure", as Lau and Fisher[24] have suggested, the dynamics of such a structure must differ significantly from the strictly axisymmetric inviscid model investigated here. 5. Possible future developments of the jet model include: (a) Addition of a co-flowing outer stream. (b) Addition of a plane boundary normal to the jet axis, to simulate an impinging jet. (c) Allowing the jet and ambient fluids to have different densities. (d) Allowing the jet circulation to decay downstream, e.g. by removing vortices whose radius falls below a limiting value. However, the present work provides a caution against transferring conclusions from axisymmetric turbulent models to real jet flows.
220
A.V.J. EDWARDSad C. L. MORFEY
Acknowledgements--The present study was initiated by Prof. P. O. A. L. Davies, whose guidance during the first author's Ph.D. studies is acknowledged. Financial support was provided by the Science Research Council and the University of Southampton. REFERENCES 1. G. L. Brown and A. Roshko, On density effects and large structure in turbulent mixing layers. J. Fluid Mech. 64, 775-816 (1974). 2. S. C. Crow and F. H. Champagne, Orderly structure in jet turbulence. J. Fluid Mech. 48, 547-591 (1971). 3. A. Roshko, Structure of turbulent shear flows: a new look. AIAA J. 14, 1349--1357 (1976). 4. G. K. Batchelor, Computation of the energy spectrum in homogeneous two-dimensional turbulence. Phys. Fluids 12, II.233-II.239 (1969). 5. C. L. Morley, Sound radiation from a simulated jet flow. In Mechanics of Sound Generation in Flows (Edited by E.-A. Mfiller), pp. 12-18. Springer-Verlag, Berlin (1979). 6. P. O. A. L. Davies and J. C. Hardin, Potential flow modelling of unsteady flow. In Numerical Methods in Fluid Dynamics (Edited by C. A. Brebbia and J. J. Connor). Pentech Press, London (1974). 7. A. J. Chorin and P. S. Bernard, Discretization of a vortex sheet with an example of roll-up. J. Comput. Physics 13, 423-429 (1973). 8. R. R. Clements, An inviscid model of two-dimensional vortex shedding. J. Fluid Mech. 57, 321-336 (1973). 9. R. R. Clements and D. J. Maull, The representation of sheets of vorticity by discrete vortices. Prog. Aerospace Sci. 16, 12%146 (1975). 10. R. R. Clements, Flow representation, including separated regions, using discrete vortices. In Cutup. Fluid Dynamics, AGARD LS-86 (1977). 11. D. J. Mautl, An Introduction to the Discrete Vortex Method. Cambridge University Engineering Daptartment. CUED/A-Aero/TR8 (1979). 12. H. Lamb, Hydrodynamics, 6th Edn, pp. 202-249. Dover, New York (1945). 13. E. Acton, A Modelling of Large Jet Eddies. Ph.D. Thesis, Cambridge University (1977); see also J. Fluid Mech. 98, 1-31
(1980). 14. C. L. Morley and A. V. J. Edwards, Energy and Momentum Constraints on the use of a Cutoff in Lumped-vortex Flow Models. University of Southampton, ISVR Memorandum 587 (1978). 15. S. M. Datums and D. Kiichemann, On a vortex-sheet model for the mixing between two parallel streams. Proc. R. Soc. Lurid. A 339, 451--461 (1974). 16. J. C. Lau, M. J. Fisher and H. V. Fuchs, The intrinsic structure of turbulent jets. J. Sound Vib. 22, 379-.406 (1972). 17. P. Bradshaw, D. H. Ferriss and R. F. Johnson, Turbulence in the noise-producing region of a circular jet. J. Fluid Mech. 19, 591-624 (1964). 18. B. Bose, Experimental Study o[ the Wave-like Characteristics of Turbulence in Shear Layers and Behind Grids. Ph.D. Thesis, Southampton University (1968). 19. M. J. Fisher and P. O. A. L. Davies, Correlation measurements in a non-frozen pattern of turbulence. J. Fluid Mech. 18, 97-116 (1963). 20. B. J. Hill, Measurement of local entrainment rate in the initial region of axisymmetric turbulent air jets. J. Fluid Mech. 51,773-779 (1972). 21. L Maestrello and E. McDaid, Acoustic characteristics of a high-subsonic jet. AIAA J. 9, 1058-1066 (1971). 22. A. V. J. Edwards, The Computer Simulation of an Unsteady Flow. Ph.D. Thesis, Southampton University (1978). 23. W.T.Ashurst, Numerical simulation of turbulent mixing layers via vortex dynamics. In Turbulent Shear Flows--l (Edited by F. Durst, B. F. Launder, F. W. Schmidt and J. H. Whitelaw), pp. 402-413. Springer-Verlag, Berlin (1979). 24. J. C. Lau and M. J. Fisher, The vortex-street structure of "turbulent" jets. J. Fluid Mech. 6'/, 299-337 (1975). APPENDIX. R E L A T I O N B E T W E E N E N T R A I N M E N T AND D I F F U S I O N OF MEAN VELOCITY IN A JET For a monotonic mean velocity profile ti(y), a round jet (Fig. 15a) has its mass flow distributed with respect to li according to . a~ - i
P(~)=yu ~yy
/fo
®
(~
(o)
y
(AI)
P(;)
//~
Increasing x
II
IN.'~I
. Omox
0
mx
Fig. 15. (a) Mean axial velocity profile in jet. (b) Amplitude density of mean velocity with respect to mass flow rate.
A computer simulation of turbulent jet flow
221
In other words, P(a) da is the proportion of the mean flow rate (across a plane x = const.) for which the mean velocity lies between (a, a + da), thus
oP(a) da = I. Let Q(x) and X(x) denote the total mean flux of mass and x-momentum across
(A2) x = const., then
dX = a dQ (pressure terms and u 2 - (a) 2 neglected) = aQ. P0i) da.
(A3)
Integration gives
f~
X(x)/Q(x) = J o aP(a) dfi.
(A4)
Momentum conservation implies X(x)= const.; diffusion of the mean velocity profile (Fig. 15b) implies a decrease in the integral on the right; thus, if the profile spreads, Q(x) must increase with x, meaning that fluid is entrained by the jet.