Numerical Simulations of Incompressible Flows

The following sections are included:

• Dedication

• Biographical Sketch of Dr. Dochan Kwak

• Publications of Dr. Dochan Kwak

• CONTENTS

This paper reviews recent progress made in incompressible Navier-Stokes simulation procedures and their application to problems of engineering interest. Discussions are focused on the methods designed for complex geometry applications in three dimensions, and thus are limited to primitive variable formulation. A summary of efforts in flow solver development is given followed by numerical studies of a few example problems of current interest. Both steady and unsteady solution algorithms and their salient features are discussed. Solvers discussed here are based on a structured-grid approach using either a finite-difference or a finite-volume frame work. As a grand- challenge application of these solvers, an unsteady turbopump flow simulation procedure has been developed which utilizes high performance computing platforms. In the paper, the progress toward the complete simulation capability of the turbo-pump for a liquid rocket engine is reported. The Space Shuttle Main Engine (SSME) turbo-pump is used as a test case for evaluation of two parallel computing algorithms that have been implemented in the INS3D code. The relative motion of the grid systems for the rotor-stator interaction was obtained using overset grid techniques. Unsteady computations for the SSME turbo-pump, which contains 114 zones with 34.5 million grid points, are carried out on SGI Origin 3000 systems at NASA Ames Research Center. The same procedure has been extended to the development of NASA-DeBakey Ventricular Assist Device (VAD) that is based on an axial blood pump. Computational, and clinical analysis of this device are presented.

We trace the evolution of artificial compressibility or preconditioning methods, originally introduced by Chorin to enable time-marching solutions of the incompressible equations and later adapted and generalized to a wide range of flow regimes. In particular, we use asymptotic expansions to investigate the theoretical underpinnings of the approach and to elucidate how it is instrumental in improving the accuracy and efficiency of a wide range of flow computations. These include low Mach compressible flows, viscous-dominated flows, unsteady flows and multi-phase flows. Representative computations are provided throughout to illustrate the efficacy of the techniques.

A family of low-order finite element solvers for incompressible flows is described. Both the advection and divergence terms are treated using consistent numerical fluxes along edges. Several techniques to accelerate convergence to steady state are explored and compared. The techniques are then used in a fully implicit time-marching scheme that solves a steady problem at every timestep. Several examples demonstrate the usefulness of the developed scehmes.

A direct primitive variable solution technique, that does not require either artificial compressibility or a pressure Poisson solver, is discussed for the solution of the incompressible potential, Euler and reduced or full Navier Stokes equations. Both steady state and transient flows are evaluated. The procedure is equally applicable for the solution of compressible flows as given in a number of the references. With the compressible code, the incompressible limit of Mach number identically equal to zero provides results identical with those contained here for the incompressible form of the equations. The convergence properties for the two systems are essentially identical. No additional boundary conditions on pressure or any other flow variable are required in either formulation.

The following sections are included:

• Introduction

• Analysis

• Numerical Evaluation

• Conclusion

• References

Numerous papers have appeared in the literature over the past thirty years discussing projection-type methods for solving the incompressible Navier-Stokes equations. A recurring difficulty encountered is the proper choice of boundary conditions for the auxiliary variables in order to obtain at least second order accuracy in the computed solution. A further issue is the formula for the pressure correction at each timestep. An overview of boundary condition choices that give second-order convergence for all solution variables is presented here based on recently published results by Brown, Cortez and Minion [2].

We describe a split-step finite-difference scheme for solving the incompressible Navier-Stokes equations on composite overlapping grids. The split-step approach decouples the solution of the velocity variables from the solution of the pressure. The scheme is based on the velocity-pressure formulation and uses a method of lines approach so that a variety of implicit or explicit time stepping schemes can be used once the equations have been discretized in space. We have implemented both second-order and fourth-order accurate spatial approximations that can be used with implicit or explicit time stepping methods. We describe how to choose appropriate boundary conditions for the pressure to make the scheme accurate and stable. A divergence damping term is added to the pressure equation to keep the numerical dilatation small. Several numerical examples are presented.

