Large-Eddy Simulation of In-Cylinder Flows

Advances in physical models, numerical methods, and computational power together have brought large-eddy simulation (LES) to the point where it warrants serious consideration for computing in-cylinder turbulent flows. This article includes: a discussion of LES and how it differs from Reynolds-averaged Navier Stokes (RANS) modeling; motivation for transitioning to LES for in-cylinder flow in IC engines; and quantitative results for two motored engine configurations. Dynamicsubgrid-scale stress models are emphasized, and results from two different numerical approaches are discussed. Compared to RANS, LES requires a low level of empirical input, provides more complete information on in-cylinder flow structure, and makes previously inaccessible phenomena (e. g. , cycle-to-cycle variability) amenable to analysis.


INTRODUCTION
Multidimensional modeling of in-cylinder flow and other thermal-fluids processes has advanced dramatically over the past decade.Computations of flow, liquid fuel sprays, mixing, and combustion now are feasible for practical engine configurations through full two-and four-stroke engine cycles.Multidimensional models are exercised in capacities that range from the fundamental exploration of in-cylinder physical processes to the design of production combustion systems.As more information is sought from these models (e.g., direct and quantitative information on combustion stability and emissions), it becomes increasingly imperative to strengthen the underlying foundation of physical models and numerical methodology.
Hydrodynamic turbulence remains among the most compelling and ubiquitous challenges.The fundamental equations describing thermal-fluids processes in the engine essentially are known.Gas motion, for example, is described by the compressible Navier-Stokes equations.But at the Reynolds numbers characteristic of in-cylinder flows, apparently chaotic three-dimensional velocity fluctuations occur on length scales that range from a fraction of the bore diameter to sub-millimeter, and on a correspondingly broad range of time scales (Haworth and Jansen, 1998).Computers do not exist, and will not become available in the foreseeable future, that are capable of dealing with this broad range of spatial and temporal scales.
The dynamic range of scales to be resolved is reduced by averaging or filtering the governing equations and modeling the additional terms (apparent ÒstressesÓ arising from fluctuations about the filtered mean) that appear.The prevailing approach is ensemble-or phase-averaging, wherein fluctuations at all scales are dealt with by the turbulence model.This is ÒReynolds-averaged Navier StokesÓ (RANS) modeling.Variants of the two-equation k-e model with wall functions prevail for multidimensional modeling of in-cylinder flows.
The choice of turbulence model has profound implications.The turbulence model drives the modeling of all other physical processes including heat transfer, fuel sprays, mixing, and combustion.Moreover, the approach adopted to reduce the range of scales imposes restrictions on the information that can be extracted from the computation.For example, cycle-to-cycle flow and combustion variability is not accessible from an ensemble-averaged formulation.
It long has been argued (e.g., Reynolds 1976) that turbulence modeling based on spatial filtering offers advantages compared to time-or phase-averaging.These advantages are particularly compelling for the IC-engine application (El Tahry and Haworth, 1992;Haworth and Jansen, 1998).Computations of turbulent flows based on spatially filtered equations are known as Òlarge-eddy simulationÓ or LES.
In this overview of LES for in-cylinder flows in reciprocating IC engines, we begin with a discussion of what is meant by LES in this context.Specific models and numerical methodologies are described next, followed by results for two simple engine configurations.We close with a discussion of outstanding issues and research directions.

