Parallel computation of unsteady three-dimensional incompressible viscous flow using an unstructured multigrid method

The development and validation of a parallel unstructured tetrahedral non-nested multigrid (MG) method for simulation of unsteady 3D incompressible viscous flow is presented. The Navier-Stokes solver is based on the artificial compressibility method (ACM) and a higher-order characteristics-based finite-volume scheme on unstructured MG. Unsteady flow is calculated with an implicit dual time stepping scheme. The parallelization of the solver is achieved by a MG domain decomposition approach (MG-DD), using the Single Program Multiple Data (SPMD) programming paradigm. The Message-Passing Interface (MPI) Library is used for communication of data and loop arrays are decomposed using the OpenMP standard. The parallel codes using single grid and MG are used to simulate steady and unsteady incompressible viscous flows for a 3D lid-driven cavity flow for validation and performance evaluation purposes. The speedups and efficiencies obtained by both the parallel single grid and MG solvers are reasonably good for all test cases, using up to 32 processors on the SGI Origin 3400. The parallel results obtained agree well with those of serial solvers and with numerical solutions obtained by other researchers, as well as experimental measurements.     

Incompressible flow computations over moving boundary using a novel upwind method

In this paper a novel method for simulating unsteady incompressible viscous flow over a moving boundary is described. The numerical model is based on a 2D Navier–Stokes incompressible flow in artificial compressibility formulation with Arbitrary Lagrangian Eulerian approach for moving grid and dual time stepping approach for time accurate discretization. A higher order unstructured finite volume scheme, based on a Harten Lax and van Leer with Contact (HLLC) type Riemann solver for convective fluxes, developed for steady incompressible flow in artificial compressibility formulation by Mandal and Iyer (AIAA paper 2009-3541), is extended to solve unsteady flows over moving boundary. Viscous fluxes are discretized in a central differencing manner based on Coirier’s diamond path. An algorithm based on interpolation with radial basis functions is used for grid movements. The present numerical scheme is validated for an unsteady channel flow with a moving indentation. The present numerical results are found to agree well with experimental results reported in literature.     

Parallel-multigrid computation of unsteady incompressible viscous flows using a matrix-free implicit method and high-resolution characteristics-based scheme

A three-dimensional parallel unstructured non-nested multigrid solver for solutions of unsteady incompressible viscous flow is developed and validated. The finite-volume Navier–Stokes solver is based on the artificial compressibility approach with a high-resolution method of characteristics-based scheme for handling convection terms. The unsteady flow is calculated with a matrix-free implicit dual time stepping scheme. The parallelization of the multigrid solver is achieved by multigrid domain decomposition approach (MG-DD), using single program multiple data (SPMD) and multiple instruction multiple data (MIMD) programming paradigm. There are two parallelization strategies proposed in this work, first strategy is a one-level parallelization strategy using geometric domain decomposition technique alone, second strategy is a two-level parallelization strategy that consists of a hybrid of both geometric domain decomposition and data decomposition techniques. Message-passing interface (MPI) and OpenMP standard are used to communicate data between processors and decompose loop iterations arrays, respectively. The parallel-multigrid code is used to simulate both steady and unsteady incompressible viscous flows over a circular cylinder and a lid-driven cavity flow. A maximum speedup of 22.5 could be achieved on 32 processors, for instance, the lid-driven cavity flow of Re = 1000. The results obtained agree well with numerical solutions obtained by other researchers as well as experimental measurements. A detailed study of the time step size and number of pseudo-sub-iterations per time step required for simulating unsteady flow are presented in this paper.     

Numerical simulation of 3D fluid–structure interaction flow using an immersed object method with overlapping grids

The newly developed immersed object method (IOM) [Tai CH, Zhao Y, Liew KM. Parallel computation of unsteady incompressible viscous flows around moving rigid bodies using an immersed object method with overlapping grids. J Comput Phys 2005; 207(1): 151–72] is extended for 3D unsteady flow simulation with fluid–structure interaction (FSI), which is made possible by combining it with a parallel unstructured multigrid Navier–Stokes solver using a matrix-free implicit dual time stepping and finite volume method [Tai CH, Zhao Y, Liew KM. Parallel computation of unsteady three-dimensional incompressible viscous flow using an unstructured multigrid method. In: The second M.I.T. conference on computational fluid and solid mechanics, June 17–20, MIT, Cambridge, MA 02139, USA, 2003; Tai CH, Zhao Y, Liew KM. Parallel computation of unsteady three-dimensional incompressible viscous flow using an unstructured multigrid method, Special issue on “Preconditioning methods: algorithms, applications and software environments. Comput Struct 2004; 82(28): 2425–36]. This uniquely combined method is then employed to perform detailed study of 3D unsteady flows with complex FSI. In the IOM, a body force term F is introduced into the momentum equations during the artificial compressibility (AC) sub-iterations so that a desired velocity distribution V₀ can be obtained on and within the object boundary, which needs not coincide with the grid, by adopting the direct forcing method. An object mesh is immersed into the flow domain to define the boundary of the object. The advantage of this is that bodies of almost arbitrary shapes can be added without grid restructuring, a procedure which is often time-consuming and computationally expensive. It has enabled us to perform complex and detailed 3D unsteady blood flow and blood–leaflets interaction in a mechanical heart valve (MHV) under physiological conditions.     