A semi-implicit form of the method of spectral deferred corrections is applied to the solution of the incompressible Navier-Stokes equations. A methodology for constructing semi-implicit projection methods with arbitrarily high order of temporal accuracy in both the velocity and pressure is presented. Three variations of projection methods are discussed which differ in the manner in which the auxiliary velocity and the pressure are calculated. The presentation will make clear that projection methods in general need not be viewed as fractional step methods as is often the practice. Two simple numerical examples are used to demonstrate fourth-order accuracy in time for an implementation of each variation of projection method.

A numerical method for solving the incompressible Navier Stokes equations based on splitting the velocity into a rotational and a potential part via a Helmholz decomposition is presented, using the Lamb form of the Navier Stokes equations. The scheme conserves both mass and momentum to machine accuracy and is not limited to two dimensions. No special handling is required for multipley connected domains, since a Poisson equation is solved for the pressure. Careful selection of farfield boundary conditions eliminate artificial sources of vorticity.

A finite difference discretization on a staggered grid is used to test the scheme. A segregated solution algorithm is adopted and shows that the scheme has good convergence properties. To demonstrate that the method is equally well suited for internal and external flow calculations, flow in a driven cavity, flow over a backward facing step, flow in a cascade and flow over a flat plate are computed.

A recently proposed Helmholz Decomposition method for solving the incompressible Navier Stokes equations is applied to three dimensional flow in a driven cavity. An alternative equation in the form of a convection diffusion equation for the total pressure is implemented to overcome undesirable properties of the Poisson's equation for total pressure. Second order upwinding is used for the rotational term of the Lamb form of the momentum equations and the convection term of the pressure equation. A finite difference discretization on a staggered grid is used and a segregated solution procedure is employed.

The main goal of this article is to discuss the numerical simulation of incompressible viscous fluid flow around moving rigid bodies of elliptical shapes, such as ellipses in two–dimensions, and ellipsoids in three–dimensions. Related simulations for circular and spherical bodies, by methods closely related to those described in this article, have been investigated by the authors in preceding publications, however, the non–sphericity introduces additional difficulties to be discussed in the article. Numerical results for two and three–dimensional test problems will be presented; they include the interaction of two ellipsoids settling in a narrow channel.

The following sections are included:

• Introduction

• Governing Equations

• Mesh Generation

• Finite Volume Solution Algorithm

• Multigrid Acceleration

• Parallelisation

• Numerical Examples

• Conclusions

• Acknowledgements

• References

We provide an overview of the interface-tracking and interface-capturing techniques we have developed in recent years for computation of flow problems with moving boundaries and interfaces. Both classes of techniques are based on stabilized formulations, and determination of the stabilization parameters used in these formulations is also highlighted here. The interface-tracking techniques are based on the Deforming-Spatial-Domain/Stabilized Space-Time formulation, where the mesh moves to track the interface. The interface-capturing techniques, which were developed for two-fluid flows, are based on the stabilized formulation, over non-moving meshes, of both the flow equations and the advection equation governing the time-evolution of an interface function marking the interface location. For interface-capturing techniques, to increase the accuracy in representing the interface, the Enhanced-Discretization Interface-Capturing Technique can be used. We also provide and overview of some of the additional ideas developed to increase the scope and accuracy of these two classes of techniques.

The interface-tracking and interface-capturing techniques we developed in recent years for computation of flow problems with moving boundaries and interfaces rely on stabilized formulations such as the streamline-upwind/Petrov-Galerkin (SUPG) and pressure-stabilizing/Petrov-Galerkin (PSPG) methods. The interface-tracking techniques are based on the Deforming-Spatial-Domain/Stabilized Space-Time formulation, where the mesh moves to track the interface. The interface-capturing techniques, typically used with non-moving meshes, are based on a stabilized semidiscrete formulation of the Navier-Stokes equations, combined with a stabilized formulation of the advection equation governing the time-evolution of an interface function marking the interface location. We provide an overview of the interface-tracking and interface-capturing techniques, and highlight how we determine the stabilization parameters used in the stabilized formulations.

