Jacobi Smoother Research Articles

A scalable physically consistent particle method for high-viscous incompressible flows

AbstractA scalable matrix solver was developed for the moving particle hydrodynamics for incompressible flows (MPH-I) method. Since the MPH-I method can calculate both incompressible and highly viscous flows while ensuring stability through physical consistency, a wide range of industrial applications is expected. However, in its implicit calculation, both the pressure and velocity must be solved simultaneously via a linear equation with a nondefinite symmetric coefficient matrix. In this study, this nondefinite linear system was converted into a symmetric positive definite (SPD) system where only the velocity is unknown. This conversion enabled us to solve the system with well-known solvers such as the conjugated gradient (CG) and conjugated residual (CR) methods. For scalability, bucket-based multigrid preconditioned CG and CR solvers were developed for the SPD system. To handle multidimensionality during preconditioning, an extended Jacobi smoother that is even applicable in a nondiagonally dominant matrix system was proposed. The numerical efficiency was confirmed via a simple high-viscosity incompressible dam break calculation, and the scalability within the presented case was confirmed. In addition, the performance under shared memory parallel computations was studied.

Computing tensor operator exponentials within low‐rank tensor formats with application to the parameter‐dependent multigrid method

AbstractUncertainties in physical models can lead to parameter‐dependent linear systems. The representation and solution of these systems are an important task in numerical mathematics. We summarize our previous results on how to represent these systems using low‐rank tensor methods and how to solve these systems using the parameter‐dependent multigrid method. We propose a new approach to compute the tensor operator exponential, by which we mean the matrix exponential applied to a tensor operator, directly within low‐rank tensor formats. This approach is based on classical matrix methods combined with low‐rank arithmetic. The tensor operator exponential within a low‐rank tensor format is used to approximate the inverse diagonal of a low‐rank operator. This approximation is then used as Jacobi smoother for the parameter‐dependent multigrid method. Using this we observe in numerical experiments a grid size independent convergence rate of the multigrid method. Instead of inverting only diagonals of tensor operators, our approach also allows for the inversion of all symmetric positive definite tensor operators.

Optimized sparse approximate inverse smoothers for solving Laplacian linear systems

In this paper we propose and analyze new efficient sparse approximate inverse (SAI) smoothers for solving the two-dimensional (2D) and three-dimensional (3D) Laplacian linear system with geometric multigrid methods. Local Fourier analysis shows that our proposed SAI smoother for 2D achieves a much smaller smoothing factor than the state-of-the-art SAI smoother studied in Bolten et al. (2016) [12]. The proposed SAI smoother for 3D cases provides smaller optimal smoothing factor than that of weighted Jacobi smoother. Numerical results validate our theoretical conclusions and illustrate the high-efficiency and high-effectiveness of our proposed SAI smoothers. Such SAI smoothers have the advantage of inherent parallelism.

A parallel space-time multigrid method for the eddy-current equation

We expand the applicabilities and capabilities of an already existing space-time parallel method based on a block Jacobi smoother. First we formulate a more detailed criterion for spatial coarsening, which enables the method to deal with unstructured meshes and varying material parameters. Further we investigate the application to the eddy-current equation, where the non-trivial kernel of the curl operator causes severe problems. This is remedied with a new nodal auxiliary space correction. We proceed to identify convergence rates by local Fourier analysis and numerical experiments. Finally, we present a numerical experiment which demonstrates its excellent scaling properties.

Optimal Complex Relaxation Parameters in Multigrid for Complex-Shifted Linear Systems

We derive optimal complex relaxation parameters minimizing smoothing factors associated with multigrid using red-black successive overrelaxation or damped Jacobi smoothing applied to a class of lin...

Numerical analysis of Krylov multigrid methods for stationary advection-diffusion equation

