Direct numerical simulation of incompressible turbulent flows

Direct numerical simulation of incompressible turbulent flows

Computers & Fluids 30 (2001) 555±579 www.elsevier.com/locate/comp¯uid Direct numerical simulation of incompressible turbulent ¯ows R. Friedrich a,*,...

261KB Sizes 0 Downloads 21 Views

Computers & Fluids 30 (2001) 555±579

www.elsevier.com/locate/comp¯uid

Direct numerical simulation of incompressible turbulent ¯ows R. Friedrich a,*, T.J. H uttl a,1, M. Manhart a, C. Wagner b a

Lehrstuhl f ur Fluidmechanik, TU Munchen, Boltzmannstrasse 15, 85748 Garching, Germany b DLR, Bunsenstrasse 10, 37073 G ottingen, Germany Accepted 10 November 2000

Abstract The paper discusses recent achievements of direct numerical simulation (DNS) of incompressible ¯ows. The various spatial discretization techniques which can be used in the case of simple or complex geometry are referred to, along with suitable time advancement schemes. The advantage of using a staggered variable arrangement and ecient Poisson solvers is stressed before initial and boundary conditions for in¯ow, entrainment and out¯ow boundaries are discussed. Comments on the spatial resolution required for DNS and on the use of zonal embedded grids to reduce the computational e€ort complete the part on numerical aspects. A few examples are then given where DNS serves to study the physics of turbulent ¯ows. The ®rst is fully developed ¯ow through solid and permeable pipes. Wall permeability turns out to be a means to model wall roughness. Striking changes in the turbulence structure are observed near rough walls. The second example deals with ¯ow separation in sudden expansions of solid and permeable pipes. It is found that the free shear layer surrounding the re-circulation zone experiences a markedly di€erent downstream development if the incoming ¯ow is con®ned by rough walls. The reattachment length is reduced. At the same time the ampli®cation rates of the maximum values of the Reynolds stress components in the mixing layer are lower than those in the solid wall con®guration. The results underline the importance of the wall blocking e€ect. The third example concerns fully developed ¯ow through curved and coiled pipes. Strong secondary ¯ow generates complex mean ¯ow ®elds and six nonzero Reynolds stress components contributing to the production of turbulent kinetic energy. This renders the ¯ow unpredictable by two-equation models and thus a challenge to turbulence modellers. Ó 2001 Elsevier Science Ltd. All rights reserved.

*

Corresponding author. Tel.: +49-89-289-16144; fax: +49-89-289-16145. E-mail addresses: [email protected] (R. Friedrich), [email protected] (T.J. H uttl). 1 Also corresponding author. Address: Fichtenstrasse 1, D-82178 Puchheim, Germany. Tel.: +49-89-890-20957; fax: +49-89-1489-95242. 0045-7930/01/$ - see front matter Ó 2001 Elsevier Science Ltd. All rights reserved. PII: S 0 0 4 5 - 7 9 3 0 ( 0 1 ) 0 0 0 0 6 - 8

556

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