LES and RANS
The premise of LES is to resolve as broad a range of spatial and temporal scales as practicable, and to model only the small scales that cannot be resolved.Classic theory suggests that the small-scale structure of turbulence should be more universal than the large-scale structure (Tennekes and Lumley, 1972); in that case, LES turbulence models are expected to be more general than RANS models.This idea is not new.As early as the 1960Õs (Smagorinsky 1993) numerical simulations of atmospheric flows incorporated these ideas.Because LES has been developed largely within the context of Eulerian grid-based numerical methods, the resulting turbulence models often are referred to as Òsubgrid-scaleÓ (SGS) or Òsubfilter-scaleÓ models.
An important distinction between RANS and LES is the sense in which the numerical solution ÒconvergesÓ with increasing spatial and temporal resolution.In RANS modeling, mesh and/or time-step refinement improves numerical accuracy but does not inherently increase the dynamic range of scales that is resolved; a consistent RANS methodology converges to an exact solution of the filtered Navier-Stokes equations.In LES, by contrast, smaller spatial/temporal scales are resolved as the grid spacing/computational time step are refined.A consistent LES converges to an exact solution of the unfiltered Navier-Stokes equations: that is, to direct numerical simulation (DNS).This suggests that high-fidelity numerical methods are an important component of LES.
A drawback of early subgrid-scale models was the explicit appearance of mesh size (filter width) as a length scale in the modeled turbulent stress.The specification of a length scale that was unrelated to local flow conditions limited the utility and generality of these models.A significant theoretical advance in LES is the ÒdynamicÓ coefficient-determination procedure of Germano et al. (1991).In dynamic models, the subgrid-scale stress is prescribed based on a local assessment of how well resolved the flow actually is.This reduces the level of empiricism, and helps to isolate the physical turbulence model from inaccuracies in the numerical method.It is theoretical contributions of this kind, together with significant advances in computational power and numerical methodology, that together have brought LES to the point where it warrants serious consideration for the in-cylinder application.

What is LES?
There is a broad and growing published literature in LES, and several ÒschoolsÓ have emerged.One approach has been to use no explicit turbulence model at all, the rationales being: Ð numerical inaccuracy provides a subgrid-scale model that is at least as reliable as any explicit model; Ð the numerical dissipation of the low-order numerical methods generally used in engineering computational fluid dynamics (CFD) dominates the modeled subgridscale stress; Ð numerical inaccuracy behaves in a manner that is qualitatively similar to an explicit SGS model.These points are not without merit; Òmodel-lessÓ LES has been legitimized by Boris et al. (1992) and Oran and Boris (1993).An example for in-cylinder flow (including an explicit flame-propagation model) is the work of Naitoh et al. (1992).In other schools, the term ÒLESÓ is reserved for simulations using extremely high-fidelity non-dissipative numerical methods.
To the extent possible, we seek to isolate the physical model from inaccuracy in the numerical method.Thus for present purposes, LES implies: Ð an explicit physical model for the turbulent subgrid-scale stress; Ð a numerical methodology adequate to distinguish the effects of the physcial model from those of numerical inaccuracy; Ð a consistent formulation (including boundary conditions) that converges to DNS.To maintain temporal resolution commensurate with the spatial scales resolved, the computational time step generally must correspond to a material-speed Courant number not exceeding unity.
It is equally important to be clear about what it is that one seeks from LES; that is, how do we know when we get it right?Here the criterion is that computed one-point statistics (first and second moments, at a minimum) agree with suitably filtered experimentally measured, DNS-derived, or theoretical values.In its application to in-cylinder flow in IC engines, this implies that computations must extend through multiple engine cycles to minimize contamination by the initial conditions.An alternative, for example, is Òvery-large-eddy simulationÓ (VLES): there a single engine cycle is computed and standard turbulence wall functions are used.VLES has been explored for in-cylinder flows by Rutland (1998).

