Geometric Multigrid Method Research Articles

An adaptive finite volume solver for steady Euler equations with non-oscillatory k-exact reconstruction

In this paper, we present an adaptive finite volume method for steady Euler equations with a non-oscillatory k-exact reconstruction on unstructured mesh. The numerical framework includes a Newton method as an outer iteration to linearize the Euler equations, and a geometrical multigrid method as an inner iteration to solve the derived linear system. A non-oscillatory k-exact reconstruction of the conservative solution in each element is proposed for the high order and non-oscillatory behavior of the numerical solutions. The importance on handling the curved boundary in an appropriate way is also studied with the numerical experiments. The h-adaptive method is introduced to enhance the efficiency of the algorithm. The numerical tests show successfully that the quality solutions can be obtained smoothly with the proposed algorithm, i.e., the expected convergence order of the numerical solution with the mesh refinement can be reached, while the non-oscillation shock structure can be obtained. Furthermore, the mesh adaptive method with the appropriate error indicators can effectively enhance the implementation efficiency of numerical method, while the steady state convergence and numerical accuracy are kept in the meantime.

Geometric multigrid method for steady Buoyancy convection in vertical cylinders

An implementation of the Geometric Multi-Grid (GMG) method, with full approximation scheme/storage (FAS) algorithm, in a numerical study of steady Buoyancy-driven convection with axis-symmetric flows in vertical cylinders is presented. The particular system examined is a cylinder with an aspect ratio a (radius divided by height) of 4. The fluid is heated from below and cooled from above, and the circular wall of the vessel is insulated. The Rayleigh (Ra) and fluid Prandtl (Pr) numbers are respectively 5 6.4 10  and 7. A non-linear system of equations was formulated in stream function-vorticity-temperature variables and discretized using a monotonic conservative finite difference scheme of second order accuracy. The steady state condition was solved for purposes of comparing two numerical methods: the GMG FAS and the GaussSeidel method with lexicographic ordering (GS-LEX). The GMG FAS method has significantly higher efficiency in CPU performance compared to the pure GS-LEX method for fine grids only if the tolerance value for stopping iteration process is chosen not too small. A procedure for the selection of adjustment parameters for GMG FAS algorithm is proposed and tested for different grid sizes. Details regarding convergence criteria are addressed.

Parallel geometric multigrid

The paper describes practical approach for minimising the parallelisation overhead for a solution to the boundary value problems by geometric multigrid methods. It is shown that proposed multiple coarse grid correction strategy makes it possible not only to create the task of the smoother least demanding, but also to avoid load imbalance and limit the communication overhead. Estimation of maximum speedup and efficiency of parallel robust multigrid technique and parallel V-cycle are given.

Resilience for Massively Parallel Multigrid Solvers

Fault tolerant massively parallel multigrid methods for elliptic partial differential equations are a step towards resilient solvers. Here, we combine domain partitioning with geometric multigrid methods to obtain fast and fault-robust solvers for three-dimensional problems. The recovery strategy is based on the redundant storage of ghost values, as they are commonly used in distributed memory parallel programs. In the case of a fault, the redundant interface values can be easily recovered, while the lost inner unknowns are recomputed approximately with recovery algorithms using multigrid cycles for solving a local Dirichlet problem. Different strategies are compared and evaluated with respect to performance, computational cost, and speedup. Especially effective are asynchronous strategies combining global solves with accelerated local recovery. By this, multiple faults can be fully compensated with respect to both the number of iterations and run-time. For illustration, we use a state-of-the-art petascale supercomputer to study failure scenarios when solving systems with up to $6 \cdot 10^{11}$ (0.6 trillion) unknowns.

A nearly optimal multigrid method for general unstructured grids

In this paper, we develop a multigrid method on unstructured shape-regular grids. For a general shape-regular unstructured grid of ${\cal O}(N)$ elements, we present a construction of an auxiliary coarse grid hierarchy on which a geometric multigrid method can be applied together with a smoothing on the original grid by using the auxiliary space preconditioning technique. Such a construction is realized by a cluster tree which can be obtained in ${\cal O}(N\log N)$ operations for a grid of $N$ elements. This tree structure in turn is used for the definition of the grid hierarchy from coarse to fine. For the constructed grid hierarchy we prove that the convergence rate of the multigrid preconditioned CG for an elliptic PDE is $1 - {\cal O}({1}/{\log N})$. Numerical experiments confirm the theoretical bounds and show that the total complexity is in ${\cal O}(N\log N)$.