1. Introduction Turbulence in ¯uid ¯ows is a phenomenon that cannot be de®ned easily. However, common characteristics of all turbulent ¯ows are: irregularity in time and space, ecient mixing, nonlinearity, large Reynolds numbers, dissipation and continuous spectra of length and time scales. Principally, there are three di€erent approaches to predict turbulent ¯ows: · Statistical turbulence modelling (STM) · Large eddy simulation (LES) · Direct numerical simulation (DNS). In this order they form a hierarchy of techniques with their potential and level of information as classifying principle. The statistical approach to turbulence (STM) is associated with the highest loss of information and with a closure problem which is still not satisfactorily solved. On the second-order level, the main problems are associated with the modelling of the pressure±strain rate and the dissipation rate tensors. Existing models bear many similarities with models for rheological ¯ows. The principle of fading memory often used there is fundamental for STM. Spectral information is completely lost, since any statistical quantity is an average over all turbulence scales. The LES technique is intermediate between STM and DNS since it allows to predict the dynamics of the large turbulent scales. Only the e€ect of the smallest scales on the resolved ones is modelled, in the hope that the small scales are less anisotropic and less in¯uenced by the geometry of the ¯ow than the large ones. The three-dimensional time-dependent low-pass ®ltered transport equations are solved numerically. The ®ltering process introduces subgrid-scale (SGS) stresses (and SGS scalar ¯uxes) which `¯uctuate' in space and time and describe the transfer of momentum (and scalar) across the cut-o€ wave number (resp. frequency) from large to small scales and back (backscatter e€ects). Most of the currently available SGS models are purely dissipative, i.e. they fail to predict backscatter phenomena which appear in wall layers, mixing layers or generally in layers with high gradients of ¯ow variables. This de®ciency does not matter if errors in the prediction of turbulent ¯ows up to 20% can be tolerated. From a scienti®c point of view, however, accurate SGS models have to be developed. LES is a technique to predict high-Reynolds-number ¯ows. DNS considers all degrees of freedom appearing in the ¯ow and is thus the only way to predict and analyse turbulent ¯ows in all details without using any approximating assumption. Since the required number of grid points increases faster than the square of the Reynolds number, DNS is at present feasible only at low or moderate Reynolds numbers. DNS is a tool to study turbulence physics, to test new statistical or SGS models and to predict complex ¯ows of technical interest. It has been applied to various ¯ows in the past. Without trying to be complete, we mention some of the more recent fundamental studies in which earlier investigations are quoted. Direct simulations of homogeneous isotropic turbulence were performed by Vincent and Meneguzzi [70]. They found that vortical structures are organized in tubes of length comparable to integral scales and diameter of a few Kolmogoro€ length scales. DNS of homogeneous shear turbulence has been reported e.g. by Lee et al. [42] for values of the mean shear rate comparable to those in wall layers. The authors could demonstrate that streaky

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

557

structures appear also in unbounded ¯ows forming an array of high- and low-speed regions which develop as a consequence of high shear rate rather than of wall blocking. Mixing layers were directly simulated by Lesieur et al. [44], Moser and Rogers [52], showing that coherent vortices undergo helical pairings if three-dimensional forcing of the large scales is applied and identifying vortex pairing as one possible route to turbulence. Many direct simulations refer to wall bounded ¯ows and have provided useful databases for analysing near-wall e€ects. The classical work of Kim et al. [37] related to fully developed channel ¯ow has initiated a bulk of further investigations and has even been used to test the accuracy of experimental techniques. Only recently Eckelmann [11] could clarify the deviation of his hot-®lm data obtained in an oil channel from DNS data. Channel ¯ow DNS has been used by Choi et al. [7] to develop techniques for active turbulence control. Since sweep and ejection events (inrush of high-speed ¯uid towards the wall and transport of low-speed ¯uid away from the wall) are mainly responsible for the high skin friction of turbulent ¯ows, suppression of these events by local wall blowing and suction (controlled by sensors and a feedback loop) turned out to be principally e€ective, but inapplicable in practice. Recently, Lim et al. [45] demonstrated by DNS the e€ectiveness of spanwise magnetic ¯uxes in inhibiting near-wall coherent motions in electrically conducting ¯uids. Injection or suction of a ¯uid normal to a wall which are also used to control turbulent shear ¯ow has been investigated by Sumitani and Kasagi [66] based on DNS of channel ¯ow with heat transfer. Using DNS data for fully developed channel ¯ow, Antonia and Kim [2] studied low-Reynolds-number e€ects on near-wall turbulence. DNS data of fully developed channel ¯ow were produced by H artel and Kleiser [22] in order to analyse and model SGS motions in wall layers which turn out to be considerably controlled by backscatter e€ects. Moser and Moin [51] used DNS to analyse mild e€ects of streamline curvature in fully developed channel ¯ow and found Taylor±G ortler vortices to be responsible for half of the di€erences in the Reynolds and the wall shear stress between convex and concave sides. Perot and Moin [56,57] have isolated the viscous and kinematic mechanisms of the wall/turbulence interaction and have developed improved near-wall models for the dissipation, pressure±strain and turbulent transport terms of the Reynolds stress equations, without resorting to additional model constants or ad hoc damping functions. DNS data of fully developed ¯ow in pipes with circular cross-section are discussed by Eggels et al. [12]. They were obtained with second-order central discretizations and compared well with LDA measurements. Secondary motions of Prandtl's second kind appearing in straight channels with square cross-section were investigated by Gavrilakis [17] using second-order central di€erencing on staggered grids and high resolution. Huser and Biringen [33] performed DNS of the same ¯ow but with ®fth-order upwind-biased ®nite di€erences for the convection terms and fourth-order central di€erences for the viscous terms. A comparison of both computations is dicult because of di€erences in Reynolds number and spatial resolution. It appears, however, that the higher-order accuracy of Huser and Biringen's scheme cannot compensate their rather poor spatial resolution. Developing turbulent ¯ows like boundary layers have been the subject of DNS since the 80s. Respectable DNS data have been presented by Spalart [64,65] for self-similar and zero-pressure gradient boundary layers. The separating turbulent boundary layer along a ¯at plate was recently directly simulated by Na and Moin [53] based on second-order central di€erence schemes. The authors developed useful methods of getting approximate but physically realistic in¯ow and freestream boundary conditions. Boundary layers with convex transverse curvature developing on

558

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

long cylinders were investigated by Neves et al. [54,55]. Curvature e€ects were identi®ed by comparison with plane channel data. Several turbulence statistics were found to scale with a curvature dependent velocity scale. Applications of DNS to ¯ow in complex geometries are less numerous. A classical example is the ¯ow over a backward-facing step. Recently a turbulent boundary layer separating at a backstep and re-attaching has been investigated by Le et al. [41] with a second-order central di€erence scheme on staggered grids providing excellent agreement with statistical data of a concurrent experiment. A related problem, namely the turbulent ¯ow in a sudden pipe expansion was studied by Wagner and Friedrich [72] and Wagner [71], using similar numerical techniques. Turbulent ¯ow over surface-mounted riblets was simulated by Choi et al. [6] in generalized coordinates in the wall-normal and spanwise directions, but cartesian coordinates in the streamwise direction. Spatial discretizations were again second-order central. The computed drag turns out to be in good agreement with existing experimental data. The analysis helped to elucidate the mechanism of drag reduction by riblets. More applications of DNS could have been listed. The reader is referred to review articles of Rogallo and Moin [61], Reynolds [59], Kasagi and Shikazono [36], Hartel [21], Moin and Mahesh [50] and to the article of Kleiser and Zang [38] about DNS of transitional ¯ows. Progress in the future will e.g. be reported in the biannual conference series `Turbulence and Shear Flow Phenomena', `European Turbulence Conference', `ERCOFTAC Workshop on Direct and LargeEddy Simulation' and `AFOSR International Conference on DNS and LES'. 2. Numerical methods 2.1. Spatial discretization Accurate numerical methods are required for reliable DNS. Since turbulent ¯ows generally have a broad range of scales which increases with the Reynolds number, all scales, even the smallest have to be treated as accurately as possible. Methods with signi®cant arti®cial viscosity should be avoided. Excellent work has been done in the past using spectral methods, which provide very accurate spatial di€erentiation. In these methods the ¯ow variables are expanded in terms of smooth, global and mostly orthogonal functions (trial functions). The word `global' refers to trial functions that are de®ned on the whole computational domain. `Smooth' means in®nitely di€erentiable. The expansion coecients are determined via a Galerkin, tau or collocation method. In the latter, the coecients are only used to perform analytical di€erentiation, hence the name pseudospectral method. The collocation is the simplest and less time consuming of the three methods. Moreover, boundary conditions can be more easily implemented since the nodal unknowns are explicitly used. For domains with periodic boundary conditions, Fourier series are the usually chosen trial functions. Jacobi polynomials (e.g. Chebyshev and Legendre polynomials) are appropriate in situations with nonperiodic boundary conditions. For unbounded domains a mapping of the domain onto ®nite intervals is necessary before spectral expansions can be used. Further details are discussed in Ref. [3]. In cases of complex geometry or when di€erent resolutions and/or di€erent numerical methods shall be used within the ¯ow ®eld, spectral domain decomposition methods can be applied

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

559

with subdomains that are either overlapping or nonoverlapping. Due to high computational costs of these methods there are only few DNS applications at present. Hartel et al. [23] used the spectral domain decomposition method in one ¯ow direction only to study intrusion fronts. Concerning general numerical aspects, see e.g. Refs. [14,35]. Besides spectral methods, ®nite di€erence methods based on second- or fourth-order accurate central di€erence formulae have been traditionally used for DNS of turbulent ¯ows. The suitability of a second-order central di€erence scheme for DNS is discussed in Ref. [5] based on DNS of channel ¯ow. The authors come to the conclusion that ®nite-di€erence schemes (in general) can be more easily applied to complex geometries than spectral methods. Moreover, for a second-order central scheme the same number of grid points as for the spectral scheme is sucient to have a `spectral' resolution of the velocity ¯uctuations. However, about twice the number of grid points in each coordinate direction is needed to achieve a `spectral' resolution of the vorticity ¯uctuations. More recently, high-order compact ®nite-di€erence schemes (in which values of the derivatives are considered in complement to the values of the function) have received more attention, see Refs. [10,43], and are now being applied in DNS of turbulence. 2.2. Time discretization Direct simulations mostly use explicit time advancement schemes or semi-implicit methods to overcome severe time-step restrictions due to locally ®ne grids. Such ®ne grids generally arise when the steep wall-normal gradients are resolved or near coordinate singularities like the axis of a cylindrical coordinate system or near a sharp corner. Linear stability criteria lead to time steps for explicit schemes which are much below the Kolmogorov time scale tk ˆ …m=e†1=2 . m (m2 /s) is the kinematic viscosity and e (m2 /s3 ) the dissipation rate of turbulent kinetic energy k. Another time scale which is of the order of tk in near-wall turbulence, is the viscous time scale m=u2s (us is the friction velocity). It should also provide an upper bound for the time step permissible to guarantee physically correct solutions. Now, fully implicit time advancement schemes as e.g. used by Choi et al. [5] allow for time steps bigger than m=u2s from the point of view of numerical stability. The question whether time steps Dt such as e.g. 2m=u2s lead to an accurate prediction of turbulence statistics has been answered by Choi and Moin [4] for the minimal ¯ow unit of Jimenez and Moin [34], with the result that Dt should be close to 0:2m=u2s to sustain turbulence. It is quite likely that this estimate is too conservative for DNS of turbulent ¯ow in a full channel. Own computations [49] showed that the minimal ¯ow unit is a very delicate ¯ow which changes from the laminar to the turbulent state sporadically on each channel side and therefore needs higher spatial and temporal resolution than fully developed turbulent ¯ow. Explicit time advancement schemes are mostly based on second-, third- or fourth-order accurate Runge±Kutta methods. They exhibit a larger stability domain and thus allow for larger time steps than the more traditional second-order accurate Adams±Bashforth and Leapfrog schemes. In some ¯ow cases the severe time step restriction of explicit schemes can be somewhat relaxed by applying a Galilean coordinate transformation. For pipe ¯ow Akselvoll and Moin [1] have shown that domain decomposition into a core and a wall region and the implicit treatment of convection and di€usion terms with derivates in only one coordinate direction (circumferential in the core region and radial in the wall region) allows for time steps which are by a factor of 2.5 larger than those of Eggels et al. [12].

560

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

2.3. Arrangement of variables and solution of the Poisson equation The non-dimensional unsteady incompressible Navier±Stokes equations without external force have the form:   o~ u 1 2 ‡N ~ u ‡ rp r~ u ˆ 0; …1† ot Re r ~ u ˆ 0:

…2†