LES and Engines
Many flows of practical interest are statistically stationary and/or spatially one-or two-dimensional; RANS modeling is especially appealing in such cases.LES, by contrast, is inherently three-dimensional and transient.Unsteady threedimensional RANS has been the approach of choice for incylinder flows for many years.Moreover, it has been argued (El Tahry and Haworth, 1992;Haworth and Jansen, 1998) that the computational meshes typically used for RANS modeling of practical in-cylinder configurations (currently 10 5 to 10 6 computational elements with second-order spatial discretization) already are sufficient to capture 80-90% of the flowÕs kinetic energy.
In-cylinder Reynolds numbers are proportional to crankshaft rotational speed.For an automotive-scale engine at a low-speed light-load operating condition, the integralscale-based turbulence Reynolds number Re T at the time of ignition (about 30 CAD before piston top-dead-centerÑ TDC) is approximately equal numerically to the engine rotational speed in revolutions-per-minute: that is, Re T is in the low thousands (Haworth and Jansen, 1998).These lowto-moderate Reynolds numbers make the engine an attractive candidate for LES.
In these respects, the transition to LES is natural.Several important considerations remain, however.The turbulent stress model (including wall treatment) needs to be reformulated.The numerical methodology needs to be assessed for suitability to LES.Initial and inflow/outflow boundary conditions need to be modified.Multiple engine cycles must be computed 1 .
(1) While a single engine cycle normally is computed in RANS, multiple engine cycles should be considered there as well.Without extreme care in the specification of initial conditions, it is unlikely that the first computed cycle will correspond to the desired ensemble mean.
And, as we move beyond consideration of the gas-phase turbulent flow alone, the other physical models (e.g., fuel sprays and combustion) need to be adjusted.
In return for this effort, LES promises a significantly improved representation of the in-cylinder turbulent flow field.Even if practical considerations limit the number of engine cycles computed to a small number, the simulated flow field contains substantially more information than a corresponding RANS: in principle, all of the information that can be resolved meaningfully using the specified numerical methodology, mesh, and computational time step is retained.Moreover, LES allows consideration of previously inaccessible physical phenomena.Cycle-to-cycle variability finally can be examined directly.And outstanding issues in turbulence structure and modeling (e.g., temporal vs spatial vs ensemble averaging) are amenable to analysis.

SUBGRID-SCALE MODELS
All models discussed here are of the Smagorinsky type.The influence of unresolved (subgrid-scale) motions on the resolved scales is treated as an additional viscosity.The subgrid-scale stress tensor t SGS ji has a form identical to that of the viscous stress: (1) The subgrid-scale viscosity m SGS is proportional to a norm of the local resolved rate-of-strain tensor and to a filter width D, (2) The single model coefficient is C S .

Constant-Coefficient Model
In the constant-coefficient Smagorinsky model, C S is a constant and D is proportional to the local mesh spacing.To accommodate non-uniform meshes, D is specified as , where V is the volume associated with a computational node or element.Modelers have adopted different values of C S to match benchmark data in different flows.The non-universality of C S and the requirement that the filter width D be specified explicitly motivate the need for a more general approach.et al. (1991) proposed a method for evaluating subgrid-scale model coefficients from information contained in the resolved fields.In their ÒdynamicÓ procedure, two m r

Germano
different filter widths are introduced, and , where . Quantities filtered at the smaller scale are denoted by the overbar notation while the hat notation denotes filtering at the larger scale.In usual practice, the smaller filter is implicit in the numerical method; that is, the quantity denotes the computed velocity delivered by the numerical method at position x and time t.The second filter, or Òtest filter,Ó corresponds to a second explicit filtering operation at the larger scale.Thus represents the LES-computed velocity field filtered at scale .An exact relationship can be derived between the subgridscale stress tensors at the two different filter widths (the ÒGermano identityÓ).Substitution of a Smagorinsky form for the subgrid-scale stress into the Germano identity, along with some additional assumptions (Lilly, 1992), leads to a closedform local expression for the quantity : (3) where L ij and M kl are second-order tensors that are computable from the LES-resolved velocity field: (4) (5) Here L ij is a measure of the turbulent stress between the two filter widths, while M kl is a function of the resolved rateof-strain.It is significant that the dynamic model returns not simply the value of the model coefficient C S , but the product of C S and the square of the filter width ; only the filterwidth ratio needs to be specified in M kl .It can be argued that the ratio should be more uniform than itself (Haworth and Jansen, 1998).
The angled brackets in Equation (3) represent averages over homogeneous directions.Unfortunately, there generally are no homogeneous directions in flows of engineering interest.A third model variant overcomes this difficulty.Meneveau et al. (1996) proposed to accumulate the averages required in the dynamic model over flow pathlines rather than over directions of statistical homogeneity.This again leads to a closed-form expression for that involves the two second-order tensors of Equations ( 4) and ( 5).The Lagrangian form of the dynamic model requires the solution of two additional transport equations for quantities that represent weighted averages of L ij M ij and M kl M kl over fluid-particle trajectories.In addition, it requires the specification of a relaxation time scale that corresponds to a Lagrangian memory time for fluid elements.The resulting model is applicable to arbitrary statistically non-homogeneous turbulent flows.