This paper outlines the problem of applying multigrid methods for the stationary scalar advection-diffusion equation for the given advection filed. We assume that all functions are smooth enough to be represented by discretization methods. Two main discretization are considered: finite difference method for Dirichlet boundary conditions and pseudo-spectral method for the periodic boundary. We test the following smoothers for multigrid methods: Jacobi smoother, Gauss-Seidel smoother, Krylov subspace smoother (GMRES method with and without preconditioners). The analysis is performed in the space formed by the cross product of discretization parameters, diffusion coefficient values, multigrid levels and smoothers.We demonstrate that the most efficient strategy depends on parameter value and given velocity field. Best variants include Gauss-Seidel smoothers which is optimal for advection-dominated problem while multigrid method is used as a preconditioner for a Krylov method. Such methods can be used for spectral or pseudo-spectral methods where explicit dense matrix storage is impossible.

Research on Aggregations for Algebraic Multigrid Preconditioning Methods

Aggregation based algebraic multigrid is one of the most efficient methods to solve sparse linear systems. In this paper, a new method based on cliques and several classical aggregations are implemented based on sparse data structures, and compared in solving sparse linear systems from model partial differential equations. The results show that with Vanek scheme, the number of levels is often less than that with others, and when Jacobi smoothing is used, the number of preconditioned iterations and the time elapsed for iteration are both the least, while the time for setup is too long. From the view point of strong connections, it is beneficial to aggregate the points in a common clique, but the experiments have not shown this privilege. Instead, the two-point aggregation based on strong connection is effective in most cases, and the four-point and the eight-point aggregations which are based on loop of the two-point algorithm with two or three times respectively are effective too. In respect to total computation time, the best scheme is always one of these three schemes. In addition, the larger the average degree the related adjacent graph has, the less time required for the Vanek scheme to setup is and the more efficient it is, so is the aggregation scheme based on cliques.

Hybrid Multi-level solvers for discontinuous Galerkin finite element discrete ordinate diffusion synthetic acceleration of radiation transport algorithms

This paper examines two established preconditioners which were developed to accelerate the solution of discontinuous Galerkin finite element method (DG-FEM) discretisations of the elliptic neutron diffusion equation. They are each presented here as a potential way to accelerate the solution of the Modified Interior Penalty (MIP) form of the discontinuous diffusion equation, for use as a diffusion synthetic acceleration (DSA) of DG-FEM discretisations of the neutron transport equation. The preconditioners are both two-level schemes, differing in the low-level space utilised. Once projected to the low-level space a selection of algebraic multigrid (AMG) preconditioners are utilised to obtain a further correction step, these are therefore “hybrid” preconditioners. The first preconditioning scheme utilises a continuous piece-wise linear finite element method (FEM) space, while the second uses a discontinuous piece-wise constant space. Both projections are used alongside an element-wise block Jacobi smoother in order to create a symmetric preconditioning scheme which may be used alongside a conjugate gradient algorithm. An eigenvalue analysis reveals that both should aid convergence but the piece-wise constant based method struggles with some of the smoother error modes. Both are applied to a range of problems including some which are strongly heterogeneous. In terms of conjugate gradient (CG) iterations needed to reach convergence and computational time required, both methods perform well. However, the piece-wise linear continuous scheme appears to be the more effective of the two. An analysis of computer memory usage found that the discontinuous piece-wise constant method had the lowest memory requirements.

Energy-efficient multigrid smoothers and grid transfer operators on multi-core and GPU clusters

We investigate time and energy to solution for the CPU- and GPU-based execution of the compute intensive smoother and grid transfer operators in a geometric multigrid linear solver. We use a hybrid parallel implementation for both shared and distributed memory multi-core host systems comprising CUDA-capable devices. Our numerical experiments are designed to assess the effect of combining an MPI-parallel multigrid framework with OpenMP host threads or CUDA accelerators instead of MPI-only CPU computations for various parallel setups. We present runtime and energy measurements from a quad-CPU test system equipped with two GPUs. We find that using an accelerated asynchronous smoother can yield substantial savings of time and energy to solution over using a host-only Jacobi smoother in small and medium sized host systems with one or two multi-core CPUs. The acceleration of the grid transfer operators also yields a benefit, yet smaller than the benefit from the smoother. For large host systems a hybrid MPI-OpenMP parallelization turns out to be most beneficial with respect to energy consumption, although it is not the fastest option.

Analysis of a New Space-Time Parallel Multigrid Algorithm for Parabolic Problems