N …~ u† represents the nonlinear convection term. Harlow and Welsh [24] have introduced the concept of staggered meshes to integrate Eqs. (1) and (2) numerically. It stores the pressure in the cell centre and each velocity component in the centre of that cell face which is perpendicular to the direction of the velocity component. The great advantage of this concept is that it avoids oscillations in the computed pressure. It is therefore applied even in computations on generalized coordinates [5,77,78]. Collocated variable arrangements lead to inherently unstable discretizations [58]. The most common solution methods of Eqs. (1) and (2) deal with a Poisson equation for the pressure (or the pressure correction) and with the momentum equations for the computation of velocity. Harlow and Welsh's marker-and-cell (MAC) method is the prototype of such a method and equivalent to explicit versions of the projection or fractional step methods independently proposed by Chorin [8] and Temam [67]. A fully implicit fractional step method is used in Ref. [5]. We shortly discuss an explicit projection method of second-order accuracy in time. In a ®rst step a provisional velocity ®eld ~ u is ®rst computed from ~ u

~ un 2 Dt

1

ˆ

N …~ un † ‡

1 2 n 1 r~ u : Re

…3†

The second step corrects the velocity ®eld by considering the equations ~ u un‡1 ~ ˆ 2 Dt r ~ un‡1 ˆ 0:

rpn‡1 ;

…4† …5†

Taking the divergence of Eq. (4) and making use of Eq. (5) results in a Poisson equation for the pressure: 1 r ~ u : …6† 2 Dt In Eq. (3) a leapfrog step has been used for time advancement of the convection term. Boundary conditions for p are obtained by projecting the correction equation (4) on the outward normal unit vector ~ n to a speci®c boundary C of the ¯ow ®eld. This provides the Neumann condition  n‡1  op 1  n‡1 ~ ~ ˆ uC uC  ~ n; …7† on C 2 Dt r2 pn‡1 ˆ

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

561

where ~ uC is the boundary value of the provisional velocity ®eld. While solving the Eqs. (4)±(7) the following compatibility condition has to be satis®ed: Z Z   ~ un‡1 ~ u  ~ n ds: …8† r ~ u dX ˆ X

C

It will be identically satis®ed at each time step, if the condition of zero mass ¯ux across all boundaries of the computational domain X, viz Z ~ n ds ˆ 0 …9† uC  ~ C

is ful®lled. It is important that space discretizations satisfy Eq. (8). Returning to the question of boundary conditions for ~ uC , Peyret and Taylor [58] show that the use of staggered grids and centred ®nite di€erences leads to solutions independent of the value of ~ uC . Therefore,  ~ un‡1 ~ n ˆ 0, which corresponds to a zero-derivative condition for Eq. (7), is feasible. Yet, it uC ~ C does not imply zero physical pressure gradient. The ecient solution of the Poisson equation (6) at each time step is of paramount importance. It depends crucially on the problem to be solved. If the turbulence ®eld is homogeneous in one or even more directions and can be discretized in these directions by cartesian coordinates, fast Fourier transform of the discretized Poisson equation leads to sets of decoupled Helmholtz problems. 2D algebraic equations of this type can be eciently solved by cyclic reduction techniques. A set of discrete 1D Helmholtz equations can be treated by a tridiagonal-matrix algorithm. Even with fast elliptic solvers the computational time spent in solving the Poisson equation is of the order of 50±60% of the total computational time of a DNS. This amount can go up to more than 90% if the ¯ow ®eld is inhomogeneous in all three directions or if generalized coordinates are used. In such cases ecient elliptic solvers are needed. At present multigrid methods are used [5,20,76]. They are iterative solvers, like the conjugate gradient (CG) methods, see Ref. [19]. While CG methods are applicable only to symmetric coecient matrices of the discrete Poisson problem (6), the CGSTAB (CG squared stabilized) methods and GMRES (generalized minimal residual)-type methods deal with nonsymmetric systems. CGSTAB algorithms have been proposed by Van den Vorst and Sonneveld [69] and Van den Vorst [68]. The ®rst GMRES method has been suggested by Saad and Schultz [62]. For further details we refer to Refs. [9,13]. 2.4. Initial and boundary conditions With respect to economical use of computational time, an ideal situation is, if the computation can be started with initial conditions that are close to the ®nal turbulent state. Only rarely this happens to be so (e.g. when a ¯ow parameter is only slightly changed). In most cases, such initial conditions are not known. They have to be generated in some way. Fortunately, turbulent ¯ows get independent of their initial conditions if sucient time has elapsed. Thus, initial velocity ¯uctuations can be quite arbitrary and in general need not match speci®c spectra. It is sucient to generate a spatially uncorrelated divergence-free ¯uctuating velocity ®eld using a random number generator which provides numbers between 1 and 1, to weight them according to expected turbulence intensities, and to allocate them to the mesh cells. These random ¯uctuations are

562

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

superimposed on computed or measured mean velocity ®elds. It may then happen that the ¯uctuations decay in the course of the computation since they have been introduced at the highest wave number causing unphysically high viscous forces. This can be overcome by arti®cially increasing the Reynolds number and then gradually decreasing it to the desired value within one or two problem times. There is then sucient transient time in which ¯uctuating energy can be distributed among all wave numbers. Five types of boundary conditions are required for DNS: homogeneous, wall, in¯ow, entrainment and out¯ow boundary conditions. Flows which are statistically homogeneous in one or more directions can be handled using periodic boundary conditions in these directions, however, care must be taken to ensure that the computational domain is suciently large to contain a representative number of large-scale structures or more speci®cally, at least twice as large as the separation at which two-point correlations vanish. Boundary conditions at smooth solid walls are straightforward. The representation of permeable or rough walls needs some modelling. In¯ow/ out¯ow boundary conditions are needed when the ¯ow is inhomogeneous in the main ¯ow direction. Their speci®cation forms one of the major problems of DNS because the ¯ow variables at these boundaries depend on the unknown ¯ow outside the domain. Conditions have to be designed in such a way that the propagation of boundary errors is minimized. Le and Moin [40] describe a method of generating in¯ow ¯uctuations for turbulent ¯ow over a backward-facing step. It provides a set of stochastic signals that satisfy prescribed second-order statistics at the inlet. A severe drawback of this method, however, lies in the fact that the ¯ow quickly loses its statistical characteristics within a few integral lengths from the inlet and recovers the proper characteristics only slowly downstream. This is because turbulent structures are missing at the inlet. A substantial ``transition'' length is thus unavoidable. It can be made very short, if the unsteady in¯ow conditions are taken from a DNS which is run parallel to the DNS of the real problem. This technique has been used in pipe expansion ¯ow by Wagner and Friedrich [72]. Recently, Lund et al. [47] have developed a method for boundary layers where the velocity ®eld at a downstream station is re-scaled and re-introduced at the inlet. A similar strategy was used earlier in LES of turbulent boundary layers by Richter et al. [60]. Convective out¯ow boundary conditions have been proposed independently by Lowery and Reynolds [46] for mixing layers and by Richter et al. [60]. Lele et al. [43] use such a condition in the form o~ u o~ u ‡ Uc ˆ0 ot @x

…10†

and conduct tests to determine the proper convection velocity Uc in the main ¯ow direction x. Finally, they keep Uc constant, equal to the mean exit velocity. Wagner and Friedrich [72] use a convection equation to specify the velocity ¯uctuations ~ u0 in the exit plane of a sudden pipe expansion, viz: o~ u0 o~ u0 ‡ hui ˆ0 ot ox

…11†

and extrapolate the mean velocity linearly. The mean axial velocity hui varies with the wallnormal distance. The idea behind condition (11) is that the ¯uctuating ®elds are frozen while they are convected through the exit plane.

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

563