We develop a difference method for the solution of the incompressible Navier-Stokes equations on a staggered grid. It is based on a Padé type fourth order approximation in space and a semi-implicit second order method in time. The method does not allow for any parasitic odd-even oscillatory solutions. We carry out a complete stability analysis for the periodic case, and construct an iterative method for efficient solution of the systems of equations to be solved for each timestep. A few numerical experiments are included for illustration of the efficiency of the iterative method.

B-splines are basis functions for piecewise polynomials having a high level of derivative continuity. They possess attractive properties for complex flow simulations: they have compact support, provide a straightforward handling of boundary conditions and grid nonuniformities, yield numerical schemes with high resolving power, and the order of accuracy is a mere input parameter. This paper discusses some basic approximation properties of B-spline schemes, and their relationship with more conventional numerical schemes. Then, we present a spline collocation method for solving the incompressible Navier-Stokes equations in velocity-pressure formulation by means of the fractional step method.

This paper describes the development and application of a high-order finite-difference approach for the simulation of viscous flows on stretched, curvilinear and dynamic grids. The solver utilizes 4^th- and 6^th -order compact-differencing schemes for the spatial discretization, coupled with both explicit and implicit time-marching methods. Up to 10^th-order, Pade-type low-pass spatial filter operators are incorporated to eliminate the spurious high-frequency modes which otherwise arise due to the lack of inherent dissipation in the spatial scheme. Special attention is given to proper metric evaluation procedures for three-dimensional moving and curvilinear meshes so that the advantages of the high-order approach are retained in practical calculations. The improvements derived from the new high-order methodology are highlighted by considering several applications, including: vortex advection, shear-layer instability, Large-Eddy-Simulations of isotropic turbulence and channel flow, as well as the fluid-structure interaction of boundary layer instabilities with an elastic panel.

High-Reynolds-number flows are simulated by solving the incompressible Navier-Stokes equations. A finite-difference method with third-order upwinding is employed without incorporating a turbulence model. In addition, a multi-directional formulation is used to improve the accuracy of the solution. The validity of this method is thoroughly discussed, and examples are presented to show the applicability of the present approach to a variety of problems using Cartesian and body-fitted coordinates systems. Some of the examples include: a) flows around bluff-bodies using Cartesian coordinates system, and b) flows past subsonic two- and three-dimensional airfoils using O-type coordinates system. The results show that typical flow mechanisms are clearly captured by the present direct numerical simulations. For the airfoil simulations, computed lift coefficients are found to agree quite well with the experimental data even for attack angles above the stall angle.

This paper deals with the lattice Boltzmann approach to CFD, especially to computation of incompressible flow problems. The approach starts from the Boltzmann equation instead of the Navier-Stokes equations on which the conventional CFD is based. Introducing BGK type simplified collision term and choosing appropriate discrete molecular velocity sets, the lattice Boltzmann BGK (LBGK) equation can be derived from the discrete velocity Boltzmann equation using suitable finite difference approximation. Calculations are carried out for a wide variety of both 2-D and 3-D, steady and unsteady incompressible flow problems. Results are compared with those of the Navier-Stokes equations. The lattice Boltzmann method is able to reproduce the dynamics of incompressible flows and could be an alternative to solving the Navier-Stokes equations.

We present a review of the CIP method, which is a kind of semi-Lagrangian scheme and has been extended to treat incompressible flow in the framework of compressible fluid. Since it uses primitive Euler representation, it is suitable for multi-phase analysis. The recent version of this method guarantees the exact mass conservation even in the framework of semi-Lagrangian scheme. Comprehensive review is given for the strategy of the CIP method that has a compact support and subcell resolution including front capturing algorithm with functional transformation.