Implementation of the Vanka-type multigrid solver for the finite element approximation of the Navier–Stokes equations on GPU

We present a complete GPU implementation of a geometric multigrid solver for the numerical solution of the Navier–Stokes equations for incompressible flow. The approximate solution is constructed on a two-dimensional unstructured triangular mesh. The problem is discretized by means of the mixed finite element method with semi-implicit timestepping. The linear saddle-point problem arising from the scheme is solved by the geometric multigrid method with a Vanka-type smoother. The parallel solver is based on the red–black coloring of the mesh triangles. We achieved a speed-up of 11 compared to a parallel (4 threads) code based on OpenMP and 19 compared to a sequential code.

Applications of MMPBSA to Membrane Proteins I: Efficient Numerical Solutions of Periodic Poisson-Boltzmann Equation.

Continuum solvent models have been widely used in biomolecular modeling applications. Recently much attention has been given to inclusion of implicit membranes into existing continuum Poisson-Boltzmann solvent models to extend their applications to membrane systems. Inclusion of an implicit membrane complicates numerical solutions of the underlining Poisson-Boltzmann equation due to the dielectric inhomogeneity on the boundary surfaces of a computation grid. This can be alleviated by the use of the periodic boundary condition, a common practice in electrostatic computations in particle simulations. The conjugate gradient and successive over-relaxation methods are relatively straightforward to be adapted to periodic calculations, but their convergence rates are quite low, limiting their applications to free energy simulations that require a large number of conformations to be processed. To accelerate convergence, the Incomplete Cholesky preconditioning and the geometric multigrid methods have been extended to incorporate periodicity for biomolecular applications. Impressive convergence behaviors were found as in the previous applications of these numerical methods to tested biomolecules and MMPBSA calculations.

A Semi-Implicit Version of the MPAS-Atmosphere Dynamical Core

Abstract An important question for atmospheric modeling is the viability of semi-implicit time integration schemes on massively parallel computing architectures. Semi-implicit schemes can provide increased stability and accuracy. However, they require the solution of an elliptic problem at each time step, creating concerns about their parallel efficiency and scalability. Here, a semi-implicit (SI) version of the Model for Prediction Across Scales (MPAS) is developed and compared with the original model version, which uses a split Runge–Kutta (SRK3) time integration scheme. The SI scheme is based on a quasi-Newton iteration toward a Crank–Nicolson scheme. Each Newton iteration requires the solution of a Helmholtz problem; here, the Helmholtz problem is derived, and its solution using a geometric multigrid method is described. On two standard test cases, a midlatitude baroclinic wave and a small-planet nonhydrostatic gravity wave, the SI and SRK3 versions produce almost identical results. On the baroclinic wave test, the SI version can use somewhat larger time steps (about 60%) than the SRK3 version before losing stability. The SI version costs 10%–20% more per step than the SRK3 version, and the weak and strong scalability characteristics of the two versions are very similar for the processor configurations the authors have been able to test (up to 1920 processors). Because of the spatial discretization of the pressure gradient in the lowest model layer, the SI version becomes unstable in the presence of realistic orography. Some further work will be needed to demonstrate the viability of the SI scheme in this case.

GPU accelerated VOF based multiphase flow solver and its application to sprays

A new volume of fluid (VOF) based two-phase flow solver has been developed that has been parallelized for general purpose graphical processor units (GPGPU). Incompressible Navier–Stokes equation is explicitly solved using finite volume method on a 2-D Cartesian grid. The pressure Poisson equation is solved on the GPU in conjunction with the geometric multigrid method. This solver performance is established by validating it against standard benchmark test cases. Furthermore, the solver acceleration is determined by comparing it with the serial version of the code. The solver was then used to study the primary jet breakup of a 2-D planar liquid jet due to a co-flowing gas stream. The effect of liquid to gas phase velocity ratio on jet breakup was investigated.

Solver preconditioning using the combinatorial multilevel method on reservoir simulation