Entrainment boundary conditions, e.g. far away from the edges of a boundary layer, jet or mixing layer require careful formulation to prevent the appearance of spurious disturbances. Mappings to in®nity can help. In any case it is necessary to ensure that the potential-¯ow behaviour of the freestream velocity ®eld is well represented [65]. For DNS of a separated boundary layer Na and Moin [53] have devised a physically realistic entrainment boundary condition in which the wall-normal velocity component is prescribed while the streamwise component adjusts itself according to a zero-vorticity condition. 2.5. Resolution requirements DNS has to resolve all turbulent scales. As an example we discuss the spatial resolution required in a DNS of fully developed ¯ow through smooth pipes. Both the largest and smallest turbulence structures are dynamically signi®cant and must be accurately predicted. The largest structures are of the order of the pipe radius. The size of the smallest is controlled by 1=4 viscosity and the dissipation rate e, and is O…lk †, where lk ˆ …m3 =e† is the Kolmogorov length scale.  peaks at the wall and reaches values of the order of 0:2u4s =m. Of course, this value is a result of computation and thus unknown a priori. It is therefore more convenient to use the viscous length scale m=us , instead of lk , as the smallest turbulent scale to be resolved. In order to express it in terms of global parameters, we use the Blasius formula for the resistance coef®cient k: k ˆ 0:3164Reb

1=4

:

…12†

Reb is the Reynolds number based on the bulk velocity ub and the pipe diameter D. For smooth pipes, k is: 2

2

k ˆ 8…us =ub † ˆ 8…Res =Reb † :

…13†

If we use a cylindrical coordinate system and resolve the pipe radius by an equidistant mesh of width 2m=us (in reality one prefers nonequidistant meshes), we need (taking Eqs. (12) and (13) into account), the following number of mesh points in r-direction: D=2 1 ˆ Res ˆ 0:049718Re7=8 b : 2m=us 4

…14†

Assuming a mesh size of 8m=us along the pipe circumference, which means that streaky structures are resolved by roughly 12 points, and the same mesh size in axial direction to discretize a pipe of length 5D, the total number of mesh points is: Nˆ

5p 3 21=8 Re ˆ 0:0004826Reb : 256 s

…15†

For Reynolds numbers Res of the order of 300±400, the length 5D suces to guarantee the complete decay of two-point velocity correlations. For Res ˆ 360, e.g. the total number of mesh points is 2:86  106 . Formula (15) can also be used to estimate the Reynolds number limitations for DNS on a speci®c computer.

564

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

Fig. 1. Rms velocity ¯uctuations in fully developed channel ¯ow at us h=m ˆ 180 for a DNS with a zonal embedded grid. ( ): ®ne grid; ( ): global grid.

2.6. Zonal embedded grids Close to a solid wall the energy-carrying turbulence structures are small compared to those far from the wall which needs grid re®nement. A way to do it is with zonal embedded grids. One designs several zones of grids in which the mesh spacing decreases from zone to zone as one approaches the wall. Thus, the total number of grid points and the computational costs can be reduced. Kravchenko et al. [39] have developed an ecient method in which a Galerkin method with B-spline basis functions is used in the inhomogeneous ¯ow direction, the direction with the zonal boundaries. Fourier spectral methods are applied in the remaining homogeneous directions. An advantage of a B-spline based method over ®nite-volume methods is high accuracy. The technique has been tested for fully developed channel ¯ow, both in DNS and LES calculations, and excellent agreement with classical DNS and with experimental data has been obtained. Manhart [48] has used zonal embedded grids in the context of a ®nite volume technique on staggered grids which is of second-order accuracy in space and time. In computations of fully developed channel ¯ow, the grid is re®ned in the near-wall layer (®ne grid zone) by a factor of two in all three coordinate directions compared to the (coarser) global grid covering the whole domain. The ®ne grid overlaps the coarse grid. At the interface where the ®ne grid ends and the coarse grid continues, intergrid boundary conditions have to be speci®ed by interpolation such that mass, momentum and energy are conserved cellwise. This is hard to achieve. It turns out that ®rst-order interpolation for the ®ne-grid velocity boundary condition and Neumann boundary conditions for the ®ne-grid pressure correction work satisfactorily. Fig. 1 shows a comparison of rms velocity ¯uctuations computed on a ®ne, a coarse and on zonal grids. The grid interface is at z‡ ˆ 55. It can be seen that the algorithm works quite well.

3. DNS results of ¯ow through pipes A few examples will be given below in which DNS is used to study the physics of turbulent ¯ows. They cover the e€ect of wall permeability or wall roughness in straight pipes, the separation

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

565

and re-attachment of turbulent ¯ow in a sudden pipe expansion and the e€ects of curvature and torsion on the turbulence structure in coiled pipes. 3.1. Fully developed ¯ow in solid and permeable pipes A major e€ect characterizing the interaction of a solid wall with a turbulent ¯ow is the so-called wall blocking e€ect. It limits the maximum turbulent kinetic energy production in the wall layer to levels much lower than in mixing layers with comparable main shear, where the increased production rate is primarily due to the enhancement of velocity ¯uctuations in the direction of main shear. This was the motivation to study turbulent ¯ow through permeable pipes where these ¯uctuations are nonzero at the wall. Six di€erent direct simulations of fully developed ¯ow through solid and permeable pipes were performed with a ®nite volume code of second-order accuracy. The pipe was ®ve diameters D long and was resolved by 486  240 equidistant grid points in … z; u†-directions. In radial …r† direction 70 grid points were used and clustered close to the wall. In the solid wall case …DNS1 † with a Reynolds number Res ˆ us D=m ˆ 360, the viscous sublayer contained six points. In terms of wall units the grid spacing was: Dz‡ ˆ 3:7;

…r‡ Du† ˆ 0:073  4:71;

Dr‡ ˆ 0:185  5:56:

…16†

Table 1 summarizes the ¯ow and computational parameters of all simulations. While the no-slip boundary conditions were applied at the wall in all cases, the impermeability condition, …ur †rˆR was relaxed at permeable walls and replaced by: …ur r†rˆR ˆ a…ur r†rˆR

Dr :

…17†

For a ˆ 1 this is the discrete counterpart of o…ur r†=or ˆ 0 at the wall. The damping factor a determines the level of permeability. a ˆ 0 represents a solid wall and a ˆ 1 a perfectly permeable wall. In DNS6 each second mesh surface along the wall was assumed perfectly permeable (a ˆ 1), while the remaining ones were kept impermeable. Reb is the Reynolds number based on the bulk velocity ub and k ˆ 8…us =ub †2 , the friction factor. The same mean pressure gradient was used in DNS1;2;5;6 to drive the ¯ow. In DNS3;4 it was increased to achieve higher Reynolds number and ®nally to approach the mean mass ¯ux of DNS1 . It turned out that partially permeable walls interact with the turbulent ¯ow in a way very similar to a rough wall. In fact, the mean velocity pro®les coincide with Nikuradse's law for the rough pipe [63]: Table 1 Parameters of fully developed turbulent ¯ow in pipes with perfectly permeable walls (DNS2;3;4 ), partially permeable walls (DNS5;6 ) and with a solid wall (DNS1 ) DNS Res a Reb ub k  102 kR

1

2

3

4

5

6

360 0.0 5317 14.77 3.67 ±

360 1.0 2077 5.77 24.0 0.70

720 1.0 4327 6.01 22.2 0.65

900 1.0 5499 6.11 21.4 0.65

360 0.975 4406 12.24 5.30 0.055

360 ± 4665 12.96 4.80 0.043

566

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

huz i=us ˆ 2:87 ln…y=kR † ‡ 8:48

…18†