Vorticity Confinement has been shown, over the last several years, to be a way to quickly and cheaply approximate incompressible flows over complex configurations at high Reynolds number. For these flows, there are many important salient features that are reproduced by the method in a simple way, on relatively coarse uniform Eulerian computational grids with fast low-order computational methods. No complex, high order schemes and no extensive refinement or body conforming grid are required. The basic features of the flow that are very effectively simulated all involve vortical regions: thin boundary layers on solid surfaces that are attached or separate, thin vortex sheets and thin filaments that can convect over long distances with no significant physical diffusion, and, finally, turbulent wakes, including small scale effects. These features are very difficult to treat with conventional schemes. The main point is that almost all of the vortical regions in these flows are either thin, so that their internal structure is not important (analogous to shocks or contact discontinuities in compressible flow), or contain small scale structures, such as turbulent wakes. All of these vortical structures are, of course, embedded in an incompressible irrotational flow field, which is also efficiently computed.

First, the salient features of Vorticity Confinement solutions will be discussed, to give the reader an understanding of the characteristics of the method. Then, a short presentation of the original formulation of the method will be given and a number of recent papers giving more details referenced. This will be followed by some recent results where the method is used to treat the small scales in turbulent wakes in LES-type simulations. Many other recent results concern free vortices and attached flow over complex configurations, as well as validation studies. Papers describing these will be referred to. Finally, a new, simpler formulation will be described that is effective for both Vorticity Confinement as well as convection of thin streams of passive scalars. Preliminary results of this new formulation for convecting vortices and scalars in 2-D will then be presented.

An approach to the optimization of propellers is presented that is based on a helicoidal vortex model and the rigorous calculation of the induced velocities to maximize the thrust for a given output torque of the engine. The mathematical model is presented for inviscid flow and for the viscous correction. The design of a two-bladed propeller at low advance ratio is carried out and the distributions of circulation, induced velocities, and chord are obtained. Comparison with published data indicates that the method is efficient and flexible and produces detailed and accurate results.

A noniterative method for nonlinear parabolic partial-differential equations is described and applied to boundary-layer equations for two-dimensional incompressible laminar flows. Comparison of calculated results indicates that the accuracy of this method is comparable to those obtained with the iterative method.

A noniterative method for nonlinear parabolic partial-differential equations is described and applied to boundary-layer equations for two-dimensional laminar and turbulent flows. Comparison of calculated results indicates that the accuracy of this method is comparable to those obtained with an iterative method.

Parallel computations of viscous turbulent flows over high-lift airfoils are performed using unsteady, incompressible and compressible Reynolds-averaged Navier-Stokes equations with three two-equation turbulence models (the standard k - ε, k - ω , and k - ω SST model). Both the incompressible and compressible flow solvers are developed and demonstrated to investigate the compressibility effects. The compressible code involves an upwind-differenced scheme for the convective terms and a lower-upper symmetric Gauss-Seidel scheme for temporal integration. The incompressible code with pseudo-compressibility method also adopts the same schemes as the compressible code. Both codes are parallel-processed using the message passing interface programming method and show good parallel speedups. The compressible and incompressible codes are validated by predicting the flow around the RAE 2822 transonic airfoil and NACA 4412 airfoil, respectively. In addition, both the incompressible and compressible code using the Chimera overlapping grid scheme are used to compute the flow over the NLR 7301 airfoil with flap and the NASA GAW-1 high-lift airfoil. Compressibility effects on surface pressure coefficients, velocity profiles, and skin friction coefficients are numerically simulated.

The effect of ice shape on aerodynamic performance of ice-contaminated airfoils is investigated by comparing several different ice shapes. The computational method used is based on Navier-Stokes flow physics and validated with experimental data. The effect of turbulence modeling is first studied for different ice shapes to select appropriate models for the iced airfoil study. Results show that ice shapes can be classified into a few geometric categories based on their common characteristics. Aerodynamic performance of iced airfoils is strongly dependent on several key ice shape features such as ice ridges or horns, and the geometry of the suction side of an iced airfoil is more sensitive in determining overall aerodynamic characteristics.

Blood flow in a severely stenotic carotid bifurcation has been analyzed. A realistic three-dimensional model of an atherosclerotic carotid bifurcation was obtained using MR images. An unstructured grid was then generated on the domain and the full Navier-Stokes equations were solved numerically. Steady flows at different Reynolds numbers as well as an unsteady pulsatile flow have been calculated. Non-Newtonian behavior of blood was considered in some simulations. The results obtained are found to be in good qualitative agreement with simulations carried out earlier for two-dimensional flows.