The purpose of this paper is to report the first preliminary study of the recently introduced combinatorial multilevel (CML) method for solver preconditioning in large-scale reservoir simulation with coupled geomechanics. The CML method is a variant of the popular algebraic multigrid (AMG) method yet with essential differences. The basic idea of this new approach is to construct a hierarchy of matrices using the discrete geometry of the graph, based on support theory for preconditioners. In this way, the CML method combines the merits of both geometric and algebraic multigrid methods. The resulting hybrid approach not only provides a simpler and faster setup phase compared to AMG, but the method can be proven to exhibit strong convergence guarantees for arbitrary symmetric diagonally dominant matrices. In addition, the underlying theoretical soundness of the CML method contrasts to the heuristic AMG approach, which often can show slow convergence for difficult problems. This new approach is implemented in a reservoir simulator for both pressure and poroelastic displacement preconditioners in the multistage preconditioning technique. We present results based on several known benchmark problems and provide a comparison of performance and complexity with the widespread preconditioning schemes used in large-scale reservoir simulation. An adaptation of CML for non-symmetric matrices is shown to exhibit excellent convergence properties for realistic cases.

New Extrapolation Multiscale Multigrid Method for Second Order Elliptic Problem

Motivated by the idea of algebraic multigrid (AMG) method, we take the transpose of bilinear interpolation operator as restriction operator, an improved geometric V-cycle multigrid (IGVMG) method is proposed. By combining the IGVMG method and Chen's extrapolation formula, a new extrapolation multiscale multigrid (EXMMG) method for a second order elliptic problem is constructed. Then the IGVMG method is used to solve the solutions for coarse grid and finest grid, respectively. Chen's extrapolation formula is applied to update the approximate solution of coarse grid, and the cubic interpolation operator is applied to provide a better initial value on fine grid from coarse grid. Some numerical experiments of the EXMMG method are presented, using fourth order compact difference scheme. It is demonstrated that the EXMMG method is effective.

An analytical coarse grid operator applied to a multiscale multigrid method

Equations describing flow through heterogeneous porous media are often difficult to solve by direct numerical simulation given their expensive computational requirement. In this paper, numerical convergence is achieved for solving the flow equation with high contrast variable coefficients using two low computational methods that incorporate analytical approximations into the geometric multigrid. The methods are inheritably fast and easily implementable with coefficients describing general non-periodic media. The first method, named here as Analytical Coarse Operator (ACO), defines the operator at each coarse-scale of the V-cycle using an analytical approximation of the upscale tensor, it works in the same way as iterative homogenization. Its efficiency and reliability are measured by comparing with similar algorithms using arithmetic and harmonic averages of the coefficient, since they are widely used in the literature. The second method, called Analytical Prolongation Operator (APO), may be termed as aggressive coarsening, because it obtains the solution by skipping few levels in the V-cycle. It is an extension of the first approach where a multiscale prolongation operator is defined using an analytical approximation of the solution of the cell-problem from homogenization theory. The efficiency and reliability of this method are measured by comparing the number of V-cycles with the full numerical implementation. Various test cases demonstrated that convergence can be achieved using typical implementation procedures of the geometric multigrid method, therefore it highlights the low cost and easy implementation features of the methods. The examples include media having separable and non separable scales, periodic or not. In the case of non periodic medium, the problem of determining coarse grids that are representative element volume (REV) is also addressed. The medium from the SPE10 benchmark project is an application of the method for non separable scale.

Evaluation of multilevel sequentially semiseparable preconditioners on computational fluid dynamics benchmark problems using Incompressible Flow and Iterative Solver Software

This paper studies a new preconditioning technique for sparse systems arising from discretized partial differential equations in computational fluid dynamics problems. This preconditioning technique exploits the multilevel sequentially semiseparable (MSSS) structure of the system matrix. MSSS matrix computations give a data‐sparse way to approximate the LU factorization of a sparse matrix from discretized partial differential equations in linear computational complexity with respect to the problem size. In contrast to the standard block diagonal and block upper‐triangular preconditioners, we exploit the global MSSS structure of the 2×2 block system from the discretized Stokes equation and linearized Navier‐Stokes equation. This avoids approximating the Schur complement explicitly, which is a big advantage over standard block preconditioners. Through numerical experiments on standard computational fluid dynamics benchmark problems in Incompressible Flow and Iterative Solver Software, we show the performance of the MSSS preconditioners. They indicate that the global MSSS preconditioner not only yields mesh size independent convergence but also gives viscosity parameter and Reynolds number independent convergence. Compared with the algebraic multigrid (AMG) method and the geometric multigrid (GMG) method for block preconditioners, the MSSS preconditioning technique is more robust than both the AMG method and GMG method, and considerably faster than the AMG method. Copyright © 2015 John Wiley & Sons, Ltd.

A coupled ordinates method for solution acceleration of rarefied gas dynamics simulations