with a slope modi®ed from 2.5 to 2.87. The roughness heights kR in Table 1 have been obtained by comparing DNS data with Eq. (18). Permeable or rough walls change the turbulence structure in a pipe dramatically. The viscous sublayer disappears. DNS6 , e.g. represents the case of a roughness height kR ˆ 0:043D or kR‡ ˆ kR us =m  15 which is three times the viscous sublayer thickness. The Reynolds shear stress grows linearly with the distance from the (permeable) wall and strongly reduces the viscous stress, i.e. the mean velocity and in turn the mass ¯ux. As an example of the changes in the turbulence structure due to wall permeability or roughness, we present pro®les of the mean velocity, the Reynolds shear stress, the ¯atness factor of the radial rms velocity ¯uctuation and the rms vorticity ¯uctuations. Fig. 2 shows mean velocity pro®les for the six DNS cases. DNS1 is con®rmed by PIV measurements. Partially permeable pipes generate pro®les that are typical for ¯ow through sandrough pipes. Therefore, wall permeability mimics wall roughness. The pro®le for the highest Res seems close to that for Re-number independence. Fig. 3 illustrates the steep rise of the Reynolds shear stress close to the wall (which is in total contrast to the solid wall case). The ®gure also indicates the limit of Res ! 1. The ¯atness factor of the radial velocity ¯uctuation in Fig. 4 re¯ects

Fig. 2. Mean velocity pro®les versus wall distance in wall units. The squares represent LDA measurements of Westerweel et al. [79], (Ð): DNS1 , (- - -): DNS2 , (--): DNS3 , (--): DNS4 , (±  ±  ±): DNS5 , (  ): DNS6 .

Fig. 3. Reynolds shear stress pro®les (lines as in Fig. 2).

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

567

Fig. 4. Flatness factors of the radial velocity ¯uctuations (lines as in Fig. 2).

Fig. 5. Rms vorticity ¯uctuations xz;rms (lines as in Fig. 2).

Fig. 6. Rms vorticity ¯uctuations xu;rms (lines as in Fig. 2).

strongly intermittent (sweep) events close to the solid wall. On the other hand, the velocity ¯uctuations do not `see' perfectly permeable walls which results in ¯atness factors F …u0r † ˆ 3. Rms vorticity ¯uctuations in Figs. 5±7 indicate further structural changes. Furthermore the spacing of streaky structures is modi®ed. Perfectly permeable walls nearly double the spacing while partially permeable or rough walls (DNS5;6 ) slightly decrease it. For more details the reader is referred to Refs. [74,75].

568

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

Fig. 7. Rms vorticity ¯uctuations xr;rms (lines as in Fig. 2).

3.2. Flow in pipes with sudden expansion Flow separating at the corner of a sudden pipe expansion and re-attaching further downstream is the prototype of a complex turbulent ¯ow. Its reliable prediction with statistical methods is still a challenge today. The free shear layer which develops downstream of the corner has its own dynamics, di€erent from the upstream near-wall turbulence dynamics. The maximum in the variance of the radial velocity ¯uctuations in the shear layer exceeds that in the wall layer by a factor of 6. There is also a 4-fold increase in the Reynolds shear stress leading to a corresponding strong enhancement of turbulence production while at the same time the mean rates of strain are only marginally modi®ed with respect to the wall layer upstream. Motivated by the idea that the rapid growth of turbulent kinetic energy in the free shear layer might be mainly due to the fact that the velocity ¯uctuations in the direction of main shear (the u0r ¯uctuations) relax from the impermeability constraint at the wall, we have conducted two direct simulations in sudden pipe expansions with solid (case S) and partially permeable (a ˆ 0:975) or rough walls (case R). The computational domain consists of an upstream pipe section of length 1D, diameter D and a downstream section of diameter De ˆ 1:2D and length 2D. The expansion ratio is thus De =D ˆ 1:2. The incoming ¯ow is fully developed at a Reynolds number of Res ˆ 360 in both cases, and is taken from two separate simulations of ¯ow in solid and permeable pipes with slightly coarser grids than those for DNS1 and DNS5 , namely 256  128  70 cells in …z; u; r†directions. The grid that discretizes the pipe expansions consists of 154  128  104 points. It is nonequidistant only in r-direction. Only two results will be discussed here which seem to underline the importance of the wall blocking e€ect or the relaxation of the free shear layer from this e€ect. Figs. 8 and 9 contrast the

2

Fig. 8. Contour lines of hu0r i for case S, normalized with the surface averaged total kinetic energy of the in¯ow plane h12ui ui iav . The maximum in the shear layer is 3:669  10 2 , compared to 6:04  10 3 upstream.

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

569

2

Fig. 9. Contour lines of hu0r i for case R, normalized with the surface averaged total kinetic energy of the in¯ow plane h12ui ui iav . The maximum in the shear layer is 3:427  10 2 , compared to 8:54  10 3 upstream.

spatial development of the radial velocity variance hu0r 2i in solid and rough pipes. The values of the contour lines have been nondimensionalized with the mean kinetic energy in the in¯ow plane as obtained from surface and time averaging the magnitude of the instantaneous velocity vector there. Thus the spatial development is independent of the energy content of the system. The ratio 2 2 of the maximum in the shear layer (SL) to that of the incoming ¯ow (UP), hu0r imax;SL =hu0r imax;UP is 6.07 for solid walls and 4.01 for rough walls. Moreover, the position of the maximum is closer to the step in case R. This is consistent with the reduction in the re-attachment length zR =h from 10.1 in case S to 8.8 in case R. Similar e€ects are observed in the Reynolds shear stress. The corresponding ratios of maximal values, hu0z u0r imax;SL =hu0z u0r imax;UP are 7.15 in case S and 4.92 in case R. This has consequences not only for the production terms but also for the pressure±strain, turbulent di€usion and dissipation terms. Turbulence models must be capable of capturing these e€ects. More details of turbulent ¯ow in solid wall pipe expansions are found in Refs. [16,71±73]. 3.3. Flow in coiled pipes Turbulent ¯ow in coiled pipes appears in many industrial applications and is hard to predict with statistical methods. The centreline of a coiled pipe has not only curvature but also torsion and forms a helix of radius ra and pitch ps . Curvature alone leads to secondary motions due to centrifugal e€ects. The parameter characterizing these e€ects is the Dean number De. Adding torsion to a curved pipe leads to a coiled (or helically coiled) pipe. Torsion alone has no ¯uid dynamic e€ect as long as the walls are smooth. The Germano number Gn speci®es the e€ect of torsion on the ¯ow structure in coiled pipes. The third important parameter is the Reynolds number. We call hpiav the mean pressure averaged over a circular cross-section R2 p of the coiled pipe and Ds the distance measured along the centreline between two cross-sections. Then the pressure di€erence which drives the ¯ow, balances the wall shear stress, sw;m , averaged over the pipe's circumference: R hDpiav : 2 Ds From Eq. (19) the corresponding averaged friction velocity us results (see Refs. [25±30]): sw;m ˆ

…19†

us ˆ …sw;m =q†1=2 :

…20†

The ¯ow parameters speci®ed in direct simulations can now be de®ned as: Res ˆ

Rus ; m

Des ˆ

p j Res ;

Gns ˆ s Res ;

…21†

570

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

where j, s denote the `curvature' and `torsion' of the pipe axis: jˆ

ra2

ra ; ‡ ps2



ra2

ps : ‡ ps2

…22†

