A~tr&--A aurner~c~ study of fairy-deve~o~d turbtdent flow past a ba~ward-fauns step is performed to analyze the e&et of mesh ~~~nerne~t on the computed results. The time averaged equations of conse~a~ou of mass and momentum are solved by a fi~te-volume method using two versions of the K-E model af turbulence. The computations are performed with both the standard and the non&ear K-E turbulence models for 166 x 73 and 332 x 145 mesh points. The comparkon of the results with available ex~~me~ta~ findings indicate that no&near terms must be ittcorporated into the K-s model to account far normal stress differences, and that even with very fine ~so~ution, the standard K-E model fails to provide accurate predictions elf the flow field.
~~pa~~~i~~ of a turbulent Sunday layer and its sub~e~~e~t reattachment over a sohd surface occurs in a variety of e~gi~eer~~g applications (e.g. di~users, Aow over ~rcr~ wings, base flow over coasting missiles, etc.). Therefore, u~dersta~~~g the ch~a~e~stics of such flow fields is crucial to the design and performance analysis of these e~g~~eer~~gdevices. hong the few selected geometries that have been studied in detail is the Bow past a backward-facing step, It is often used as a test-case for analyzing the efficacy of ~~rnp~tat~~nal algorithms and turbulence modeis since it embodies several crucial aspects of turbulent separated flows-namely flow separation, reattachment, and the occurrence of secondary separation regiuns. Xn addition, for the case of Bow over a backyard-faei~g step, the separation point is tied at the sharp corner of the step, thus al~~w~~gthe exarn~~atjon of the physicai processes ~thout the additional comp~e~~t~esresulting from the movement of the separation point. A schematic of the basic ffow field is sbuwn in Fig. 1. The flow separates at the step corner and reattaches at a distance, L, from the inlet. This distance is a function of the remold number ideated here based on the upstream cc~ter~ine mean velocity and downstream channel width) and expansion ratio (defined as the ratio of downstream to upstream channel heights, ~~~~~~).The analysis of the backwards facing step t~~ca~~yinvolves the modeling of the flow field for a given expansion ratio as a function of the fiow Reynolds number. For very large vdues of the Reynolds number, the flow in the channel is fuily turbulent, and is within the realm of e~gi~ee~ng app~~cati~~s~Several recent works [l--3] have considered the fiow past a backward~fac~~gstep and provide an exceflcnt review of the current state of the art from both ex~e~me~ta~ and computational paint of view. Among these, the recent work by Speziale and Nga [I] deals with an e@icient ~pp~oa~b to ~~rnp~t~tio~al analysis of fully tnrbu~ent flows past a backwards-facing step with particular emphasis on t~rbui~n~ modeii~g and the predictive capabiijt~es of a new nonlinear K-E model. The work by ~e~en~igi~and Nellor [Z] also consider the flow past a backward~fac~ng step from a computational point of view but with a full Reyndds stress closure model, The experimentat aspects of fully turbulent Bow past a ba~k~rd~fa~i~g step is addressed by Kim et a5. f3f. In additive to the case of fully turbulent flow, the transition from Iaminar to turbn~~~~flow past backward-facing steps have also been of considerable interest to the researchers, and the interested readers are referred to the comprehensive work on the laminar and turbulent flow past a backward-facing step by ~m&y
me
ef ai. f41. Q% present work addresses another facet of this interesting problem, earnest, the effect of grid r~s~~~lion on the computed results. FOTthis purposet the work by Speriate and Ngo [If which considers the funny-developed turbulent flow past a channel with a backw~rd-~a~~g step with an expansion ratio of 32 at a Reynolds number of 132, is chosen as the burns 607 ES ?Q:i;-F
608
S. THANGAM
and N. HUR
Fig. I. Physical configuration and coordinate system
reference. In this work the computations were performed with 166 x 73 nonuniform finite volumes using both the standard K-E model and a nonlinear K-E model [.5]. Their results clearly demonstrate the efficacy of the nonlinear K-E model in the prediction of turbulent separated flows. In the present study, the adequacy of the grid resolution for both turbulence models is examined based on computations with a 332 x 146 mesh. It is shown that the effect of even quadrupling the number of mesh points is minimal for both the standard and the nonlinear models. In the sections to follow, the formulation of the physical problem and the method of solution are given along with a discussion of the results and the conclusions.
2. FORMULATION
OF
THE
PHYSICAL
PROBLEM
AND
THE
METHOD
OF
SOLUTION
The physical problem under consideration is the fully-developed turbulent flow of an incompressible viscous fluid past a backward-facing step. A schematic of the physical configuration and the coordinate system are specified in Fig. 1. The governing equations are the incompressible, time-averaged Navier-Stokes and continuity equations:
where p is the fluid density, U and ii are the mean velocities in the x and y directions, p is the Based on the mean pressure, and TX,, t,, and rYYare the Reynolds stress components. nonlinear K-E model proposed by Speziale [5], the Reynolds stress tensor from which Z,,, 7,, . and ?YYare determined is given by
where, iii = (ii, 6) is the mean velocity vector and
are the mean rate of strain tensor and its frame indifferent Oldroyd derivative, respectively. In (4), K is the turbulent Kinetic energy, E is the turbulent dissipation rate, and C, and Cn are
Turbulent separated Bow
609
empirical constants that assume the values of 0.09 and 1.68, respectively [S]. The turbulent kinetic energy and the dissipation rates are obtained separately from the modeled versions of the transport equations which assume the following simplified form for the two-dimensional problem at hand.
where, LIDand aK are the turbulent Prandtl numbers which are taken to be 1.3 and 1.0, respectively, while C,, and C,, are dimensionless constants which are taken to be 1.43 and 1.92, respectively [6]. It should be noted here that (5) and (6f are essentially the same as those used by Speziale and Ngo [l] except for the inclusion of the turbulent diffusion terms. The mean velocity field at the channel inlet is specified as that for fully-developed turbulent channel flow, while the outflow conditions are extrapolated. The channel length is selected to be sufficiently far upstream and downstream of the step to ensure fully-developed turbulent flow conditions. The law of the wall is applied to all solid boundaries (a single log-waI1 layer starting at Y+ = 11.6 is used; for details about the actual implementation of the boundary conditions, see Lilley and Rhode [9] and Speziale and Ngo [l]). The governing equations (l)-(3) are discretized along with the transport equations for the turbulence kinetic energy and dissipation rate (5) and (6) using a finite-volume method [?I. In this scheme the velocities are defined at the boundaries of the control volume, while the pressure, turbulent kinetic energy and dissipation rate are defined at the centroid of the control volume. The discretized equations for all the field variables are obtained by integration over the control volume. The resulting system of nonhomogeneous algebraic equations are solved by a successive line under relaxation (SLUR) procedure with the repeated application of the tridiagonal matrix solution [8]. In the present work, a finite-volume computational code developed by Lilley and Rhode [9] which uses the standard K-E turbulence model was modified to incorporate the nonlinear K-E model given by (4); the standard K-E model is recovered in the limit as CD -+ 0. The details of the modifications are given elsewhere [I], and the following steps briefly explain the solution procedure: (4 A-set of pseudo-velocity values are obtained by solving the momentum equations by assuming uniform pressure in the computational domain. values are then used to obtain pressure by solving the @I These pseudo-vel~ity discretized equation for the conservation of mass with the SLUR method. (cl Based on the pressure obtained in (b), the velocity components are obtained from the discretized momentum equations with the SLUR method. (4 The correction for velocity components are obtained by evaluating the correction for the mass flux in each cell. (4 Using the corrected velocities from (d), the transport equations for the turbulence kinetic energy and dissipation are solved using the SLUR method. The stress values are then updated, and based on the velocity values obtained from (d), the procedure is repeated until convergence is obtained. The resulting mean velocity field can then be integrated to obtain the form of the streamlines. The imputations reported in this work were performed using nonuniform mesh and consist of two categories: a 166 X 73 mesh to emulate the results of Speziale and Ngo [l] and a 332 x 146 mesh to evaluate the effect of mesh refinement. Both the standard and the nonlinear K-E model were considered to determine the improvements from further mesh refinements. The distribution of the mesh was selected based on that proposed by earlier works [l, 21, with
610
S. THANGAM
and N. HUR
the finer mesh distributed near the step corner and near walls. For simplicity, the 332 x 146 mesh was generated directly by doubling in each direction on the mesh used by Speziale and Ngo [l]. In addition, to improve the accuracy, all the second derivatives were computed using a fourth order scheme inside the recirculation zone since this was found necessary for the elimination of numerically induced oscillations in the flow field (Speziale and Ngo [l]). The computed solution was assumed to have converged to its steady state when the root mean square of the average difference between successive iterations was less than 1O-4 for the mass source [7,9] and the reattachment point which in turn was obtained by interpolation. Additiona details regarding the computational scheme may be found elsewhere [l, 9]. A DE/VAX-8700 with a processing speed of 6 VUP (VAX Units of Performance with the VAX 11/780 being 1 VUP) was used for the calculations with 64 bit arithmetic. The computations required approximately 240 minutes of CPU time for 1000 iterations for the 166 x ‘73 mesh using the standard K-E model and about 480 minutes for the same number of iterations involving the nonlinear K-E model. The time requirement for the 332 x 146 mesh was about four times that for the 166 x 73 mesh. To meet the specified convergence criteria, at Ieast 1800 iterations were needed for the standard K-E model, while about 2700 iterations were needed for the nonlinear K-E model.
3. RESULTS
AND
DISCUSSION
The main aim of the present work is to examine the effect of grid refinement on the computed solutions for fully turbulent flow past a backward-facing step for a particular expansion ratio and Reynotds number. As indicated earlier, a 3:2 expansion ratio (defined as the ratio of downstream to upstream channel heights, HZ/H,) and a Reynolds number (based on the upstream centerline mean velocity and downstream channel width) of 132,000 were selected to facilitate comparison with the results already available in the literature [I, 31. Since the computational work by Speziale and Ngo [l] found that the a significant improvement in the computed results could be obtained by the use of nonlinear K-E turbulence model in comparison to the standard model, in the present work the effect of grid refinement on both the standard and the nonlinear models are considered. The computed results are presented in Figs 2-5 for the following four cases: a 166 x 73 non-uniform mesh with the standard 1y-E model; a 332 x 146 nonuniform mesh with the standard K-E model; a 166 x 73 nonuniform mesh with the nonlinear K-E model; and a 332 x 146 nonuniform mesh with the nonlinear K-E model. For brevity, only the streamlines, mean velocity profiles and the streamwise turbulence intensity profiles will be presented for each case and compared with the experimental results of Kim et al. [3]. It should be noted here that the rest&s for the 166 x 73 mesh with the standard and nonlinear K-E models arc essentially the same as that given in [l], except that the influence of the turbulence diffusion terms have been included, and a slightly modified algorithm with an improved convergence criteria has been used for the computations. They are presented here primarily for comparison purposes. In Figs 2 (a)-(c), the results for the computations with a 166 x 73 mesh based on the standard K-E model are given. In Fig. 2 (a), computer generated streamlines are shown, and the mean reattachment point for this case is computed to be 5.57. As discussed in Speziale and Ngo (l), this compares qualitatively with the experimental value of about 7.0 [3]. In Fig. 2 (b) profiles of the dimensionless mean velocity, ii/U, are shown at selected streamwise locations along with the available experimental data from Kim et al. 131. This is followed in Fig. 2 (c) with the profiles of turbulence intensity (%)““/L& at selected axial locations. As can be seen for both these cases, the comparison of the computed and experimental results are, at best. qualitative. In Figs 3 (a)-(c), the results for the computations with 332 x 146 mesh based on the standard K-E model are presented. In Fig. 3 (a), computer generated streamlines are shown, and the mean reattachment point for this case is computed to be 5.58. This should again be compared with the experimental value of about 7.0 131. In Figs 3 (b) and 3 (c) profiles of the
611
Turbulent separated flow
1.a
1.0
0.0
1.0
0.0
1.0
0.0
0.0
1.0
* Experimental Data (Kim, Kline and Johnston, 1980); -Computational
1.0
0.0
Results
(W
o
X/AH-l.333
X/At4
n.tJo
* 2.867
D.fS
x/An
0.00
-
ats
5.333
X/&i
o.uo
-
O.lf
7.887
%/At4 = 8.533
0.00
0.15
x&&i
0.00
-
10.3%
0.15
(iiii)‘“lu* G Experimental Data (Kim, Kline and Johnston, 1980); -Computational
Results (note: z
= -z&p)
Fig. 2. Standard K-E model: 166 X 73 mesh. (a) Computer generated streamlines (reattachment length, L = 5.57). (b) Dimensionless mean velocity proties at selected locations. (c) Dimensionless turbulent intensity profiles at seiected locations.
dimensionless mean axial velocity, tz/UO, and profiles of the turbulence intensity (~)~~/U~ at selected axial locations are shown along with the available experimental data from Kim et al, [3]. As can be seen, there is little improvement even though the number of mesh points is four times that used before, We now consider the results based on the nonlinear K-E model. In Figs 4 (a)-(c), the results for the computations with a 166 X 73 mesh are given; Figs 5 (a)-5 (c) contain the results for 332 x 146 mesh, As for the standard K-E model, the computer generated streamlines are shown
S. THANGAM
0
x/An
-
1.333
X/AH
-
2.867
X/AH
-
and
5.333
N. HUR
X/AH
.
~~~~~~~
0.0
1.0
0.0
1.0
0.0
1.0
0.0
1.0
0.0
1.0
0.0
1.0
U/U, O Experimental
0
x/AH
-
1.33
3
X/AH
Data (Kim, Kline and Johnston,
-
2.887
X/AH
-
5.33
1980); -Computational
x/An
-
7.007
x/AH
Ii
gR
.
2
.
.
. i
uo.oo
1 . .
1 .
IF 0 ,: AK
0.15
0.00
0.15
0.00
LIL 0.00
0.15
8.55
. .
0
-
Results
0.15
0.00
0.15
0.00
0.1s
(ii.y2/U, q
Experimental
Data (Kim, Kline and Johnston,
Fig. 3. Standard K-E model: 332 length, L = 5.58). (b) Dimensionless turbulent
x
1980); -
Computational
Results (note: uu = --T.&)
146 mesh. (a) Computer generated streamlines (reattachment mean velocity profiles at selected locations. (c) Dimensionless intensity profiles at selected locations.
in Fig. 4 (a). The mean reattachment point for this case is computed to be 6.28 which compares more favorably with the experimental value of 7.0 [3]. In Figs 4(b) and 4(c) profiles of the dimensionless mean velocity and the dimensionless streamwise turbulence intensity are shown at selected locations along with the available experimental data from Kim ef al. [3]. As can be seen, much better comparison with the experimental results are evident for both the mean velocity profiles and the turbulence intensity in contrast to the standard K-E model. This is attributed primarily to the more accurate representation of the normal stress differences in the
613
Turbulent separated flow
X/AH 64
I
d.0
i.0
ii.0
Lo
i.0
6.0
i.0
d.0
d.0
l.0
I;,U, D Experimental Data (Kim, Kline and Johnston, 1980); -Computational
Results
@I
X/AH - 5.333
2l7iTTdiL-i 0.00
0.15
0.00
x/An - 7.067
0.00
0.15
0.16
x/An - 0.5:
0.00
0. 1s
0.00
0.1s
(iiii)‘“/U, a Experimental Data (Kim, Kline and Johnston, 1980); -
Computationai Results (note: G = -z,,/pf
(Cl Fig. 4. Nonlinear K-E model; 166 x 73 mesh. (a) Computer generated streamlines (reattachment Length, L = 6.28). (b) Dimensionless mean velocity profiles at selected locations. (c) Dimensionless turbulent intensity profiles at selected locations.
nonlinear model (4), as pointed out by Speziale [5]. The computer generated streamlines based on the finer 332 x 146 mesh are shown in Fig. 5 (a), and the corresponding mean reattachment point is 6.30. The dimensionless profiles of the mean axial velocity and the turbulence intensity are shown in Figs 5 (b) and 5 (c). Again the agreement between the computational and the experimental results are quite good, though the improvement due to the quadrupling of the mesh points is observed to be minimal. In summa~, the results obtained with the 166 x 73 mesh are found to be adequate for a
614
S. THANGAM
and N. HUR
0
10
X/AH
0
fn X/Art-l.
ti
0
Pi
0
A
0
d
0.0
0
0.0
i. 0
d. 0
i.0
i.0
d-0
d.o
i.0
d.0
i.
0
;;,U,
0 Experimental Data (Kim, Kline and Johnston, 1980); -Computational
Resutts
(hl
Xl&i
-
5.x
- 7.86
X/AH
I3
-I’
. !?
h IT----
l
uo.oo
0.15
O.OG
0.15
0.00
n.
15
J ‘,
n.on
~~)“*/u*
n Experimental Data (Kim, Kline and Johnston, 1980); -
0.15
0.00
0. 15
i3
e X/AH -
10.33.
1’ l
LI
I
i
0.00
0.15
Computational Results (note: E = -qp)
(C) Fig. 5. Nonlinear K-E model: 332 x 146 mesh. (a) Computer generated streamlines (reattachment length, L = 6.30). (b) Dimensionless mean velocity profiles at selected locations. (c) ~imensionIess turbulent intensity profiles at selected locations.
channel with a 3:2 expansion ratio for both the standard as well as the nonlinear K-E turbulence models. In both cases, additional mesh refinement have shown to provide only minimal improvements in the computed solutions. However, as shown by Speziale and Ngo [lf , a substantial improvement can be obtained by the use of the nonlinear K-E turbulence mode1 instead of the standard model. It should also be noted that while the computational time is increased nearly by a factor of 1.5 for the nonlinear model over the standard model, it is still within the realm of mini-computers with average computing speeds.
Turbulent separated flow
615
4. CONCLUSIONS A highly resolved numerical study of fully-developed turbulent flow past a backward-facing step has been performed using the standard and the nonlinear K-E models of turbulence. The incompressible, time-averaged Navier-Stokes and continuity equations along with the transport equations for the turbulence kinetic energy and dissipation rate were derived to include the turbulence diffusion terms and solved by a finite-volume method. The computations were performed using a 166 X 73 and a 332 X 146 mesh resolution for each of the two turbulence models under consideration. The following conclusions were arrived at: (9 The addition of the turbulence diffusion terms and the implementation of substantial mesh refinement de~nitively establishes that the nonlinear K-E model yields considerable improvement in the prediction of the mean reattachment point to within 10% of the average experimental value. (ii) The comparison of the results with available experimental data confirm the findings of Speziale and Ngo [l] that both the mean velocity profiles and the turbulence intensities are also better predicted by the nonlinear K-E model. (iii) The results indicate that the 166 x 73 mesh used by Speziale and Ngo [l] is quite adequate, and that further mesh refinement yields only a minimal improvement. In summary, these results based on a single log-wall layer model (with yf = 11.6) definitively establish the overall superiority of the predictive capabilities of the nonlinear K-E turbulence model [5] over the standard K-E model in so far as the backward-facing step problem is concerned. However, it should be noted that different boundary conditions, use of multiple wall /ayer models, or the inclusion of low Reynolds number near wall corrections can significantly alter these findings. A more detailed analysis of these effects will be the subject of a future study. Ac/cnowledgemen&--The authors are indebted to Dr Charles G. Speziale of ICASE, NASA Langley Research Center for suggesting this problem, and for his thoughtful criticism and encouragement during the course of this study. The authors also acknowledge the support of the Stevens Institute of Technology and its Computer Center.
REFERENCES [l] C. G. SPEZIALE and T. NGG, ht. J. Enging Sci. 26, 1099 (1988). [2] M. C. CELENLIGIL and G. L. MELLOR, ASME J. Flus Enging 107,467 (1978). 131J. KIM, S. J. KLINE and J. P. JOHNSTON, ASh4.E J, Fluids Enaina 102.302 (1980). [4j B. F. ARMALY, F. DURST, J. C. F. PEREIRA and B. SCHGNUNG, i Fluid M&h. 127,473 (1983). ‘51 C. G. SPEZIALE, J. Fluid Mech. 178.459 (1987). j6] B. E. LAUNDER and D. B. SPALDING, kath;matical Models of Turbulence. Academic Press, London (1972). 71 S. V. PATANKAR, Numerical Heat transfer and Fluid Flow. Hemisphere, New York (1980). ;S] E. ISAACSON and H. B. KELLER, Analysis of Numerical Methodk Wiley, New York (1970). ,9] D. G. LILLEY and D. L. RHODE, A computer code for swirling turbulent axisymmetric recirculating flows in practical isothermal ~mbustor geometries. NASA Contractor Report CR-3442 (1982). (Received 24 July 1990; accepted 21 August 1990)