Comments on Dynamic Models
Dynamic determination of subgrid-scale model coefficients offers several advantages.The subgrid-scale stress increases locally in areas of low grid resolution in response to the high activity found between the two filter scales (large L ij ), and decreases to zero in case all scales of motion are fully resolved locally L ij ® 0. The dynamic model as outlined here has the correct asymptotic behavior close to solid walls; explicit wall functions or damping functions are not needed, provided that mesh resolution is sufficient near walls.And the effective local filter width is an output of the model.
Dynamic Smagorinsky models require the formation of two second-order tensors at every computational node.These tensors are functions of the resolved velocity field and the resolved rate-of-strain field in a neighborhood of the node.The extent of the neighborhood depends on the nature of the second (test) filter adopted.Usually only nearest neighbors are used, for convenience.Jansen (1994) discusses several alternatives for test filtering on unstructured grids.Finally, averaging over homogeneous directions or over flow pathlines is used to compute the local value of .In the Lagrangian form of the dynamic model, two new scalar transport equations are solved.High numerical accuracy is not required here, as these quantities essentially serve to establish an upstream-averaging Òdomain of influenceÓ; a simple local time-explicit upwinding procedure can be employed (Meneveau, 1996).Thus even in its Lagrangian form, the computational overhead for the dynamic Smagorinsky models is not substantial.In fact, the operations required are computationally less intensive than solving the usual k-e equations.The SGS models are implemented simply by modifying the local fluid viscosity.
Although much of the formal theory and development has been done in the context of non-dissipative spectral CFD methods, the dynamic models are well suited to grid-based computation as well (structured or unstructured; dissipative or non-dissipative; finite-difference, finite-volume, or finiteelement).Numerical dissipation is ÒsensedÓ as a lack of activity between the two filter scales, and the level of subgrid-scale stress is reduced accordingly.In the case of dissipative methods, the dynamic procedure should be interpreted as providing subgrid-scale viscosity over and above that implicit in the discretization scheme.

NUMERICAL METHODOLOGY
The first-order time, second-order space discretizations typically employed in engineering RANS computations for incylinder flows are highly dissipative and not well suited to maintaining the high temporal fidelity required for LES.Initial explorations in canonical flows (Section 4.1) yielded poor results, and these methods were abandoned.Here discussion is limited to two relatively new alternative numerical approaches.