We use a helical, orthogonal …s; r; h†-coordinate system proposed by Germano [18] in our DNS of fully developed turbulent ¯ow in coiled pipes. It relies essentially on the Frenet triad of unit vectors associated with the pipe centreline and provides a local ¯ow description which is valid for R < 1=j. s is the arclength measured along the centreline and …r; h† are polar coordinates in the circular cross-section. h is measured clockwise. h ˆ 0 de®nes a vertical line from the centre to the upper side of the pipe wall, and h ˆ p=2 indicates the outer side of the pipe where the ¯uid experiences concave wall curvature in main ¯ow direction. The horizontal axis h ˆ p=2 is a symmetry axis for mean ¯ow through curved pipes. There is no such symmetry axis for ¯ow a€ected by torsion. The computational domains are pipes of 15:23R in length and diameter D ˆ 2R. They are resolved by 256  70  180 cells in …s; r; h†-directions. The grid spacing is nonequidistant only in r-direction. Three ¯ow cases are discussed: ¯ow through a straight pipe (DP), a curved pipe (DT) and a helically coiled pipe (DXXH). The geometrical and ¯ow parameters are listed in Table 2. Cases DP and DT are also discussed by H uttl and Friedrich [26], further investigations on these and other cases can be found in Refs. [15,25±32]. Case DXXH has the highest torsion of all cases investigated. The fact that Res is kept constant means that the same pressure di€erence is used to drive these ¯ows. It should be noted that parameters in Table 2 involving the bulk velocity ub are results of the computation and not a priori known. k ˆ 8…us =ub †2 is the friction factor. A ®rst observation to be made is that the secondary ¯ow induced by curvature reduces the volume rate (hence Reb ) and increases the friction factor. Computational parameters are usually discussed in wall units. The pipe length is 3503 wall units, i.e roughly three times the length of a streaky structure. The length along the inner side is 90% of that and therefore still suciently large. The axial grid spacing Ds‡ is 15 along the outer side and 12.3 along the inner side (case DT). The radial grid spacing varies between 1.03 and 5.38 (in the core region). Close to the wall the circumferential grid spacing is 8.01. If one takes the variation of the friction velocity along the pipe circumference into account (see Fig. 10), Ds‡ , based on local values, would have values of 18 along the outer side of the curved pipe and 6.6 along the inner side. It will, therefore, be interesting to see the e€ect of axial grid re®nement in succeeding computations. Statistical results presented below have been obtained by time averaging and spatial averaging in streamwise direction over those points which have the same …r; h†position as the point of the in¯ow plane (the s-coordinate in case DXXH does not connect these

Table 2 Geometrical and ¯ow parameters for straight (DP), curved (DT) and coiled pipe (DXXH) Case

jR

sR

ra =R

ps =R

Res

Reb

Des

Gns

ub =us

k  102

DP DT DXXH

0.0 0.1 0.1

0.0 0.0 0.165

1 10.0 2.68

0.0 0.0 4.43

230 230 230

6812 5624 5576

0.0 72.7 72.7

0.0 0.0 37.9

14.80 12.23 12.12

3.64 5.34 5.44

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

571

Fig. 10. Pro®les of the friction velocity us : (Ð) DT, (- - -) DXXH.

Fig. 11. Pro®les of the mean wall friction in axial direction sw;s : (Ð) DT, (- - -) DXXH.

Fig. 12. Pro®les of the mean wall friction in circumferential direction sw;h : (Ð) DT, (- - -) DXXH.

points). On the average 200 statistically independent time samples of the ¯ow ®eld have been used to obtain stable statistics. Figs. 11 and 12 show the axial and circumferential components of the wall shear stress normalized with us (from Eq. (20)) and plotted as functions of h. The highest axial shear stress is on the outer part of the bend. The secondary ¯ow leads to stress components sw;h which reach more than 60% of the mean axial wall stress sw;m . The e€ect of strong torsion increases the amplitude of

572

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

sw;h near h ˆ p by nearly 15% and decreases it even more near h ˆ 0. Fig. 12 already indicates that the main reason for the secondary ¯ow is the centrifugal e€ect due to pipe curvature, but strong torsion modulates this e€ect noticeably. We recall that the elevation of the pipe after one revolution is about 10 times the radius ra of the helix. The mean axial velocity hus i in Fig. 13 is only marginally modi®ed by torsion. The centrifugal e€ect forces the rapidly ¯owing ¯uid particles to the outside of the curved and coiled pipes with the consequence that hus i is low on the inner side and high on the outer side. The vertical pro®les have maxima near the wall due to the secondary motion. Note that hus i (like all other velocity components, even ¯uctuations) has been nondimensionalized with us . The mean circumferential and radial velocity components in Figs. 14 and 15 have remarkable amplitudes of the order of the rms velocity ¯uctuations. The e€ect of torsion on these components is pronounced. While turbulent ¯ow through straight pipes is characterized by four axisymmetric components of the Reynolds stress tensor, there are six non-zero components in the case of curvature and torsion. Their amplitudes reach values of the order of u2s locally, as demonstrated in Figs. 16±21. The strong spatial variations of these quantities are due to the complex structure of the mean ¯ow, as seen in the contour plots (Figs. 8 and 9) presented by H uttl and Friedrich [26] for case DT. The two vortices observed on the inner side of the curved pipe are regular Dean vortices (see the mean

Fig. 13. Pro®les of the mean axial velocity component hus i: (Ð) DP, (- - -) DT (h ˆ 0), (  ) DT (h ˆ p=2), (--) DXXH (h ˆ 0), (- - - -) DXXH (h ˆ p=2).

Fig. 14. Pro®les of the mean circumferential velocity component huh i: (Ð) DP, (- - -) DT (h ˆ 0), (  ) DT (h ˆ p=2), (±  ±  ±) DXXH (h ˆ 0), (- - - -) DXXH (h ˆ p=2).

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

573

Fig. 15. Pro®les of the mean radial velocity component hur i: (Ð) DP, (- - -) DT (h ˆ 0), (  ) DT (h ˆ p=2), (±  ±  ±) DXXH (h ˆ 0), (- - - -) DXXH (h ˆ p=2).

Fig. 16. Pro®les of the Reynolds stress hu00s u00s i: (Ð) DP, (- - -) DT (h ˆ 0), (  ) DT (h ˆ p=2), (±  ±  ±) DXXH (h ˆ 0), (- - - -) DXXH (h ˆ p=2).

Fig. 17. Pro®les of the Reynolds stress hu00h u00h i: (Ð) DP, (- - -) DT (h ˆ 0), (  ) DT (h ˆ p=2), (±  ±  ±) DXXH (h ˆ 0), (- - - -) DXXH (h ˆ p=2).

574

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

Fig. 18. Pro®les of the Reynolds stress hu00r u00r i: (Ð) DP, (- - -) DT (h ˆ 0), (  ) DT (h ˆ p=2), (±  ±  ±) DXXH (h ˆ 0), (- - - -) DXXH (h ˆ p=2).

Fig. 19. Pro®les of the Reynolds shear stress hu00s u00h i: (Ð) DP, (- - -) DT (h ˆ 0), (  ) DT (h ˆ p=2), (±  ±  ±) DXXH (h ˆ 0), (- - - -) DXXH (h ˆ p=2).

Fig. 20. Pro®les of the Reynolds shear stress hu00s u00r i: (Ð) DP, (- - -) DT (h ˆ 0), (  ) DT (h ˆ p=2), (±  ±  ±) DXXH (h ˆ 0), (- - - -) DXXH (h ˆ p=2).

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

575

Fig. 21. Pro®les of the Reynolds shear stress hu00h u00r i: (Ð) DP, (- - -) DT (h ˆ 0), (  ) DT (h ˆ p=2), (±  ±  ±) DXXH (h ˆ 0), (- - - -) DXXH (h ˆ p=2).

streamfunction W or the vector plot). These vortices are a€ected by torsion. They are slightly turned clockwise and the lower one is appreciably damped. It is also important to note that the overall level of the turbulence ¯uctuations is reduced due to curvature, especially in the near-wall regions. Torsion on the other hand enhances the turbulence activity in the core region compared to pure curvature. These results show that reliable statistical predictions of ¯ow through curved or coiled pipes are a great challenge to turbulence modellers since already the production of turbulent kinetic energy depends on all components of the Reynolds stress tensor. Two-equation turbulence models must, therefore, fail.