We present and analyze a new space-time parallel multigrid method for parabolic equations. The method is based on arbitrarily high order discontinuous Galerkin discretizations in time and a finite element discretization in space. The key ingredient of the new algorithm is a block Jacobi smoother. We present a detailed convergence analysis when the algorithm is applied to the heat equation and determine asymptotically optimal smoothing parameters, a precise criterion for semi-coarsening in time or full coarsening, and give an asymptotic two grid contraction factor estimate. We then explain how to implement the new multigrid algorithm in parallel and show with numerical experiments its excellent strong and weak scalability properties.

Comparison of multigrid algorithms for high‐order continuous finite element discretizations

SummaryWe present a comparison of different multigrid approaches for the solution of systems arising from high‐order continuous finite element discretizations of elliptic partial differential equations on complex geometries. We consider the pointwise Jacobi, the Chebyshev‐accelerated Jacobi, and the symmetric successive over‐relaxation smoothers, as well as elementwise block Jacobi smoothing. Three approaches for the multigrid hierarchy are compared: (1) high‐order h‐multigrid, which uses high‐order interpolation and restriction between geometrically coarsened meshes; (2) p‐multigrid, in which the polynomial order is reduced while the mesh remains unchanged, and the interpolation and restriction incorporate the different‐order basis functions; and (3) a first‐order approximation multigrid preconditioner constructed using the nodes of the high‐order discretization. This latter approach is often combined with algebraic multigrid for the low‐order operator and is attractive for high‐order discretizations on unstructured meshes, where geometric coarsening is difficult. Based on a simple performance model, we compare the computational cost of the different approaches. Using scalar test problems in two and three dimensions with constant and varying coefficients, we compare the performance of the different multigrid approaches for polynomial orders up to 16. Overall, both h‐multigrid and p‐multigrid work well; the first‐order approximation is less efficient. For constant coefficients, all smoothers work well. For variable coefficients, Chebyshev and symmetric successive over‐relaxation smoothing outperform Jacobi smoothing. While all of the tested methods converge in a mesh‐independent number of iterations, none of them behaves completely independent of the polynomial order. When multigrid is used as a preconditioner in a Krylov method, the iteration number decreases significantly compared with using multigrid as a solver. Copyright © 2015 John Wiley & Sons, Ltd.

Domain-Specific Optimization of Two Jacobi Smoother Kernels and Their Evaluation in the ECM Performance Model

Our aim is to apply program transformations to stencil codes in order to yield the highest possible performance. We recognize memory bandwidth as a major limitation in stencil code performance. We conducted a study in which we applied optimizing transformations to two Jacobi smoother kernels: one 3D 1st-order 7-point stencil and one 3D 3rd-order 19-point stencil. To obtain high performance, the optimizations have to be customized for the execution platform at hand. We illustrate this by experiments on two consumer and two server architectures. We also verified the need for complex optimizations with the help of the Execution-Cache-Memory performance model. A code generator with knowledge about stencil codes and execution platforms should be able to apply our transformations automatically. We are working towards such a generator in project ExaStencils.

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.

Towards a complete FEM-based simulation toolkit on GPUs: Unstructured grid finite element geometric multigrid solvers with strong smoothers based on sparse approximate inverses

We describe our FE-gMG solver, a finite element geometric multigrid approach for problems relying on unstructured grids. We augment our GPU- and multicore-oriented implementation technique based on cascades of sparse matrix–vector multiplication by applying strong smoothers. In particular, we employ Sparse Approximate Inverse (SPAI) and Stabilised Approximate Inverse (SAINV) techniques. We focus on presenting the numerical efficiency of our smoothers in combination with low- and high-order finite element spaces as well as the hardware efficiency of the FE-gMG. For a representative problem and computational grids in 2D and 3D, we achieve a speedup of an average of 5 on a single GPU over a multithreaded CPU code in our benchmarks. In addition, our strong smoothers can deliver a speedup of 3.5 depending on the element space, compared to simple Jacobi smoothing. This can even be enhanced to a factor of 7 when combining the usage of approximate inverse-based smoothers with clever sorting of the degrees of freedom. In total the FE-gMG solver can outperform a simple (multicore-) CPU-based multigrid by a total factor of over 40.

Multigrid Method With Adaptive IDR-Based Jacobi Smoother