The simulation of 3D unsteady incompressible flows with moving boundaries on unstructured meshes

The development of a computational model for the simulation of three-dimensional unsteady incompressible viscous fluid flows with moving boundaries is presented. The numerical model is based upon the solution of the Navier–Stokes equations on unstructured meshes using the artificial compressibility approach. An ALE formulation is adopted and the equations are discretized using a cell vertex finite volume method. The formulation ensures the satisfaction of the geometric conservation law when the mesh is allowed to move. An implicit time discretization is adopted and a dual time approach is employed. Explicit relaxation is used for the sub-iterations, with multigrid acceleration. For moving geometries, the mesh is deformed by adopting a spring analogy, combined with a wall distance function approach. The numerical procedure is validated on a standard problem and is then used for the simulation of flow over a flexible fish-like body.     

Numerical simulation of incompressible flows using adaptive unstructured meshes and the pseudo-compressibility hypothesis

The numerical simulation of incompressible viscous flows, using finite elements with automatic adaptive unstructured meshes and the pseudo-compressibility hypothesis, is presented in this work. Special emphasis is given to the automatic adaptive process of unstructured meshes with linear tetrahedral elements in order to get more accurate solutions at relatively low computational costs. The behaviour of the numerical solution is analyzed using error indicators to detect regions where some important physical phenomena occur. An adaptive scheme, consisting in a mesh refinement process followed by a nodal re-allocation technique, is applied to the regions in order to improve the quality of the numerical solution. The error indicators, the refinement and nodal re-allocation processes as well as the corresponding data structure (to manage the connectivity among the different entities of a mesh, such as elements, faces, edges and nodes) are described. Then, the formulation and application of a mesh adaptation strategy, which includes a refinement scheme, a mesh smoothing technique, very simple error indicators and an adaptation criterion based in statistical theory, integrated with an algorithm to simulate complex two and three dimensional incompressible viscous flows, are the main contributions of this work. Two numerical examples are presented and their results are compared with those obtained by other authors.     

A p-multigrid spectral difference method for two-dimensional unsteady incompressible Navier-Stokes equations

This paper presents the development of a 2D high-order solver with spectral difference method for unsteady incompressible Navier-Stokes equations accelerated by a p-multigrid method. This solver is designed for unstructured quadrilateral elements. Time-marching methods cannot be applied directly to incompressible flows because the governing equations are not hyperbolic. An artificial compressibility method (ACM) is employed in order to treat the inviscid fluxes using the traditional characteristics-based schemes. The viscous fluxes are computed using the averaging approach (Sun et al., 2007; Kopriva, 1998) [29] and [12]. A dual time stepping scheme is implemented to deal with physical time marching. A p-multigrid method is implemented (Liang et al., 2009) [16] in conjunction with the dual time stepping method for convergence acceleration. The incompressible SD (ISD) method added with the ACM (SD-ACM) is able to accurately simulate 2D steady and unsteady viscous flows.     

Computational technologies for simulation of complex near-wall turbulent flows using unstructured meshes

The work is devoted to the peculiarities of the implementation of hybrid Reynolds’ Averaged Navier-Stokes equations-Large Eddy Simulation (RANS-LES) approaches of the Detached Eddy Simulation (DES) family for simulation of complex near-wall turbulent flows using unstructured meshes. The problems of determining required geometric characteristics in the mesh nodes and adaptation of hybrid approaches to the used accurate numerical approximation scheme in space are considered. The classic benchmark problem of the decay of homogeneous isotropic turbulence and the results of the computation of a complex turbulent flow near the wall with the presence of flow separation and reattachment are considered to verify the implemented technique and to demonstrate its efficiency.     