Non-equilibrium rarefied flows are frequently encountered in a wide range of applications, including atmospheric re-entry vehicles, vacuum technology, and microscale devices. Rarefied flows at the microscale can be effectively modeled using the ellipsoidal statistical Bhatnagar–Gross–Krook (ESBGK) form of the Boltzmann kinetic equation. Numerical solutions of these equations are often based on the finite volume method (FVM) in physical space and the discrete ordinates method in velocity space. However, existing solvers use a sequential solution procedure wherein the velocity distribution functions are implicitly coupled in physical space, but are solved sequentially in velocity space. This leads to explicit coupling of the distribution function values in velocity space and slows down convergence in systems with low Knudsen numbers. Furthermore, this also makes it difficult to solve multiscale problems or problems in which there is a large range of Knudsen numbers. In this paper, we extend the coupled ordinates method (COMET), previously developed to study participating radiative heat transfer, to solve the ESBGK equations. In this method, at each cell in the physical domain, distribution function values for all velocity ordinates are solved simultaneously. This coupled solution is used as a relaxation sweep in a geometric multigrid method in the spatial domain. Enhancements to COMET to account for the non-linearity of the ESBGK equations, as well as the coupled implementation of boundary conditions, are presented. The methodology works well with arbitrary convex polyhedral meshes, and is shown to give significantly faster solutions than the conventional sequential solution procedure. Acceleration factors of 5–9 are obtained for low to moderate Knudsen numbers on single processor platforms.

A Coupled Ordinates Method for Convergence Acceleration of the Phonon Boltzmann Transport Equation

Sequential numerical solution methods are commonly used for solving the phonon Boltzmann transport equation (BTE) because of simplicity of implementation and low storage requirements. However, they exhibit poor convergence for low Knudsen numbers. This is because sequential solution procedures couple the phonon BTEs in physical space efficiently but the coupling is inefficient in wave vector (K) space. As the Knudsen number decreases, coupling in K space becomes dominant and convergence rates fall. Since materials like silicon have K-resolved Knudsen numbers that span two to five orders of magnitude at room temperature, diffuse-limit solutions are not feasible for all K vectors. Consequently, nongray solutions of the BTE experience extremely slow convergence. In this paper, we develop a coupled-ordinates method for numerically solving the phonon BTE in the relaxation time approximation. Here, interequation coupling is treated implicitly through a point-coupled direct solution of the K-resolved BTEs at each control volume. This implicit solution is used as a relaxation sweep in a geometric multigrid method which promotes coupling in physical space. The solution procedure is benchmarked against a traditional sequential solution procedure for thermal transport in silicon. Significant acceleration in computational time, between 10 and 300 times, over the sequential procedure is found for heat conduction problems.

Nonsymmetric Preconditioning for Conjugate Gradient and Steepest Descent Methods 1

We numerically analyze the possibility of turning off post-smoothing (relaxation) in geometric multigrid when used as a preconditioner in conjugate gradient linear and eigenvalue solvers for the 3D Laplacian. The geometric Semicoarsening Multigrid (SMG) method is provided by the hypre parallel software package. We solve linear systems using two variants (standard and flexible) of the preconditioned conjugate gradient (PCG) and preconditioned steepest descent (PSD) methods. The eigenvalue problems are solved using the locally optimal block preconditioned conjugate gradient (LOBPCG) method available in hypre through BLOPEX software. We observe that turning off the post-smoothing in SMG dramatically slows down the standard PCG-SMG. For flexible PCG and LOBPCG, our numerical results show that post-smoothing can be avoided, resulting in overall acceleration, due to the high costs of smoothing and relatively insignificant decrease in convergence speed. We numerically demonstrate for linear systems that PSD-SMG and flexible PCG-SMG converge similarly if SMG post-smoothing is off. We experimentally show that the effect of acceleration is independent of memory interconnection. A theoretical justification is provided.

The Closest Point Method and Multigrid Solvers for Elliptic Equations on Surfaces

Elliptic partial differential equations are important from both application and analysis points of view. In this paper we apply the closest point method to solve elliptic equations on general curved surfaces. Based on the closest point representation of the underlying surface, we formulate an embedding equation for the surface elliptic problem, then discretize it using standard finite differences and interpolation schemes on banded but uniform Cartesian grids. We prove the convergence of the difference scheme for the Poisson's equation on a smooth closed curve. In order to solve the resulting large sparse linear systems, we propose a specific geometric multigrid method in the setting of the closest point method. Convergence studies in both the accuracy of the difference scheme and the speed of the multigrid algorithm show that our approaches are effective.