Multigrid (MG) methods are known to be fast linear solvers for large-scale finite-element analyses. The Gauss-Seidel method is usually adopted as the smoother for MG methods. However, recently, considerable attention has focused on induced dimension reduction (IDR)-based solvers because they are faster. In this paper, we investigate the convergence of IDR-based solvers and evaluate the performance of an MG method with an adaptive IDR-based Jacobi smoother. Numerical results show that this method has good convergence and good efficiency in parallel computations for finite-element analysis of electromagnetic fields.

A multigrid method for weakly nonlinear elliptic equations of the second order

We consider the Dirichlet problem for a weakly nonlinear elliptic equation of the second order in the divergence form. We prove the convergence of a multigrid method for solving this problem. The method under consideration is based on the application of conformal finite elements and the Jacobi smoothing procedure.

Peano—A Traversal and Storage Scheme for Octree-Like Adaptive Cartesian Multiscale Grids

Almost all approaches to solving partial differential equations (PDEs) are based upon a spatial discretization of the computational domain—a grid. This paper presents an algorithm to generate, store, and traverse a hierarchy of d-dimensional Cartesian grids represented by a $(k=3)$-spacetree, a generalization of the well-known octree concept, and it also shows the correctness of the approach. These grids may change their adaptive structure throughout the traversal. The algorithm uses $2d+4$ stacks as data structures for both cells and vertices, and the storage requirements for the pure grid reduce to one bit per vertex for both the complete grid connectivity structure and the multilevel grid relations. Since the traversal algorithm uses only stacks, the algorithm's cache hit rate is continually higher than 99.9 percent, and the runtime per vertex remains almost constant; i.e., it does not depend on the overall number of vertices or the adaptivity pattern. We use the algorithmic approach as the fundamental concept for a mesh management for d-dimensional PDEs and for a matrix-free PDE solver represented by a compact discrete $3^d$-point operator. In the latter case, one can implement a Jacobi smoother, a Krylov solver, or a geometric multigrid scheme within the presented traversal scheme which inherits the low memory requirements and the good memory access characteristics directly.

A multigrid‐based shifted Laplacian preconditioner for a fourth‐order Helmholtz discretization

AbstractIn this paper, an iterative solution method for a fourth‐order accurate discretization of the Helmholtz equation is presented. The method is a generalization of that presented in (SIAM J. Sci. Comput. 2006; 27:1471–1492), where multigrid was employed as a preconditioner for a Krylov subspace iterative method. The multigrid preconditioner is based on the solution of a second Helmholtz operator with a complex‐valued shift. In particular, we compare preconditioners based on a point‐wise Jacobi smoother with those using an ILU(0) smoother, we compare using the prolongation operator developed by de Zeeuw in (J. Comput. Appl. Math. 1990; 33:1–27) with interpolation operators based on algebraic multigrid principles, and we compare the performance of the Krylov subspace method Bi‐conjugate gradient stabilized with the recently introduced induced dimension reduction method, IDR(s). These three improvements are combined to yield an efficient solver for heterogeneous problems. Copyright © 2009 John Wiley & Sons, Ltd.

Integration of hp-adaptivity and a two grid solver for electromagnetic problems

We present implementation details and analyze convergence of a two grid solver forming the core of a fully automatic hp-adaptive strategy for electromagnetic problems. The solver delivers a solution for a fine grid obtained from an arbitrary coarse hp grid by a global hp-refinement. The classical V-cycle algorithm combines an overlapping block Jacobi smoother with optimal relaxation, and a direct solve on the coarse grid. A theoretical analysis of the two grid solver is illustrated with numerical experiments. Several electromagnetic applications show the efficiency of combining the fully automatic hp-adaptive strategy with the two grid solver.

Integration of hp-adaptivity and a two-grid solver for elliptic problems

We present implementation details and analyze convergence of a two-grid solver forming the core of the fully automatic hp-adaptive strategy for elliptic problems. The solver delivers a solution for a fine grid obtained from an arbitrary coarse hp-grid by a global hp-refinement. The classical V-cycle algorithm combines an overlapping block Jacobi smoother with optimal relaxation, and a direct solve on the coarse grid. A simple theoretical analysis is illustrated with extensive numerical experimentation.