A block LU-SGS implicit unsteady incompressible flow solver on hybrid dynamic grids for 2D external bio-fluid simulations

A hybrid dynamic grid generation technique for two-dimensional (2D) morphing bodies and a block lower-upper symmetric Gauss-Seidel (BLU-SGS) implicit dual-time-stepping method for unsteady incompressible flows are presented for external bio-fluid simulations. To discretize the complicated computational domain around 2D morphing configurations such as fishes and insect/bird wings, the initial grids are generated by a hybrid grid strategy firstly. Body-fitted quadrilateral (quad) grids are generated first near solid bodies. An adaptive Cartesian mesh is then generated to cover the entire computational domain. Cartesian cells which overlap the quad grids are removed from the computational domain, and a gap is produced between the quad grids and the adaptive Cartesian grid. Finally triangular grids are used to fill this gap. During the unsteady movement of morphing bodies, the dynamic grids are generated by a coupling strategy of the interpolation method based on ‘Delaunay graph’ and local remeshing technique. With the motion of moving/morphing bodies, the grids are deformed according to the motion of morphing body boundaries firstly with the interpolation strategy based on ‘Delaunay graph’ proposed by Liu and Qin. Then the quality of deformed grids is checked. If the grids become too skewed, or even intersect each other, the grids are regenerated locally. After the local remeshing, the flow solution is interpolated from the old to the new grid. Based on the hybrid dynamic grid technique, an efficient implicit finite volume solver is set up also to solve the unsteady incompressible flows for external bio-fluid dynamics. The fully implicit equation is solved using a dual-time-stepping approach, coupling with the artificial compressibility method (ACM) for incompressible flows. In order to accelerate the convergence history in each sub-iteration, a block lower-upper symmetric Gauss-Seidel implicit method is introduced also into the solver. The hybrid dynamic grid generator is tested by a group of cases of morphing bodies, while the implicit unsteady solver is validated by typical unsteady incompressible flow case, and the results demonstrate the accuracy and efficiency of present solver. Finally, some applications for fish swimming and insect wing flapping are carried out to demonstrate the ability for 2D external bio-fluid simulations.     

Performance and Comparison of Cell-Centered and Node-Centered Unstructured Finite Volume Discretizations for Shallow Water Free Surface Flows

Finite volume (FV) methods for solving the two-dimensional (2D) nonlinear shallow water equations (NSWE) with source terms on unstructured, mostly triangular, meshes are known for some time now. There are mainly two basic formulations of the FV method: node-centered (NCFV) and cell-centered (CCFV). In the NCFV formulation the finite volumes, used to satisfy the integral form of the equations, are elements of the mesh dual to the computational mesh, while for the CCFV approach the finite volumes are the mesh elements themselves. For both formulations, details are given of the development and application of a second-order well-balanced Godunov-type scheme, developed for the simulation of unsteady 2D flows over arbitrary topography with wetting and drying. The popular approximate Riemann solver of Roe is utilized to compute the numerical fluxes, while second-order spatial accuracy is achieved with a MUSCL-type reconstruction technique. The Green-Gauss (G-G) formulation for gradient computations is implemented for both formulations, in order to maintain a common framework. Two different stencils for the G-G gradient computations in the CCFV formulation are implemented and tested. An edge-based limiting procedure is applied for the control of the total variation of the reconstructed field. This limiting procedure is proved to be effective for the NCFV scheme but inadequate for the CCFV approach. As such, a simple but very effective modification to the reconstruction procedure is introduced that takes into account geometrical characteristics of the computational mesh. In addition, consistent well-balanced second-order discretizations for the topography source term treatment and the wet/dry front treatment are presented for both FV formulations, ensuring absolute mass conservation, along with a stable friction term treatment.     

A virtual test facility for the efficient simulation of solid material response under strong shock and detonation wave loading

A virtual test facility (VTF) for studying the three-dimensional dynamic response of solid materials subject to strong shock and detonation waves has been constructed as part of the research program of the Center for Simulating the Dynamic Response of Materials at the California Institute of Technology. The compressible fluid flow is simulated with a Cartesian finite volume method and treating the solid as an embedded moving body, while a Lagrangian finite element scheme is employed to describe the structural response to the hydrodynamic pressure loading. A temporal splitting method is applied to update the position and velocity of the boundary between time steps. The boundary is represented implicitly in the fluid solver with a level set function that is constructed on-the-fly from the unstructured solid surface mesh. Block-structured mesh adaptation with time step refinement in the fluid allows for the efficient consideration of disparate fluid and solid time scales. We detail the design of the employed object-oriented mesh refinement framework AMROC and outline its effective extension for fluid–structure interaction problems. Further, we describe the parallelization of the most important algorithmic components for distributed memory machines and discuss the applied partitioning strategies. As computational examples for typical VTF applications, we present the dynamic deformation of a tantalum cylinder due to the detonation of an interior solid explosive and the impact of an explosion-induced shock wave on a multi-material soft tissue body.     