The objective of this work is to determine the behavior of a group of droplets interacting at an intermediate Reynolds flow. The algorithm CHIMERA has been used to accomplish the direct simulation of the vaporization of water droplets.

Results for an isolated drop have been compared with available results to validate the model. We have analyzed the global transport coefficient for groups of droplets in order to determine the agreement with experimental correlations that are commonly used in classic droplet vaporization models.

This paper presents unsteady computations on the flow field around a square cylinder with a gap between the body and the ground plane. Two-dimensional unsteady, incompressible Navier-Stokes codes are developed for the computation of the viscous turbulent flows. The incompressible code with pseudo-compressibility and dual-time stepping method adopts a conventional upwind differencing scheme for the convective terms and DP-SGS scheme for efficient temporal integration. The grid system around a two-dimensional square body is efficiently generated by adopting Chimera hole concept using IB array. By computing the flow around a square cylinder without ground effect, three two-equation turbulence models are evaluated and the developed code is validated. The results show a good agreement with experimental values and other computational results. Critical gap height at which the formation of Kàrman vortex streets is interrupted, is demonstrated by the present numerical investigations.

Simulations of an underwater missile hot launch exhausting plume from a Concentric Canister Launcher (CCL) using a computer code, BUB2D, based on a generalized formulation of hydrodynamic free surface problems, are described. In this paper the missile motion is not considered. Three cases were calculated: (1) a prescribed pressure history at the CCL exit, (2) considering the CCL as a chamber with prescribed rocket motor exhaust at one end, and (3) an interacting model, which consists of a gas dynamic model in the CCL chamber as a subsystem and the gas plume and water above the CCL device as another subsystem with the two subsystems interacting. The computed results are described.

In this paper, we will show the ideas and procedure how the AUSM-family schemes are extended for solving two phase flows equations. Specifically, we shall consider three mathematical models for describing multiphase flows, the single phase model with homogeneous phase transition, two phase mixture model, and two phase gas-solid model. Details of the numerical formulation and issues involved are discussed and the effectiveness of the method are demonstrated for several practical examples.

The aim of this research is to provide a physical complete and numerical efficient simulation method to predict developed cavitation in hydrodynamic turbomachinery as well as in micro fluid dynamic applications, e.g. in high pressure injection nozzles of combustion engines. Cavitating two-phase flows are always very unstable, highly unsteady, 3-D and turbulent. To separate and to understand the inertia controlled cavitation dynamics and its interaction with viscous effects like boundary layers and separation, we start with Euler simulations by neglecting the viscosity of the fluid. Next we introduce the single-phase turbulence k - ω model of Wilcox [1] without modifications with respect to dispersed structures of bubbly liquids, which overestimates viscous effects in the transitional regime between the vapor and liquid phase and tends to suppress typical cavitation instabilities. Consequently our third approach consists of modifications of the single-phase Wilcox model to account for the strong nonlinear variation of the turbulent viscosity µ_t, depending in the local void fraction α.

The key issue of all numerical methods for simulation of cavitating flows is the treatment of the sudden density change of the fluid, in cold water up to 40.000:1, embedded in a global incompressible liquid flow. Here the two-phase fluid is modeled as dispersed mixture of an incompressible liquid and tiny vapor bubbles which grow or collapse, accordingly to the local static pressure and their convective transport. Therefore, the standard VOF method for capturing distinct interfaces without phase transition, e.g. free surface flow or single bubbles, is extended to include phase transition of dispersed mixtures. For simulation of bubble dynamics we apply the Rayleigh equation, which is completed by an energy balance to account for thermal effects, if hot water or if technical fluids others than water, e.g. refrigerants, with high vapor densities are considered.

By using our CFD tool CAVKA we present examples of cavitating flow around hydrofoils and through single hole injection nozzles. Comparing Euler and single-phase turbulence simulations, results based on the inviscid approach are closer to experiments, which indicates an interesting option to reduce computational time.