4. Conclusions DNS of incompressible turbulent ¯ows with second-order central schemes is shown to give reliable results, provided the smallest turbulent scales are properly resolved. The viscous length scale is a suitable measure of the smallest scales. Although higher-order schemes are desirable, the speci®cation of proper boundary conditions, compatible with the order of accuracy of these schemes is a delicate task which has to be treated with care. Second-order central schemes have been used more extensively in the past than higher-order schemes to predict ¯ows with complex geometry. Based on these schemes, DNS results are presented and interpreted physically for ¯ow through solid and permeable pipes, through sudden pipe expansions and through coiled pipes. It is found that permeable pipes act like rough pipes. Roughness heights of the order of the bu€er layer thickness modify the near-wall ¯ow mainly, whereas roughness heights of the order of the pipe radius change the complete ¯ow. The rapid growth of components of the Reynolds stress tensor observed within the free shear layer leaving the corner of a sudden pipe expansion, is considerably diminished when the pipe walls are permeable (or rough). At the same time, the re-attachment length of the ¯ow is reduced. These ®ndings point towards considerable structural changes of attached and separated turbulent ¯ows along permeable walls. Finally, the simulated fully developed ¯ows in coiled pipes are dominated by centrifugal e€ects which dampen the turbulence activity, but lead to six non-zero Reynolds stress components which have to be modelled in

576

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

statistical simulations. Computations of ¯ows in coiled pipes at Reynolds numbers considerably higher than those used here are urgently needed in order to guarantee fully turbulent ¯ow. References [1] Akselvoll K, Moin P. An ecient method for temporal integration of the Navier±Stokes equations in con®ned axisymmetric geometries. J Computat Phys 1996;125:454±63. [2] Antonia RA, Kim J. Low-Reynolds-number e€ects on near-wall turbulence. J Fluid Mech 1994;276:61±80. [3] Canuto C, Hussaini MY, Quarteroni A, Zang TA. Spectral methods in ¯uid mechanics. New York: Springer; 1988. [4] Choi H, Moin P. E€ects of the computational time step on numerical solutions of turbulent ¯ow. J Computat Phys 1994;113:1±4. [5] Choi H, Moin P, Kim J. Turbulent drag reduction: studies of feedback control and ¯ow over riblets. Rep. TF-55, Department of Mechanical Engineering, Stanford University, Stanford, CA, 1992. [6] Choi H, Moin P, Kim J. Direct numerical simulation of turbulent ¯ow over riblets. J Fluid Mech 1993;255:503±39. [7] Choi H, Moin P, Kim J. Active turbulence control for drag reduction in wall-bounded ¯ows. J Fluid Mech 1994;262:75±110. [8] Chorin AJ. Numerical solution of the Navier±Stokes equations. Math Computat 1968;22:745±62. [9] Deng GB, Piquet J, Quentey P, Visonneau M. Navier±Stokes equations for incompressible ¯ows: ®nite di€erence and ®nite volume methods. In: Peyret R, editor. Handbook of computational ¯uid mechanics. London: Academic Press; 1996. p. 25±97. [10] Dervieux A. About the basic numerical methods. In: Peyret R, editor. Handbook of computational ¯uid mechanics. London: Academic Press; 1996. p. 1±23. [11] Eckelmann H. Einf uhrung in die Str omungsmechanik. Stuttgart: Teubner; 1997. p. 270±2. [12] Eggels JGM, Unger F, Weiss MH, Westerweel J, Adrian RJ, Friedrich R, Nieuwstadt FTM. Fully developed turbulent pipe ¯ow: a comparison between direct numerical simulation and experiment. J Fluid Mech 1994;268:175±209. [13] Ferziger JH, Peric M. Computational methods for ¯uid dynamics. Berlin: Springer; 1996. [14] Fischer PF. An overlapping Schwarz method for spectral element solution of the incompressible Navier±Stokes equations. J Computat Phys 1997;133:84±101. [15] Friedrich R, Lechner R, Sesterhenn J, H uttl TJ. Direct numerical simulation of turbulent compressible and incompressible wall-bounded shear ¯ows. In: Knight D, Sakell L, editors. Recent advances in DNS and LES. Proceedings of the Second AFOSR International Conference on Direct Numerical Simulation and Large Eddy Simulation, Rutgers ± The State University of New Jersey, New Brunswick, USA. 7±9 June, 1999, Kluwer Academic Publishers; 1999. p. 13±26. [16] Friedrich R, Wagner C. On turbulent ¯ow in a sudden pipe expansion and its reverse transfer of subgrid-scale energy. Lecture notes in physics 453, 1995. p. 330±6. [17] Gavrilakis S. Numerical simulation of low-Reynolds-number turbulent ¯ow through a straight square duct. J Fluid Mech 1992;244:101±29. [18] Germano M. On the e€ect of torsion on a helical pipe ¯ow. J Fluid Mech 1982;125:1±8. [19] Golub GH, van Loan C. Matrix computations. Baltimore: Johns Hopkins University Press; 1990. [20] Hackbusch W. Multi-grid methods and applications. New York: Springer; 1985. [21] Hartel C. Turbulent ¯ows: direct numerical simulation and large-eddy simulation. In: Peyret R, editor. Handbook of computational ¯uid mechanics. London: Academic Press; 1996. p. 283±338. [22] Hartel C, Kleiser L. Analysis and modelling of subgrid-scale motions in near-wall turbulence. J Fluid Mech 1998;356:327±52. [23] Hartel C, Kleiser L, Michaud M, Stein CF. A direct numerical simulation approach to the study of intrusion fronts. J Eng Math 1997;32:103±20. [24] Harlow FH, Welsh JE. Numerical calculation of time-dependent viscous incompressible ¯ow with free surface. Phys Fluids 1965;8:2182±9.

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

577

[25] H uttl TJ, Friedrich R. Fully developed laminar ¯ow in curved or helically coiled pipes. DGLR-Jahrbuch, DGLRJT97-181, 1997. S. 1203±10. [26] H uttl TJ, Friedrich R. Direct numerical simulation of turbulent ¯ows in curved and helically coiled pipes. Proceedings of the Third Asian Computational Fluid Dynamics Conference, vol. 2. Bangalore, India, 7±11 December, 1998. p. 183±8. [27] H uttl TJ, Friedrich R. High performance computing of turbulent ¯ow in complex pipe geometries. In: Krause E, Jager W, editors. High Performance Computing in Science and Engineering '98. Transactions of the High Performance Computing Center Stuttgart (HLRS), 1998. Springer; 1999. p. 236±51. [28] H uttl TJ, Friedrich R. Fully developed turbulent ¯ow in conduits with circular cross section. In: Nitsche W, et al., editors. New results in numerical and experimental ¯uid mechanics II. Contributions to the 11th AG STAB/DGLR Symposium, Berlin; Notes Numer Fluid Mech 1999;72:242±9. [29] H uttl TJ, Friedrich R. In¯uence of curvature and torsion on turbulent ¯ow in helically coiled pipes. In: Rodi W, Laurence D, editors. Engineering turbulence modelling and experiments 4. Amsterdam: Elsevier; 1999. p. 247± 56. [30] H uttl TJ, Friedrich R. Turbulent ¯ow in coiled pipes. In: Proceedings of the Joint INI/ERCOFTAC Workshop on Direct and Large-Eddy Simulation, vol. 3, paper 3.6. Isaac Newton Institute, Cambridge. 12±14 May, 1999. [31] H uttl TJ, Wagner C, Friedrich R. Navier±Stokes solutions of laminar ¯ows based on orthogonal helical coordinates. In: Taylor C, Cross J, editors. Numerical methods in laminar and turbulent ¯ow, vol. 10. Swansea, UK: Pineridge Press; 1997. S. 191±202. [32] H uttl TJ, Wagner C, Friedrich R. Navier Stokes solutions of laminar ¯ows based on orthogonal helical coordinates. Int J Numer Meth Fluids 1999;29:749±63. [33] Huser A, Biringen S. Direct numerical simulation of turbulent ¯ow in a square duct. J Fluid Mech 1993;257:65±95. [34] Jimenez J, Moin P. The minimal ¯ow unit in near-wall turbulence. J Fluid Mech 1991;225:213±40. [35] Karniadakis GE, Orszag SA, Ronquist EM, Patera AT. In: Gunzberger MD, Nicolaides RA, editors. Incompressible computational ¯uid dynamics. Cambridge: Cambridge University Press; 1993. p. 203±66. [36] Kasagi N, Shikazono N. Contribution of direct numerical simulation to understanding and modelling turbulent transport. Proc R Soc London Ser A 1995;451:257±92. [37] Kim J, Moin P, Moser R. Turbulence statistics in fully developed channel ¯ow at low Reynolds number. J Fluid Mech 1987;177:133±66. [38] Kleiser L, Zang ThA. Numerical simulation of transition in wall-bounded shear ¯ows. Annu Rev Fluid Mech 1991;23:495±537. [39] Kravchenko AG, Moin P, Moser R. Zonal embedded grids for numerical simulations of wall bounded turbulent ¯ows. J Computat Phys 1996;127:412±23. [40] Le H, Moin P. Direct numerical simulation of turbulent ¯ow over a backward-facing step. Rep. TF-58, Thermosciences Division, Department of Mechanical Engineering, Stanford University, 1994. [41] Le H, Moin P, Kim J. Direct numerical simulation of turbulent ¯ow over a backward-facing step. J Fluid Mech 1997;330:349±74. [42] Lee MJ, Kim J, Moin P. Turbulence structure at high shear rate. Proceedings of the Sixth Symposium Turbulent Shear Flows, Toulouse, 1987. p. 22-6-1±6. [43] Lele SK. Compact ®nite di€erence schemes with spectral-like resolution. J Computat Phys 1992;103:16±42. [44] Lesieur M, Comte P, Metais O. Direct numerical simulations of turbulence. In: Hirsch Ch, editor. Computational methods in applied sciences. Amsterdam: Elsevier; 1992. p. 37±43. [45] Lim J, Choi H, Kim J. Control of streamwise vortices with uniform magnetic ¯uxes. Phys Fluids 1998;10:1997± 2005. [46] Lowery PS, Reynolds WC. Numerical simulation of a spatially-developing, forced, plane mixing layer. Rep. TF-26, Thermosciences Division, Department of Mechanical Engineering, Stanford University, 1986. [47] Lund TS, Wu X, Squires KD. Generation of turbulent in¯ow data for spatially-developing boundary layer simulations. J Computat Phys 1998;140:233±58. [48] Manhart M. Zonal direct numerical simulation of turbulent plane channel ¯ow. Notes on numerical ¯uid mechanics, vol. 64. Braunschweig: Vieweg; 1998.

