GPU implementation of a Helmholtz Krylov solver preconditioned by a shifted Laplace multigrid method

A Helmholtz equation in two dimensions discretized by a second order finite difference scheme is considered. Krylov methods such as Bi-CGSTAB and IDR(s) have been chosen as solvers. Since the convergence of the Krylov solvers deteriorates with increasing wave number, a shifted Laplace multigrid preconditioner is used to improve the convergence. The implementation of the preconditioned solver on CPU (Central Processing Unit) is compared to an implementation on GPU (Graphics Processing Units or graphics card) using CUDA (Compute Unified Device Architecture). The results show that preconditioned Bi-CGSTAB on GPU as well as preconditioned IDR(s) on GPU is about 30 times faster than on CPU for the same stopping criterion.     

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

In 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.     

Vekua theory for the Helmholtz operator

Vekua operators map harmonic functions defined on domain in \mathbb R²{\mathbb R^{2}} to solutions of elliptic partial differential equations on the same domain and vice versa. In this paper, following the original work of I. Vekua (Ilja Vekua (1907–1977), Soviet-Georgian mathematician), we define Vekua operators in the case of the Helmholtz equation in a completely explicit fashion, in any space dimension N ≥ 2. We prove (i) that they actually transform harmonic functions and Helmholtz solutions into each other; (ii) that they are inverse to each other; and (iii) that they are continuous in any Sobolev norm in star-shaped Lipschitz domains. Finally, we define and compute the generalized harmonic polynomials as the Vekua transforms of harmonic polynomials. These results are instrumental in proving approximation estimates for solutions of the Helmholtz equation in spaces of circular, spherical, and plane waves.     

Vekua theory for the Helmholtz operator

Vekua operators map harmonic functions defined on domain in \({\mathbb R^{2}}\) to solutions of elliptic partial differential equations on the same domain and vice versa. In this paper, following the original work of I. Vekua (Ilja Vekua (1907–1977), Soviet-Georgian mathematician), we define Vekua operators in the case of the Helmholtz equation in a completely explicit fashion, in any space dimension N ≥ 2. We prove (i) that they actually transform harmonic functions and Helmholtz solutions into each other; (ii) that they are inverse to each other; and (iii) that they are continuous in any Sobolev norm in star-shaped Lipschitz domains. Finally, we define and compute the generalized harmonic polynomials as the Vekua transforms of harmonic polynomials. These results are instrumental in proving approximation estimates for solutions of the Helmholtz equation in spaces of circular, spherical, and plane waves.     

A multigrid-based preconditioned Krylov subspace method for the Helmholtz equation with PML

In this paper, we generalize the complex shifted Laplacian preconditioner to the complex shifted Laplacian-PML preconditioner for the Helmholtz equation with perfectly matched layer (Helmholtz-PML equation). The Helmholtz-PML equation is discretized by an optimal 9-point difference scheme, and the preconditioned linear system is solved by the Krylov subspace method, especially by the biconjugate gradient stabilized method (Bi-CGSTAB). The spectral analysis of the linear system is given, and a new matrix-based interpolation operator is proposed for the multigrid method, which is used to approximately invert the preconditioner. The numerical experiments are presented to illustrate the efficiency of the preconditioned Bi-CGSTAB method with the multigrid based on the new interpolation operator, also, numerical results are given for comparing the performance of the new interpolation operator with that of classic bilinear interpolation operator and the one suggested in Erlangga et al. (2006) [10].     

Harmonic analysis associated with a discrete Laplacian

It is well known that the fundamental solution of

$${u_t}\left( {n,t} \right) = u\left( {n + 1,t} \right) - 2u\left( {n,t} \right) + u\left( {n - 1,t} \right),n \in \mathbb{Z},$$

with u(n, 0) = δ _nm for every fixed m ∈ Z is given by u(n, t) = e ^?2t I _n?m (2t), where I _k (t) is the Bessel function of imaginary argument. In other words, the heat semigroup of the discrete Laplacian is described by the formal series W _t f(n) = Σ _m∈Z e ^?2t I _n?m (2t)f(m). This formula allows us to analyze some operators associated with the discrete Laplacian using semigroup theory. In particular, we obtain the maximum principle for the discrete fractional Laplacian, weighted ? ^p (Z)-boundedness of conjugate harmonic functions, Riesz transforms and square functions of Littlewood-Paley. We also show that the Riesz transforms essentially coincide with the so-called discrete Hilbert transform defined by D. Hilbert at the beginning of the twentieth century. We also see that these Riesz transforms are limits of the conjugate harmonic functions. The results rely on a careful use of several properties of Bessel functions.

     

Liouville-type theorem for the drifting Laplacian operator

Given manifolds with a smooth measure (M, g, e ^?f dV), we consider gradient estimates for positive harmonic functions of the drifting Laplacian. If the ∞-Bakry-Emery Ricci tensor is bounded from below and \({|\nabla f|}\) is bounded, we obtain a Liouville-type theorem. This extends a classical result of Cheng and Yau.     

An integral operator connected with the Helmholtz equation