Recent studies suggest that the wave / vortex eigenmodes' coupling occurs in the linear stage of disturbance development under the influence of centrifugal forces. The boundary layers of two different types provide pertinent examples. The first one relates to a three-dimensional motion with crossflow on a flat plate. The curvature of streamsurfaces naturally warping in the crossflow direction maintains centrifugal forces. The two-dimensional boundary layer on a concave cylindrical surface illustrates the second example. Here centrifugal forces are supported by the fixed curvature of a solid surface bending in the streamwise direction. In both cases centrifugal forces are balanced out by the normal-to-wall pressure gradient. The inclusion of the vortex eigenmodes brought about by centrifugal forces makes asymptotic models self-consistent, with the Cauchy problems well posed in the limit of large Reynolds numbers. An inference of conceptual importance from the asymptotic models is that the boundary layers under consideration suffer absolute instability in the streamwise direction which leads to earlier transition.

The convergence behavior and solution accuracy of a preconditioned characteristic–based viscous flow solver are evaluated for steady laminar compressible flow past a flate plate, over a range of Mach number and wall temperature ratio (0.1 ≤ M_∞ ≤ 3.0 , 0.1 ≤ T_w/T_∞ ≤ 10.0). The flat–plate is useful as a test problem since easily computed similarity solutions for arbitrary Mach number and wall temperature ratio are available that themselves have been validated against experimental measurements. The Navier–Stokes solutions are evaluated for asymptotic convergence behavior, spatial discretization error, and are assessed for conformance to similarity and agreement with the similar solutions. The algorithm converges well for all Mach numbers tested when the wall temperature ratio is between 0.5 and 2.0, and a reduced time step gave solutions for temperatue ratios of 0.1 and 10.0. The validation comparisons indicate consistent and acceptable accuracy over a broad range of flow speeds and temperature ratios, subject to the theoretical thin–layer similarity assumption and to a lesser extent the grid resolution.

Topological methods are used to establish global and to extract local structure properties of vector fields in axi-symmetric and 3-d flows as function of time. The notion of topological skeleton is applied to the interpretation of vector fields generated numerically by the Navier-Stokes equations. The flows considered are swirling jets with super-critical swirl numbers that show low Reynolds number turbulence in the break-up region.

Data compression for the incompressible flow solutions based on Supercompact multi-wavelets is presented. The proposed multiresolution method with supercompact wavelets offers high data compression for fluid simulation and experimental data including incompressible flow solutions. Supercompact wavelets provide advantageous benefits that it allows higher order accurate representation with compact support and therefore, it avoids unnecessary interaction with remotely located data such as across vortices. Thresholding for data compression is applied based on a covariance vector structure of multi-wavelets. Several numerical tests demonstrate large data compression ratios for the outputs of flow field simulation with various levels of fidelities.

The 'Drifting Cup in a Meandering Stream' was initiated in the ancient China by a famous poet and calligrapher called Wang Xhi Zhi. He summoned 42 social celebrities for a spring party held at Ranting (Orchid Pavilion) where a meandering stream was pulled from a nearby creek. They seated themselves along the stream to compose Chinese poems. Cups filled with rice wine were drifted on the stream and the rule was to finish the poem before the cup arrived at the seat of each poet. If a poet had failed, he had to pardon himself by drinking three cups of rice wine. This witty archaic event has gained popularity afterwards and spread to neighbor countries to become a fabulous culture long enjoyed by the rulers in this region. Many Meandering Streams have been built even up to relatively modern times. The author has looked for a list of the Meandering Streams in China, Korea and Japan, and then visited on leisure time many of those historic remains to examine their designs and functions. The author has investigated by experimental model and computer simulation the flow of some Meandering Streams. The author found that there had been a crucial change in the configuration of the Meandering Streams in the history in order to control motion of the drifting cup. Curiously enough, Posuk-Chung Pavilion located in Kyongju, Korea, is in the midst of such morphological transformation.