578

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

[49] Manhart M, Deng GB, H uttl TJ, Tremblay F, Segal A, Friedrich R, Piquet J, Wesseling P. The minimal turbulent ¯ow unit as a test case for three di€erent computer codes. Notes on numerical ¯uid mechanics, vol. 66. Braunschweig: Vieweg; 1998. [50] Moin P, Mahesh K. Direct numerical simulation: a tool in turbulence research. Annu Rev Fluid Mech 1998;30:539±78. [51] Moser R, Moin P. The e€ects of curvature in wall-bounded turbulent ¯ows. J Fluid Mech 1987;175:479±510. [52] Moser RD, Rogers M. The three-dimensional evolution of a plane mixing layer: pairing and transition to turbulence. J Fluid Mech 1993;247:275±320. [53] Na Y, Moin P. Direct numerical simulation of a separated turbulent boundary layer. J Fluid Mech 1998;370:175± 201. [54] Neves JC, Moin P, Moser RD. E€ects of convex transverse curvature on wall-bounded turbulence. Part 1: The velocity and vorticity. J Fluid Mech 1994;272:349±81. [55] Neves JC, Moin P, Moser RD. E€ects of convex transverse curvature on wall-bounded turbulence. Part 2: The pressure ¯uctuations. J Fluid Mech 1994;272:383±406. [56] Perot B, Moin P. Shear-free turbulent boundary layers. Part 1: Physical insights into near-wall turbulence. J Fluid Mech 1995;295:199±227. [57] Perot B, Moin P. Shear-free turbulent boundary layers. Part 2: New concepts for Reynolds stress transport equation modelling of inhomogeneous ¯ows. J Fluid Mech 1995;295:229±45. [58] Peyret R, Taylor TD. Computational methods for ¯uid ¯ow. New York: Springer; 1983. [59] Reynolds WC. The potential and limitations of direct and large eddy simulations. In: Gatski TB, Sarkar S, Speziale CG, editors. Advances in turbulence. New York: Springer; 1990. p. 88±104. [60] Richter K, Friedrich R, Schmitt L. Large-eddy simulation of turbulent wall boundary layers with pressure gradient. In: Durst, et al., editors. Proceedings of the 6th Symposium Turbulent Shear Flows, Toulouse, France, 1987. p. 22-3-1±7. [61] Rogallo RS, Moin P. Numerical simulation of turbulent ¯ows. Annu Rev Fluid Mech 1984;16:99±137. [62] Saad Y, Schultz MH. GMRES: a generalized minimal residual algorithm for solving non-symmetric linear systems. SIAM J Sci Stat Computat 1986;7:856±69. [63] Schlichting H. Boundary-layer theory. New York: McGraw-Hill; 1968. [64] Spalart PR. Numerical study of sink-¯ow boundary layers. J Fluid Mech 1986;172:307±28. [65] Spalart PR. Direct simulation of a turbulent boundary layer up to Reh ˆ 1410. J Fluid Mech 1988;187:61±98. [66] Sumitani Y, Kasagi N. Direct numerical simulation of turbulent transport with uniform wall injection and suction. AIAA J 1995;33:1220±8. [67] Temam R. On an approximate solution of the Navier±Stokes equations by the method of fractional steps. Part 1. Archiv Ration Mech Anal 1969;32:135±53. [68] Van den Vorst HA. BI-CGSTAB: a fast and smoothly converging variant of BI-CG for the solution of nonsymmetric linear systems. SIAM J Sci Stat Computat 1992;13:631±44. [69] Van den Vorst HA, Sonneveld P. CGSTAB, a more smoothly converging variant of CGS. Tech. Report 90-50, Delft University of Technology, 1990. [70] Vincent A, Meneguzzi M. The spatial structure and statistical properties of homogeneous turbulence. J Fluid Mech 1991;225:1±20. [71] Wagner C. Direkte numerische Simulation turbulenter Str omungen in einer Rohrerweiterung. VDI-Fortschrittsberichte, Reihe 7, Nr. 283. D usseldorf: VDI; 1996. [72] Wagner C, Friedrich R. Direct numerical simulation of turbulent ¯ow in a sudden pipe expansion. Application of direct and large Eddy simulation to tansition and turbulence. AGARD-CP-551. Nevilly-sur-Seine: Agard; 1994. [73] Wagner C, Friedrich R. Turbulent ¯ow in a sudden pipe expansion. In: Benzi R, editor. Advances in turbulence V. Dordrecht: Kluwer Academic Publishers; 1995. p. 544±8. [74] Wagner C, Friedrich R. On the turbulence structure in solid and permeable pipes. Int J Heat Fluid Flow 1998;19:459±69. [75] Wagner C, Friedrich R. Direct numerical simulation of turbulent ¯ow through permeable or rough pipes. In: Papailiou KD, et al., editors. Computational Fluid Dynamics '98, vol. 1. New York: Wiley; 1998. p. 238±43. [76] Wesseling P. Introduction to multigrid methods. Chichester: Wiley; 1992.

R. Friedrich et al. / Computers & Fluids 30 (2001) 555±579

579

[77] Wesseling P, Segal A, Kassels CGM, Bijl H. Computing ¯ows on general two-dimensional nonsmooth staggered grids. J Engng Math, in press. [78] Wesseling P, Segal A, van Kan JJIM, Oosterlee CW, Kassels CGM. Finite volume discretization of the incompressible Navier±Stokes equations in general coordinates on staggered grids. Computat Fluid Dyn J 1992;1:27±33. [79] Westerweel J, Adrian RJ, Eggels JGM, Nieuwstadt FTM. Measurements with PIV in fully developed pipe ¯ow at low Reynolds number. Proceedings of the Sixth International Symposium on Applications of Laser Technique to Fluid Mechanics, Lisbon, Portugal. 20±23 July, 1992.