Let Lu be the integral operator defined by $\text{(L}{}_{\mathrm{k?}}\text{)(x, y) = ∝ s ∝ ?(x′, y′)(}\text{e}{}^{\mathrm{ik?}}\text{?}\text{) dx′ dy′, (x, y) ? S}$ where S is the interior of a smooth, closed Jordan curve in the plane, k is a complex number with Re k ? 0, Im k ? 0, and ?² = (x ?x′)² + (y ? y′)². We define $\text{q(x, y) = [}\text{dist}\text{((x, y), ?S)]}\text{1}\text{2}\text{, (x, y) ? S; L}{}_2\text{(q, S) = {? : ∝ s ∝ ¦ ?(x, y)¦}{}^2\text{q(x, y) dx dy < ∞}; W}{}_2{}^1\text{(q, S) = {? : ? ? L}{}_2\text{(q, S),}\text{??}\text{?x}\text{,}\text{?f}\text{?y}\text{? L}{}_2\text{(q, S)}}$, where in the definition of W₂¹(q, S) the derivatives are taken in the sense of distributions. We prove that L_k is a continuous 1-l mapping of L₂(q, S) onto W₂¹(q, S).     

Estimates for eigenvalues of Laplacian operator with any order

Let D be a bounded domain in an n-dimensional Euclidean space Rn. Assume that 0 ＜ λ1 ≤λ2 ≤ … ≤ λκ ≤ … are the eigenvalues of the Dirichlet Laplacian operator with any order l{(-△)lu=λu, in D u=(δ)u/(δ)(→n)=…(δ)l-1u/(δ)(→n)l-1=0,on (δ)D.Then we obtain an upper bound of the (k 1)-th eigenvalue λκ 1 in terms of the first k eigenvalues.k∑i=1(λκ 1-λi) ≤ 1/n[4l(n 2l-2)]1/2{k∑i=1(λκ 1-λi)1/2λil-1/l k∑i=1(λκ 1-λi)1/2λ1/li}1/2.This ineguality is independent of the domain D. Furthermore, for any l ≥ 3 the above inequality is better than all the known results. Our rusults are the natural generalization of inequalities corresponding to the case l = 2 considered by Qing-Ming Cheng and Hong-Cang Yang. When l = 1, our inequalities imply a weaker form of Yang inequalities. We aslo reprove an implication claimed by Cheng and Yang.     

The Goursat problem for a generalized Helmholtz operator in the plane

We consider the Goursat problem in the plane for partial differential operators whose principal part is the pth power of the standard Laplace operator. The data is posed on a union of 2p distinct lines through the origin. We show that the solvability of this Goursat problem depends on Diophantine properties of the geometry of lines on which the data is posed. The first author is supported in part by DMS-0401215. The second author is supported in part by Grant MTM2006-13000-C03-03 of the D.G.I. of Spain.     

The form boundedness criterion for the Laplacian operator

In this paper, we characterize the class of measurable functions (or, more generally, real- or complex-valued distributions) f such that the commutator operator

Cf=[f,Δ]     

On shooting methods for the discrete Helmholtz equation with constant coefficients

Summary The paper is concerned with shooting solvers for the Helmholtz equation with constant coefficients in two dimensions using finite differences for the discretization. Dirichlet boundary conditions are treated though other conditions are possible. Beginning with a single shooting method some recursive multiple shooting methods are developed. It will be shown that the performance of the algorithms may be improved considerably by a redundance-free recursion. The number of operations required for one solution will be computed, but without preparing some matrices which do not depend on the boundary conditions and the inhomogenity. For a square withn×n points the number is of the orderO(n ²⁺⁽ⁿ⁾) with . The method will be compared with a multi-grid program and finally — as an example—a Stokes-solver and some numerical results with the shooting method are given.     

Use of Shifted Laplacian Operators for Solving Indefinite Helmholtz Equations

A shifted Laplacian operator is obtained from the Helmholtz operator by adding a complex damping. It serves as a basic tool in the most successful multigrid approach for solving highly indefinite Helmholtz equations — a Shifted Laplacian preconditioner for Krylov-type methods. Such preconditioning significantly accelerates Krylov iterations, much more so than the multigrid based on original Helmholtz equations. In this paper, we compare approximation and relaxation properties of the Helmholtz operator with and without the complex shift, and, based on our observations, propose a new hybrid approach that combines the two. Our analytical conclusions are supported by two-dimensional numerical results.     

Discontinuous Galerkin approximation with discrete variational principle for the nonlinear Laplacian

A discontinuous Galerkin method is analyzed to approximate the nonlinear Laplacian model problem. The salient feature of the proposed scheme is that it is endowed with a discrete variational principle. The convergence of the discrete approximations to the exact solution is proven. To cite this article: E. Burman, A. Ern, C. R. Acad. Sci. Paris, Ser. I 346 (2008).     

Singular value expansion for the Green function of Helmholtz operator

We consider the integral operator defined on a circular disk, and with kernel the Green function of the Helmholtz operator. We present an analytic framework for the explicit computation of the singular system of this kernel. In particular, the main formulas of this framework are given by a characteristic equation for the singular values and explicit expressions for the corresponding singular functions. We provide also a property of the singular values, that gives an important information for the numerical evaluation of the singular system. Finally, we present a simple numerical experiment, where the singular system computed by a simple implementation of these analytic formulas is compared with the singular system obtained by a discretization of the Green function of the Helmholtz operator.     

An addition theorem for the manifolds with the Laplacian having discrete spectrum

The question of the preservation of discreteness of the spectrum of the Laplacian acting in a space of differential forms under the cutting and gluing of manifolds reduces to the same problem for compact solvability of the operator of exterior derivation. Along these lines, we give some conditions on a cut Y dividing a Riemannian manifold X into two parts X ₊ and X _? under which the spectrum of the Laplacian on X is discrete if and only if so are the spectra of the Laplacians on X ₊ and X _?.