Parallel-adaptive simulation with the multigrid-based software framework UG

In this paper we present design aspects and concepts of the unstructured grids (UG) software framework that are relevant for parallel-adaptive simulation of time-dependent, nonlinear partial differential equations. The architectural design is discussed on system, subsystem and component level for distributed mesh management and local adaptation capabilities. Parallelization is founded on top of the innovative programming model dynamic distributed data (DDD). Newly introduced modules and extensions of DDD are discussed. Local multigrid methods are introduced as optimal linear solvers in the solution process. The demands of local parallel mesh adaptation are further described: Beside a mesh manipulation module further steps dynamic load balancing and migration have to be introduced. Their realization in the context of local multigrid methods is significantly non-trivial and makes the major contribution to the paper presented here. Parallel I/O provides an efficient mechanism for restart, postprocessing and long-term, large-scale computations. The UG approach is verified through a considerable code-reuse fraction of nearly 90% for simulations of complicated phenomena like porous media flow and transport as well as elastoplasticity. Parallel simulations with up to 10⁸ unknowns are shown for the Couplex benchmark. Therefore a grid convergence study to verify the reliability of the computed results is possible. For an parallel-adaptive elastoplasticity computation the speedup of the multigrid solver, which is the most scalability critical simulation part, exceeds on 512 processor a value of 300. The overhead introduced by the parallel-adaptive scheme turns out to be below 10% of the whole simulation time.     

Application of single-level,pointwise algebraic,and smoothed aggregation multigrid methods to direct numerical simulations of incompressible turbulent flows

Single- and multi-level iterative methods for sparse linear systems are applied to unsteady flow simulations via implementation into a direct numerical simulation solver for incompressible turbulent flows on unstructured meshes. The performance of these solution methods, implemented in the well-established SAMG and ML packages, are quantified in terms of computational speed and memory consumption, with a direct sparse LU solver (SuperLU) used as a reference. The classical test case of unsteady flow over a circular cylinder at low Reynolds numbers is considered, employing a series of increasingly fine anisotropic meshes. As expected, the memory consumption increases dramatically with the considered problem size for the direct solver. Surprisingly, however, the computation times remain reasonable. The speed and memory usage of pointwise algebraic and smoothed aggregation multigrid solvers are found to exhibit near-linear scaling. As an alternative to multi-level solvers, a single-level ILUT-preconditioned GMRES solver with low drop tolerance is also considered. This solver is found to perform sufficiently well only on small meshes. Even then, it is outperformed by pointwise algebraic multigrid on all counts. Finally, the effectiveness of pointwise algebraic multigrid is illustrated by considering a large three-dimensional direct numerical simulation case using a novel parallelization approach on a large distributed memory computing cluster.     

Numerical simulation of flow-induced cavity noise in self-sustained oscillations

The generation and near-field radiation of aerodynamic sound from a low-speed unsteady flow over a two-dimensional automobile door cavity is simulated by using a source-extraction-based coupling method. In the coupling procedure, the unsteady cavity flow field is first computed solving the Reynolds- averaged Navier–Stokes (RANS) equations. The radiated sound is then calculated by using a set of acoustic perturbation equations with acoustic source terms which are extracted from the time-dependent solutions of the unsteady flow. The aerodynamic and its resulting acoustic field are computed for the Reynolds number of 53,266 based on the base length of the cavity. The free stream flow velocity is taken to be 50.9 m/s. As first stage of the numerical investigation of flow-induced cavity noise, laminar flow is assumed. The CFD solver is based on a cell-centered finite volume method. A dispersion-relation-preserving (DRP), optimized, fourth-order finite difference scheme with fully staggered-grid implementation is used in the acoustic solver.     

Local penalty methods for flows interacting with moving solids at high Reynolds numbers