Performance of geometric multigrid method for coupled two-dimensional systems in CFD

The performance of a geometric multigrid method is analyzed for two-dimensional Laplace, Navier, Burgers and two formulations of Navier–Stokes (streamfunction–vorticity and streamfunction–velocity) equations. These equations are discretized with the Finite Difference Method on uniform grids with numerical approximations of first- and second-orders of accuracy. The systems of equations are solved with a Modified Strongly Implicit (MSI) and a Successive Over Relaxation (SOR) solver associated with the multigrid method with a V-cycle and a Correction Scheme (CS) and a Full Approximation Scheme (FAS). The effect of the number of inner iterations of the solver, the number of grid levels problems with grid sizes of 1025×1025 points, the influence of differential equations numbers and Reynolds number up to 1000 on Central Processing Unit (CPU) time are investigated. The results show that (1) a solution of two coupled equations (Navier or Burgers) is obtained with the same efficiency multigrid textbook that occurs in the solution of only one equation (Laplace), (2) the efficiency of the multigrid method in the solution of two coupled equations (Navier–Stokes streamfunction–vorticity formulation) or only one equation (the Navier–Stokes streamfunction–velocity formulation) decreases with increasing Reynolds numbers, and (3) the poor performance of the multigrid method for solving the Navier–Stokes seems to be related to the physics of the problem and not to the type of formulation or coupling between the equations.

Geometric multigrid for an implicit-time immersed boundary method

The immersed boundary (IB) method is an approach to fluid-structure interaction that uses Lagrangian variables to describe the structure and Eulerian variables to describe the fluid. Explicit time stepping schemes for the IB method require solvers only for Eulerian equations, for which fast Cartesian grid solution methods are available. Such methods are relatively straightforward to develop and are widely used in practice but often require very small time steps to maintain stability. Implicit-time IB methods permit the stable use of large time steps, but efficient implementations of such methods require significantly more complex solvers that effectively treat both Lagrangian and Eulerian variables simultaneously. Several different approaches to solving the coupled Lagrangian-Eulerian equations have been proposed, but a complete understanding of this problem is still emerging. This paper presents a geometric multigrid method for an implicit-time discretization of the IB equations. This multigrid scheme uses a generalization of box relaxation that is shown to handle problems in which the physical stiffness of the structure is very large. Numerical examples are provided to illustrate the effectiveness and efficiency of the algorithms described herein. These tests show that using multigrid as a preconditioner for a Krylov method yields improvements in both robustness and efficiency as compared to using multigrid as a solver. They also demonstrate that with a time step 100--1000 times larger than that permitted by an explicit IB method, the multigrid-preconditioned implicit IB method is approximately 50--200 times more efficient than the explicit method.

On the multigrid cycle strategy with approximate inverse smoothing

Purpose– The purpose of this paper is to propose multigrid methods in conjunction with explicit approximate inverses with various cycles strategies and comparison with the other smoothers.Design/methodology/approach– The main motive for the derivation of the various multigrid schemes lies in the efficiency of the multigrid methods as well as the explicit approximate inverses. The combination of the various multigrid cycles with the explicit approximate inverses as smoothers in conjunction with the dynamic over/under relaxation (DOUR) algorithm results in efficient schemes for solving large sparse linear systems derived from the discretization of partial differential equations (PDE).Findings– Application of the proposed multigrid methods on two-dimensional boundary value problems is discussed and numerical results are given concerning the convergence behavior and the convergence factors. The results are comparatively better than the V-cycle multigrid schemes presented in a recent report (Filelis-Papadopoulos and Gravvanis).Research limitations/implications– The limitations of the proposed scheme lie in the fact that the explicit finite difference approximate inverse matrix used as smoother in the multigrid method is a preconditioner for specific sparsity pattern. Further research is carried out in order to derive a generic explicit approximate inverse for any type of sparsity pattern.Originality/value– A novel smoother for the geometric multigrid method is proposed, based on optimized banded approximate inverse matrix preconditioner, the Richardson method in conjunction with the DOUR scheme, for solving large sparse linear systems derived from finite difference discretization of PDEs. Moreover, the applicability and convergence behavior of the proposed scheme is examined based on various cycles and comparative results are given against the damped Jacobi smoother.