Boundary-Fitted Method
The Chad (Computational Hydrodynamics for Advanced Design) CFD code employs a variable explicit/implicit advection scheme that differences along characteristic directions (O'Rourke and Sahota, 1996).Second-order spatial and temporal accuracy are maintained, provided that the material-speed Courant number is less than unity.Because upwinding is used, the method is dissipative.Node-based variables are used with an edge-based data structure, allowing fully unstructured meshes.The formulation is fully compressible, and allows arbitrary mesh motion to accommodate moving piston and valves.This code was designed for distributed-memory scalable parallelism; computational meshes have been partitioned using Metis (Karypis and Kumar, 1995).Implementation of SGS models is detailed in Haworth and Jansen (1998).

Body-Force Method
A quite different computational approach is being pursued by Verzicco et al. (1998).Here a structured orthogonal mesh is used.Boundary conditions are imposed by assigning body forces in a manner that enforces a prescribed velocity along surfaces that need not coincide with grid coordinate lines (Mohd-Yusof, 1996).In the case of a reciprocating IC engine, for example, the piston and valves move through a fixed mesh.Spatial discretization is via a non-dissipative staggered second-order centered finite-difference scheme (Verzicco and Orlandi, 1996); explicit third-order Runge-Kutta time integration is used; and the formulation is incompressible.The combination of structured mesh and incompressibility make this approach extremely efficient computationally: simulations are run on PCÕs or small workstations.At the same time, strict incompressibility limits its utility for practical IC engines.A low-Machnumber variable-density version (r = r(t)) is being explored (Verzicco, 1998).
While body-force methods are not new, the Mohd-Yusof formulation has the significant advantage of not reducing the stability limit of the time integration scheme.This renders practical the simulation of three-dimensional high-Reynolds-number flows.It is well suited to dynamic subgrid-scale models.Dynamic models automatically account for the presence of walls (provided that mesh resolution is adequate) without requiring explicit wall-normal-distance information, which is not readily accessible.Compared to an unstructured boundary-fitted deforming-mesh, a more uniform spatial filtering results; temporal and spatial variations in the LES filter width are not an issue.

Validation for Canonical Flows
Before attempting engine-like configurations, the unstructured LES methodology (Section 3.1) was assessed for hierarchy of canonical flow configurations (Haworth and Jansen, 1998).Evolution of turbulence kinetic energy and energy spectra was compared to experimental measurements for decaying homogeneous isotropic turbulence (Comte-Bellot and Corrsin, 1971).Evolution of the normalized anisotropy tensor was compared to results from rapiddistortion theory for homogeneous linear expansions and compressions of initially isotropic turbulence (Kassinos and Reynolds, 1994).And profiles of mean velocity, Reynolds stresses, normalized third and fourth moments, and rms pressure fluctuations were compared with DNS results for a planar channel flow (Kim et al., 1987).
The salient results of this exercise are as follows: Ð Generally satisfactory evolution of first and second moments was found in all cases.Ð For the planar channel flow, significant differences were found between a constant-coefficient model and the dynamic models, with the dynamic models showing better agreement with DNS data.Ð The effects of the physical subgrid-scale model were discernable; but numerical inaccuracy emerged as the dominant contributor to remaining discrepancies between model and measurements.

A Reciprocating Piston-Cylinder Assembly with Fixed Central Valve
Also reported in Haworth and Jansen (1998) were initial comparisons between model and measurement for the axisymmetric piston-cylinder assembly of Morse et al.
There is no compression.The piston is driven in simple harmonic motion at a speed of 200 tr/min (mean piston speed V p = 0.4 m/s).Laser-Doppler anemometry has been The axisymmetric reciprocating piston-cylinder assembly of Morse et al. (1978).
used to obtain ensemble-mean (phase-averaged) radial profiles of mean and rms axial velocity at 10 mm axial increments starting from the head for crank positions of 36¡, 90¡, 144¡, and 270¡ after piston TDC.This configuration has been the subject of exhaustive numerical studies using a variety of RANS models and numerical methods (El Tahry, 1985;El Tahry and Haworth, 1992).Initial LES results for a mesh of 151.620 nodes using a fixed-coefficient Smagorinsky model (C S = 0.1 2 ) and the numerical methodology of Section 3.1 were reported by Haworth and Jansen (1998).Computations were carried through four engine cycles; the first cycle was discarded to avoid contamination by initial conditions.Mean quantities then were extracted by azimuthal averaging and ensemble averaging over the remaining three engine cycles.These results are indicated by the dashed lines in Figure 2.  Two additional sets of calculations are shown in Figure 2. Solid lines correspond to a second unstructured boundaryfitted computation with a somewhat finer mesh (258 264 nodes) and the Lagrangian dynamic Smagorinsky model; here five engine cycles have been used in the ensemble averaging (after discarding the first cycle).And the chain-dashed lines are results obtained using the body-force method (Section 3.2; Verzicco et al., 1998), a mesh of about 1.6×10 6 nodes, and an azimuthally averaged dynamic Smagorinsky model; many of the mesh points lie outside of the flow domain.
Several observations are made from Figure 2. First, results are in all cases better than any published results obtained using RANS models (not shown; see Haworth and El Tahry, 1992).The agreement between computed and measured rms velocities is particularly compelling.At measurement stations corresponding to the tip of the flapping annular intake jet, for example (e.g., 36¡, 20 mm below the head), the rms velocity is comparable to the local mean velocity magnitude.Second, boundary-fitted-mesh results obtained using the Lagrangian dynamic model show improvement over those obtained using the fixed-coefficient model and a coarser mesh.And finally, the rms levels obtained using the body-force method and a very fine mesh are consistently higher than those of the boundary-fitted method; and the former agree better with measurements for stations close to the head at 144¡.Since comparable subgrid-scale models have been used for the body-force calculations and for the finer-mesh boundary-fitted calculations, it probably is the combination of the non-dissipative numerical method and the significantly finer mesh for the former that is responsible for the better agreement.
A final plot is included for this configuration to make some points concerning the magnitude of the subgrid-scale viscosity in LES. Figure 3 shows the evolution over (almost) five full engine cycles of the mean and maximum subgridscale viscosity.The subgrid-scale viscosity from the dynamic models tends to be strongly intermittent, with the highest values very localized in space and in time.In Figure 3, the maximum viscosity in the computational domain is one-to-two orders of magnitude higher than the volume-averaged mean viscosity.
The mean subgrid-scale viscosity varies by a factor of three over the engine cycle.Here the volume-averaged mean viscosity is only one-to-three times the fluid viscosity; the ratio is somewhat higher for the the in-cylinder region alone.The computational domain extends upstream of the chamber and valve, and includes a large plenum where there is little activity.In the computations for the canonical flows of Section 4.1, the volume-averaged subgrid-scale viscosity typically was between five and ten times the fluid viscosity.By comparison, the global in-cylinder effective viscosity near TDC for a k-e RANS model can be estimated as m k-e = C m k 2 / e » 0.2 hV p = 0.0024 m 2 /s more than two orders of magnitude higher than the fluid viscosity.Both the maximum and the mean show cycle-to-cycle variation in Figure 3.  (1978).Computations started at Ð360 degrees.Solid line: 0.1* maximum value in the computational domain.Dashed line: volume-weighted mean value over the computational domain.
Further results including cycle-to-cycle and azimuthal variability can be found in Haworth Jansen (1998) and Verzicco et al. (1998).

A Two-Valve Pancake-Chamber Motored Four-Stroke-Cycle Engine
Boundary-fitted, Lagrangian dynamic subgrid-scale model LES is in progress for a transparent-combustion-chamber (TCC) motored two-valve four-stroke-cycle engine where extensive particle-image velocimetry (PIV) measurements are available.The geometric configuration is identical to that published by Reuss et al. (1995), with the exception that the intake valve is not shrouded (Fig. 4).The engine is motored at 1200 tr/min while intake-and exhaust-plenum pressures are held at 40 kPa and 105 kPa, respectively.Examples of measured velocity fields on a cutting plane normal to the cylinder axis 6.1 mm below the head at piston TDC are given in Figure 5 (Reuss, 1998).A RANS computation (e.g., k-e) is expected to reproduce the ensemble average, which shows only a coherent large-scale swirl (see Reuss et al., 1995, for example).Figure 5 illustrates dramatically that the instantaneous flow structure on any individual engine cycle bears little resemblance to the ensemble average.
The nature of cycle-to-cycle flow variations is of particular interest.Fansler and French (1992), for example, reported bimodal probability density functions (pdfÕs) of velocity at Outer surfaces of the computational mesh for the TCC engine (Reuss et al., 1995).some measurement locations near port/cylinder interfaces in a loop-scavenged two-stroke-cycle engine.And Reuss (1998) has found that for some TCC intake configurations, flow visualizations from individual cycles can be ÒsortedÓ into a relatively small number of patterns.These are phenemena that we wish to explore and quantify using LES.
Computations were initiated at intake bottom-dead-center (BDC, Ð180¡) using an extremely coarse mesh of 25 176 nodes.Results were interpolated onto a finer mesh of 178 515 nodes at 630¡; and both the coarse-and fine-mesh computations continue from that point on.This is a challenging flow computationally.Because the engine is motored, there are large mismatches between the in-cylinder and intake/exhaust port pressures as the intake and exhaust valves open: the flow is choked for 20 to 30 CAD following each valve-opening event, including the 20¡ duration of valve overlap.
Only a few early results ar reported here.Figure 6 shows the computed evolution of the in-cylinder swirl and tumble ratios through four engine cycles 2 .For this motored operating condition, the in-cylinder motion induced by back-flow through the exhaust valve following exhaust-valve opening is nearly as strong as the intake-generated flow.The TDCcompression swirl ratio in the measurement plane of Figure 5 determined from the ensemble-averaged PIV is 0.6 (ranging from 0.3 to 2.0 on individual cycles); the steady-flow-rig measured swirl ratio is 0.7 (Reuss, 1998).These values are in good agreement with the computed TDC-compression values in Figure 6.Even the global quantities of Figure 6 exhibit variation from one engine cycle to another.It is not yet clear whether there are specific long-term trends or whether the system has reached a statistically steady state.
(2) Swirl ratio is defined as the angular momentum of the in-cylinder fluid about the cylinder (z) axis, divided by the moment of inertia about that axis, and normalized by the crankshaft rotational speed.Tumble ratios about the x and y axes are defined analogously.The origin is the instantaneous center-of-mass of the in-cylinder fluid.Figure 6 Computed (coarse mesh) in-cylinder swirl and tumble ratios for the motored TCC engine at 1.200 tr/min.Calculations began at Ð180¡.
Figure 7 shows instantaneous coarse-mesh computed TDC-compression flow fields on the same cutting plane as Figure 5. On such a coarse mesh, one cannot expect to see the richness of flow structure evident in the experimental measurements.Nevertheless, several features are noteworthy.First, the computed velocity magnitudes are comparable to measured values.Second, there clearly is structure that is not present in the measured ensemble average; features of 5 to 10 mm in extent are evident.And third, the computations exhibit cycleto-cycle variation.In remains to quantify these variations and to determine to what extent they are consistent with those observed experimentally.

NEXT STEPS
At the time of this writing, coarse-and fine-mesh computations continue for the TCC engine.These are being carried through multiple engine cycles, and the results Computed (coarse mesh) instantaneous velocity fields on a cutting plane normal to the cylinder axis, 6.1 mm below the head, at TDCcompression, on four consecutive engine cycles.Computed velocity vectors have been interpolated onto a uniform 30 x 30 mesh.The full 92 mm bore is shown.
processed to extract information on in-cylinder turbulence structure, including differences among ensemble-, time-, and spatial-averaging.Quantitative comparisons with PIV measurementsÑincluding cycle-to-cycle variationÑare in progress.Assuming a satisfactory outcome for the TCC, the next step will be motored calculations for a production engine configuration.Care must be exercised when making quantitative comparisons between LES and experimental measurements, particularly when the LES is limited to a small number of engine cycles: suitable conditioning variables for conditional averaging need to be established.
Both of the numerical methodologies mentioned here have, in conjunction with dynamic subgrid-scale turbulence models, demonstrated good quantitative agreement with measurements for an axisymmetric engine configuration.Second-order space and time appears to be the minimum that allows the effect of the subgrid-scale model to be discerned from that of numerical inaccuracy; higher-order would be preferred.Non-dissipative schemes may have an advantage for LES over methods that use upwind biasing (Mittal and Moin, 1997); a dissipative method will tend to contaminate the energy spectrum at high wave numbers, for example (Haworth and Jansen, 1998).On the other hand, skew-upwind schemes (El Tahry, 1985) and the characteristic-based advection scheme of Chad tend to advect flow ÒstructuresÓ (vortices) with high phase accuracy.The non-dissipative body-force method described here offers promise, but remains to be demonstrated for variable-density compressible flows.Clearly, numerical methodology remains one of the outstanding issues in LES for practical flow configurations.
As the viability of LES for gas-phase turbulent flow in engines becomes established, attention will turn to SGS models for fuel sprays and combustion.The issues that arise in filtering and modeling highly nonlinear source terms in LES are essentially the same as those that must be dealt with in three-dimensional time-dependent RANS.Provided that one remains in a regime where the key reaction and/or liquid spray processes to be modeled occur at unresolved scales (which is likely to remain the case for the foreseeable future), the essential difference lies in the specification of turbulence length and time scales.This implies a different standard for LES for two-phase reacting flows: namely, a relaxation of the requirement of convergence to DNS with mesh and time-step refinement.Even with essentially the same combustion and spray models, improvements are anticipated from the (presumably) better representation of the turbulent flow field.

CONCLUSION
Significant progress has been made in subgrid-scale models and numerical methodologies for LES of in-cylinder flows in reciprocating IC engines.Quantitatively accurate results can be obtained using practical mesh densities and readily implemented dynamic subgrid-scale turbulence models.Even when practical considerations limit the number of engine cycles that can be computed to a small number, LES offers a more realistic representation of the in-cylinder turbulent flow compared to RANS.
The nature of cycle-to-cycle variability of in-cylinder flow and combustion has long been a contentious issue in the engine community.The present results support the view that cycle-to-cycle flow variability is nothing more than largelength-scale long-time-scale turbulence.That is, no special forcing or perturbation is needed to generate cycle-to-cycle variations whose character is at least qualitatively similar to that observed experimentally; quantitative agreement has been demonstrated for one simple configuration (Morse et al., 1978).Cycle-to-cycle variability is implicit in the (filtered) Navier-Stokes equations at engine Reynolds numbers.The two prerequites to capturing this variability computationally appear to be geometric fidelity and an acceptable level of numerical accuracy.The role of specific phenomena such as acoustic waves in the induction system (Poinsot, 1998) remains to be established.
Mesh requirements for LES are not substantially greater than those for RANS by virtue of the low-to-moderate incylinder Reynolds numbers; the computational time step limitation (material Courant number of unity) is considerably more restrictive.This, in combination with the requirement to compute multiple engine cycles, results in long computational times for LES.
Transition from RANS to LES probably is inevitable in the IC-engine application.Initially LES will be used to address fundamental aspects of in-cylinder turbulent flow structure, but eventually it is expected to become the engineering turbulence model of choice.The transition will be paced by the availability of robust, computationally efficient, and high-fidelity numerical methods capable of handling the geometric configurations and flow regimes encountered in reciprocating IC engines.

Figure 3
Figure 3 Computed subgrid-scale viscosity, normalized by fluid viscosity, for the axisymmetric piston-cylinder of Morse et al.(1978).Computations started at Ð360 degrees.Solid line: 0.1* maximum value in the computational domain.Dashed line: volume-weighted mean value over the computational domain.

Figure 5 PIV
Figure 5 PIV velocity fields on a cutting plane normal to the cylinder axis, z = Ð6.1 mm, at piston TDC (Reuss, 1998).A 70 mm circle centered in the 92 mm bore is shown.Top: ensemble average, 90 engine cycles.Middle/bottom: instantaneous fields on two different engine cycles.