An original numerical modelling of multiphase flows interacting with solids in unsteady regimes is presented. Based on the generalized Navier-Stokes equations for multiphase flows and Volume of Fluid (VOF) formulations, an Uzawa minimization algorithm is implemented for the treatment of incompressibility and solid constraints. Augmented Lagrangian terms are added in the momentum equations to speed the convergence of the iterative solver. Defining a priori the penalty parameters which are dedicated to incompressibility and solid constraints is difficult, or impossible, as soon as the flow involves more than one phase and inertia becomes predominant compared to viscous and gravity forces. The Lagrangian penalty terms are calculated automatically according to an original local estimate of the various physical contributions. Numerical validations have been carried out for single particle settling in confined media and viscous flow through a fixed Cubic Faced Centered array. A very good agreement is obtained between experimental, theoretical and numerical results. Extension to unsteady free surface flow interacting with particles is illustrated with the simulation of a dam break flow over moving obstacles.     

A Parallel Domain Decomposition Method for 3D Unsteady Incompressible Flows at High Reynolds Number

Numerical simulation of three-dimensional incompressible flows at high Reynolds number using the unsteady Navier–Stokes equations is challenging. In order to obtain accurate simulations, very fine meshes are necessary, and such simulations are increasingly important for modern engineering practices, such as understanding the flow behavior around high speed trains, which is the target application of this research. To avoid the time step size constraint imposed by the CFL number and the fine spacial mesh size, we investigate some fully implicit methods, and focus on how to solve the large nonlinear system of equations at each time step on large scale parallel computers. In most of the existing implicit Navier–Stokes solvers, segregated velocity and pressure treatment is employed. In this paper, we focus on the Newton–Krylov–Schwarz method for solving the monolithic nonlinear system arising from the fully coupled finite element discretization of the Navier–Stokes equations on unstructured meshes. In the subdomain, LU or point-block ILU is used as the local solver. We test the algorithm for some three-dimensional complex unsteady flows, including flows passing a high speed train, on a supercomputer with thousands of processors. Numerical experiments show that the algorithm has superlinear scalability with over three thousand processors for problems with tens of millions of unknowns.     

Parallel preconditioned WENO scheme for three-dimensional flow simulation of NREL Phase VI Rotor

The preconditioned weighted essentially non-oscillatory (P-WENO) solver for viscous flows (Huang et al. (2009) [9]) is extended to non-inertial reference frames. In the present scheme, patched multi-block grid system is employed and parallel computing is adopted as well. With the present parallel P-WENO solver, three-dimensional flows of the Phase VI Rotor from National Renewable Energy Laboratory (NREL) can be simulated and analyzed for different wind speeds. Our simulation results show good agreement with the numerical predictions based on incompressible Navier–Stokes (N–S) equations as well as the available wind tunnel data from NREL. The flow phenomena, including separation and attachment line, can be captured by the present scheme. The parallel strategy adopted is a block-domain decomposition method for the patched multi-block grid system. To balance the load among different computing nodes, a Tabu search algorithm is adopted for the parallelization. The parallel efficiency of the parallel P-WENO scheme is examined for node numbers ranging from 1 to 64. It is found that the parallel efficiency is monotonically decreased as the node number adopted is increased; the parallel efficiency is retained over 90% for all cases of different node numbers. Due to the high parallel efficiency, our parallel P-WENO solver is validated for applying to practical fluid problems from compressible to incompressible limits.     

Applications of stencil-adaptive finite difference method to incompressible viscous flows with curved boundary

This paper investigates the applicability of the stencil-adaptive finite difference method for the simulation of two-dimensional unsteady incompressible viscous flows with curved boundary. The adaptive stencil refinement algorithm has been proven to be able to continuously adapt the stencil resolution according to the gradient of flow parameter of interest [Ding H, Shu C. A stencil adaptive algorithm for finite difference solution of incompressible viscous flows. J Comput Phys 2006;214:397-420], which facilitates the saving of the computational efforts. On the other hand, the capability of the domain-free discretization technique in dealing with the curved boundary provides a great flexibility for the finite difference scheme on the Cartesian grid. Here, we show that their combination makes it possible to simulate the unsteady incompressible flow with curved boundary on a dynamically changed grid. The methods are validated by simulating steady and unsteady incompressible viscous flows over a stationary circular cylinder.     

Error control and mesh optimization for high-order finite element approximation of incompressible viscous flow

We discuss the use of a posteriori error estimates for high-order finite element methods during simulation of the flow of incompressible viscous fluids. The correlation between the error estimator and actual error is used as a criterion for the error analysis efficiency. We show how to use the error estimator for mesh optimization which improves computational efficiency for both steady-state and unsteady flows. The method is applied to two-dimensional problems with known analytical solutions (Jeffrey-Hamel flow) and more complex flows around a body, both in a channel and in an open domain.