A spectrally accurate method for overlapping grid solution of incompressible Navier–Stokes equations
Introduction
Finding numerical solutions to partial differential equations (PDEs) by decomposing the computational domain into smaller subdomains is an idea that has been around for well over a century. Domain decomposition methods have been utilized for several different purposes, including straightforward parallelization [1], [2], [3], [4], simplified mesh generation for complex geometries [5], [6], [7], [8], and the ability to use different parameters or methods in different subdomains [9], [10], [11], [12]. These techniques exist in many forms, and each has its strengths. Some decompose the global domain into overlapping subdomains [13], [14], [15], [16], [17], [18], while others employ non-overlapping subdomains [19], [9], [7], [20], [21], [22]. Some use explicit interpolation techniques for values at interface boundaries [9], [23], [20], [24], and others carry out implicit interpolation [25], [26], [22], [27], [13], [14], [28], [29], [30]. Domain decomposition techniques have been developed for use with several numerical methods including finite difference [14], [5], [9], [23], finite element [25], finite volume [17], [31], and spectral methods [21], [22].
The earliest known research in domain decomposition methods was performed by H.A. Schwarz whose work was published in 1870 [13]. The original Schwarz Alternating Method, initially proposed for analytical calculations [32], was developed for the global solution of boundary value problems for harmonic functions [15] decomposed into overlapping subdomains, . The solution in the first subdomain ( with boundaries and ) is found using the global boundary conditions on and corresponding values from at the previous iteration on . The solution of is then found by using values from the solution in on . These two steps are iterated until sufficient convergence is reached (see [33], [32], [34]).
In the 1960's, Volkov generalized the original Schwarz Alternating Method into a numerical domain decomposition technique, in a form of the Composite Mesh Difference Method (CMDM) [14]. CMDM used finite difference methods to solve the 2-dimensional Poisson equation numerically on overlapping grids. His research laid the foundation for subsequent techniques that extended the use of CMDM to other elliptical and hyperbolic PDEs, and boundary value and initial value problems, with the ability to use curvilinear meshes (see [35], [36], [37], [38], [6]). Overlapping domain decomposition methods have also been developed to model complex equations and handle various difficulties in solving practical problems. The Chimera Grid Scheme, introduced in [39], employs overset (overlapping) grids for modeling flows in complex geometries. Shortly after initial development it was enhanced for use with three dimensional flows modeled by the Euler equations [40] and later with the addition of the thin-layer Navier–Stokes equations. More recently, Chimera Grid techniques have been used to model a variety of problems with complex geometries [41], [42]. Subsequently, Henshaw and Schwendeman [17], [4] developed a method for using overlapping mesh techniques in modeling high-speed reactive flows, in two and three dimensions.
In addition, techniques that employ non-overlapping grids (sometimes called patched grids) were developed. Examples include a zonal approach that uses a flux-vector splitting technique for the determination of interface values in Euler equations [43], [44], [45], [46], Lions method [19] that uses an iterative technique to arrive at the correct values to be passed between non-overlapping subdomains in solving Laplace's equation and more general second-order elliptic problems, Dawson's approach [9] that solves the heat equation using an explicit finite difference formula to determine the interface values and allows for different time stepping to be used in different subdomains. Non-overlapping grid techniques have also been extended and employed in solving the advection–diffusion equation [20] and the Navier–Stokes equations [45], [22]. Some of the more recently developed non-overlapping domain decomposition methods achieve high finite global accuracy [47] and some spectral accuracy [21], [7], [48], [22].
While non-overlapping mesh methods allow some flexibility in mesh generation, the constraints in these techniques inhibit additional flexibility that is seen in overlapping mesh methods. By allowing variable overlap size, a broad range of potential mesh configurations are supported with overlapping methods, thus allowing for more simplified mesh generation. Additionally, overlapping methods provide a convenient framework for further extension towards moving domain methods, allowing for general and unconstrained motion of rigid body parts through the background stationary meshes [49], [50], [51], [52].
So far, existing overlapping grid methods for the time-dependent PDE coupling have been traditionally relying on low-order, finite-difference or finite-volume schemes. Although some of the methods have been extended to achieve higher-order spatial convergence, using extended stencil finite-difference or compact schemes, the upper bound of the global accuracy has been usually limited to four [23], [24], [53], and at most six [54], [55], [56], [27], [57]. Recently, Brazell, Sitaraman and Mavriplis developed a high-order overlapping Discontinuous Galerkin solver for compressible equations, that uses Lagrangian interpolation at interface boundaries, and documented a polynomial convergence up to fourth order [58]. In the current paper, we introduce a spectrally-accurate overlapping mesh method for incompressible equations, that is based on a spectral-element method. The Spectral Element Method, which can be perceived as a high-order extension of the Finite Element Method, divides a domain into several conforming and adjacent elements [59], [3], [60]. The volume within each element is discretized using Nth-order tensor-product Lagrange interpolating polynomials on Gauss–Lobatto–Legendre nodal points. Approximations in all elements are coupled at the boundaries to form a global solution [60], which achieves spectral convergence with p- (polynomial order) refinement. In this work, we combine a Spectral Element Method solver with the overlapping grid approach, to arrive at a globally spectrally-accurate method for solution of the incompressible Navier–Stokes equations on overlapping domains.
One of the inherent challenges with overlapping grid methods is to minimize the errors that are introduced due to the coupling of the individual subdomain solutions into the global solution. The coupling errors consist of spatial errors and temporal errors, and have to be treated separately. Spatial errors are introduced by the spatial interpolation stencil employed to obtain a function value at the interface points of one domain from the gridpoint values in the adjacent domains at the same time level. Some overlapping mesh methods circumvent the spatial error by requiring that the gridpoints in overlapping domains exactly coincide [61], [62], thus fully conserving communicated information, with the drawback of decreased flexibility in mesh generation. Other methods that do not require the exact match of the gridpoints and thus are more flexible, use finite order interpolation schemes to determine values from adjacent subdomains. Although simple linear interpolation techniques have been popular [6], [9], [63], [29], [64], it was shown by Chesshire and Henshaw [6] that an interpolation scheme should be consistent with the accuracy of the underlying solver and higher-order interpolation is required to maintain the accuracy of the coupled solution with high-order methods, generally, of the order of one higher than the underlying solver accuracy if the overlap width scales with the grid resolution, and the same if it stays constant during grid refinement. Thus, in fourth- and sixth-order methods [58], [23], [27], [24], a generalized Lagrangian interpolation method of the corresponding order was employed. In the present work, we use Nth-order Lagrangian interpolation approach, on non-uniform Gauss–Lobatto–Legendre grids, to arrive at a spectrally-accurate interpolation scheme (with p-refinement) that allows us to maintain the global spectral accuracy of the coupled solution.
To avoid temporal errors due to coupling, the subdomain coupling methods have often been formulated implicitly, when either the interpolation dependencies are written into the global matrix [6], [26], [27], or independent subdomain solutions are coupled through the use of Schwarz-like iterations [13] to ensure interface value convergence [28], [29], [30], [31]. Global convergence of the variables in adjacent subdomains often requires many, sometimes hundreds [29], or even thousands [31] iterations, the number generally being dependent upon the PDE, the overlap size, and timestep [65], [66]. Although Henshaw [67] showed that the convergence rates when solving elliptic PDEs on overlapping subdomains can be reduced by using a multigrid solver with a smoothing algorithm near interface boundaries, this approach is difficult to implement and requires some global interventions to the solver. However, if the coupling scheme, even without iterations, does not introduce temporal errors larger than those of the global timestepping scheme, the global temporal accuracy will be preserved [68], particularly the interface values need not strictly match, but must merely be “consistent”. The notion of consistency has been previously used when constructing and analyzing interpolation and iteration schemes for the overlapping grid methods [30], where consistency was interpreted in a spatial sense. In [68], we have proposed an explicit temporal interface extrapolation algorithm for the overlapping grid methods that preserves the overall temporal accuracy, and analyzed its stability. The proposed coupling scheme is essentially “temporally consistent” (the differences between the interface values due to temporal coupling converge with the same rate as the temporal errors of the underlying flow solver).
In the present work, we implement this temporally consistent coupling scheme for coupling the solutions of incompressible Navier–Stokes equations on overlapping grids discretized with the SEM method. This novel scheme achieves a specified order of accuracy with a small number of iterations. It ensures that the temporal error between the values in the subdomains at the interfaces is equal to the temporal accuracy of the underlying SEM solver. Assuming that this holds true, the interface values need not be an exact match, while still maintaining the global accuracy of the solution. While the accuracy of the extrapolation scheme is maintained in the absence of iterations, the stability of the formulation can generally be improved by implementing a low number of Schwarz-like iterations [68].
The current paper introduces an original method for coupling Navier–Stokes equations on the overlapping domains that is spectrally accurate in space and high-order accurate in time up to a specified order. Although an overlapping Schwarz method has been previously employed as a preconditioner for the solution of the Poisson equation in a spectral element formulation of Navier–Stokes equations on a single domain [3], [69], a spectrally-accurate approach that achieves a full unsteady coupling of independent solutions to Navier–Stokes equations on overlapping grids appears to be novel. In addition, our method achieves a specified temporal accuracy with the minimum number of iterations (usually two to three) while maintaining stability of the formulation. We demonstrate the spectral spatial and high-order temporal convergence of the developed solver, as well as the overall consistency of the coupled solution, on a number of two and three-dimensional examples with incompressible Navier–Stokes equations, including laminar and turbulent test cases.
The remainder of the paper is organized in the following manner: section 2 describes the computational methodology; section 3 discusses the techniques used to parallelize the overlapping grid solver; section 4 presents the results for verification and validation; and section 5 shows the scalability tests. In section 6 some final conclusions are drawn, including a brief discussion of future extensions and applications of this work.
Section snippets
Mathematical formulation
In the present formulation, a two or three dimensional global solution domain is decomposed into two overlapping subdomains, , an example of which is seen in Fig. 1. Values at interface boundaries, are equal to contiguous values in the adjacent subdomain, and the solution in each subdomain is governed by the incompressible Navier–Stokes equations which are represented in non-dimensional form in space , as shown below and
Parallelization
Communication among separate subdomains is currently accomplished through a dual-session communication framework, which allows for two sessions to perform computations separately while enabling real-time communication with each other in a single computational platform. Each session runs its own independent instance of the SEM solver for its subdomain where MPI intracommunicators are established for communication among the local processors. An intercommunicator is established for global
Results
In this section, we demonstrate performance of the developed overlapping grid method, on two- and three-dimensional laminar and turbulent problems, with the focus on spatial and temporal orders of convergence, accuracy in the presence of outflow disturbances and long integration times, and the ability of the method to reproduce correct turbulence characteristics on overlapping domains.
Scalability
In this section, we perform the scalability tests for the developed overlapping grid method. We utilize the previous test case of the turbulent pipe flow simulations as a platform for the scalability studies. However, the refined computational meshes were constructed, with higher element counts than in the DNS simulations described above (see label to Fig. 29, Fig. 30) in order to enable testing on high processors counts (up to 1024). The refined outside mesh consisted of three times the number
Conclusions
We have developed an overlapping grid methodology for two and three dimensional applications in a Spectral Element Method incompressible flow solver. Spectral accuracy is maintained for spatial discretization due to a spectral interpolation at interface boundaries using solution approximations in an Nth-order polynomial space on GL points. The global high-order (up to a third) temporal accuracy of the flow solver is also maintained with few, or no, iterations using solutions at previous
Acknowledgements
This work has supported by the NSF grants #1250124 and #1025359. Support by the Office of Advanced Scientific Computing Research, Office of Science, U.S. Department of Energy, under Contract DE-AC02-06CH11357 and by the NSF RTG grant at Northwestern University is also acknowledged. We thank the NSF XSEDE program for providing the allocations on the SDSC Gordon and TACC Stampede clusters where the scalability and timing tests were performed.
References (86)
An overlapping Schwarz method for spectral element solution of the incompressible Navier–Stokes equations
J. Comput. Phys.
(1997)- et al.
Parallel computation of three-dimensional flows using overlapping grids with adaptive mesh refinement
J. Comput. Phys.
(2008) - et al.
On the use of composite grid schemes in computational aerodynamics
Comput. Methods Appl. Mech. Eng.
(1987) - et al.
Composite overlapping meshes for the solution of partial differential equations
J. Comput. Phys.
(1990) A fully conservative interface algorithm for overlapped grids
J. Comput. Phys.
(1995)- et al.
An adaptive numerical scheme for high-speed reactive flow on overlapping grids
J. Comput. Phys.
(2003) - et al.
A stable and conservative interface treatment of arbitrary spatial accuracy
J. Comput. Phys.
(1999) - et al.
Preconditioned spectral multi-domain discretization of the incompressible Navier–Stokes equations
J. Comput. Phys.
(2004) A fourth-order accurate method for the incompressible Navier–Stokes equations on overlapping grids
J. Comput. Phys.
(1994)- et al.
On the use of a high order overlapping grid method for coupling in CFD/CAA
J. Comput. Phys.
(2006)
Analysis of a fractional-step method on overset grids
J. Comput. Phys.
High-order compact finite-difference methods on general overset grids
J. Comput. Phys.
Adaptive local overlapping grid methods for parabolic systems in two space dimensions
J. Comput. Phys.
An overset-grid method for 3D unsteady incompressible flows
J. Comput. Phys.
Metric-discontinuous zonal grid calculations using the Osher scheme
Comput. Fluids
Applications of a conservative zonal scheme to transient and geometrically complex problems
Comput. Fluids
A staggered-grid multidomain spectral method for the compressible Navier–Stokes equations
J. Comput. Phys.
Parallel computing of overset grids for aerodynamic problems with moving objects
Prog. Aerosp. Sci.
A Chimera method based on a Dirichlet/Neumann (Robin) coupling for the Navier–Stokes equations
Comput. Methods Appl. Mech. Eng.
Moving overlapping grids with adaptive mesh refinement for high-speed reactive and non-reactive flow
J. Comput. Phys.
A wavenumber based extrapolation and interpolation method for use in conjunction with high-order finite difference schemes
J. Comput. Phys.
High order interpolation and differentiation using B-splines
J. Comput. Phys.
A spectral element method for fluid dynamics: laminar flow in a channel expansion
J. Comput. Phys.
Multiplicative Schwarz methods for parabolic problems
Appl. Math. Comput.
Chebyshev collocation method and multi-domain decomposition for Navier–Stokes equations in complex curved geometries
J. Comput. Phys.
Domain decomposition algorithms
Acta Numer.
Domain Decomposition: Parallel Multilevel Methods for Elliptic Partial Differential Equations
Domain decomposition by the mortar element method
Composite overset structured grids
A finite difference domain decomposition algorithm for numerical solution of the heat equation
Math. Comput.
A framework for coupling Reynolds-averaged with large-eddy simulations for gas turbine applications
J. Fluids Eng.
Coupled high-fidelity URANS simulation for helicopter applications
Center for Turbulence Research, Annual Research Briefs
Computational framework for coupling compressible and low Mach number codes
AIAA J.
Ueber einige Abbildungsaufgaben
J. Reine Angew. Math.
The method of composite meshes for finite and infinite regions with piecewise smooth boundary
Tr. Mat. Inst. Steklova
On the Schwarz alternating method. I
Heat transfer LES simulations in application to wire-wrapped fuel pins
On the Schwarz alternating method. III: a variant for nonoverlapping subdomains
Nonconforming Mortar Element Methods: Application to Spectral Discretizations
A polylogarithmic bound for an iterative substructuring method for spectral elements in three dimensions
SIAM J. Numer. Anal.
A domain decomposition method for incompressible viscous flow
SIAM J. Sci. Comput.
A composite multigrid method for calculating unsteady incompressible flows in geometrically complex domains
Int. J. Numer. Methods Fluids
Domain Decomposition: Parallel Multilevel Methods for Elliptic Partial Differential Equations
Cited by (39)
Cost vs Accuracy: DNS of turbulent flow over a sphere using structured immersed-boundary, unstructured finite-volume, and spectral-element methods
2023, European Journal of Mechanics, B/FluidsShock capturing using discontinuous Galerkin method and overset grids for two-dimensional Euler equations
2023, Journal of Computational PhysicsOn the theoretical foundation of overset grid methods for hyperbolic problems II: Entropy bounded formulations for nonlinear conservation laws
2022, Journal of Computational PhysicsCitation Excerpt :Overset grids have been used with all the major approximation methods, finite difference, finite element, finite volume and spectral element [20],[10],[8],[3],[12],[4],[13],[16].
On the theoretical foundation of overset grid methods for hyperbolic problems: Well-posedness and conservation
2022, Journal of Computational PhysicsVerification and convergence study of a spectral-element numerical methodology for fluid-structure interaction
2021, Journal of Computational Physics: XOverset meshes for incompressible flows: On preserving accuracy of underlying discretizations
2021, Journal of Computational PhysicsCitation Excerpt :This transfer of solution can happen via: 1) OSS [14,31,45,66], wherein interpolation equations for fringe points constitute constraint equations in the global system of equations. 2) OS [22,33,50,56,67,68], wherein the linear system for each mesh partition is solved independently and the solution at fringe points is updated in an iterative manner using interpolation equations. The iterations can be Gauss-Seidel-like (alternating Schwarz) or Jacobi-like (additive Schwarz).