A coercive heterogeneous media Helmholtz model: formulation, wavenumber-explicit analysis, and preconditioned high-order FEM

Abstract

We consider a frequency-domain heterogeneous wave propagation model governed by the Helmholtz partial differential equation (PDE) and an impedance boundary condition. The celebrated standard (H¹) variational formulation of the model is non-coercive. It is an open problem to establish a coercive variational formulation of the heterogeneous model. The main focus of this article is on solving this continuous model formulation and analysis problem, and hence establishing an efficient preconditioned numerical algorithm for simulating our novel coercive variational formulation. We develop the variational formulation for the heterogeneous model (in a Hilbert space V equipped with a stronger norm than the H¹-norm) and prove that the associated sesquilinear form is coercive, with a wavenumber-independent coercivity constant. We use this result to derive a wavenumber-independent bound for solutions of the heterogeneous media wave propagation model in the V -norm. Additionally, we prove continuity of the sesquilinear form, with a wavenumber-explicit continuity constant. Using our analysis-supported coercive formulation, we develop a high-order frequency robust-preconditioned finite element method (FEM)-based heterogeneous media discrete wave model. For demonstrating efficiency and convergence of the coercive high-order FEM model, we use non-convex media comprising curved and non-smooth boundaries and low- to high-frequency input incident waves. For the heterogeneous media, with size varying from tens to hundreds of wavelengths, we demonstrate that our new preconditioned-FEM model requires a very low number of GMRES iterations, and the number of iterations is independent of the wavenumber of the model. We also use a class of additive Schwarz domain decomposition (DD) algorithms to implement the preconditioned-FEM model. The DD-based high-order preconditioned-FEM results and comparisons further demonstrate efficiency of the coercive formulation to simulate wave propagation in heterogenous media.

Appendix: High-order WEB-splines FEM

We omit the details of the WEB-method implemented in this article and refer to [22, 29].

1.1 A.1 Basis implementation example

For the star-shaped geometry used as an example in this paper, we illustrate the construction of the WEB-spline basis using h = (1/2)². We begin by placing a grid with uniform width h = (1/2)² over the domain, and we make sure that the grid is large enough to contain all of the supports of the relevant splines. Then, the grid cells are classified, identifying each grid cell by the lower left vertex. The grid choice and classification are shown in Fig. 5a. We take advantage of the rectangular portions of the geometry, and choose a grid for the FEM basis which fits this portion well. Then, extension of the basis splines is only needed near the curved portion of the geometry. While the initial grid has a uniform width h, the resulting WEB-splines have non-uniform support sizes.

Fig. 5

Basis spline and grid cell classification for varying spline degree p when h = (1/2)²

The grid-cell classification does not change with choice of spline degree p, but the classification of basis splines and the extension process do because the size of the supports of the splines increases as p increases. We implement splines of degree p = 2, 3, and 4 which are $\mathcal {C}^{p-1}({\Omega })$ continuous, and similar to cells we identify splines by the cell and vertex in the lower left corner of their supports. The classification and extension for p = 2, 3, and 4 are shown in Fig. 5b, c, and d respectively.

1.2 A.2 Remark on Degrees of Freedom

We choose to implement the new formulations with WEB-splines for many reasons. One advantage of this method is that for a fixed grid width h, we can increase the degree p of the splines used without a large increase in the degrees of freedom (DoF) for the problem.

Traditional triangular elements with polynomial basis splines of degree p require DoF for each element equal to (1/2)(p + 1)(p + 2), and hence, for a fixed h, the DoF grow rapidly as p is increased. This is not the case with the WEB-spines basis construction.

We demonstrate the DoF for various choices of h when increasing p in Table 28. Thus, we can achieve high-order convergence without adding much computational cost. There is a small increase in the density of the matrices as p is increased which can be seen in Table 28. Additionally, the change in density for different p for a coarse grid with h = (1/2)² and a fine grid with h = (1/2)⁵ can be seen in the sparsity pattern plots given in Fig. 6.

Table 28 Degrees of freedom (DoF) and density given in the ratio of nonzero to total entries in the system for the WEB-spline basis with various grid width h and spline degree p for the example geometry shown in Fig. 4

Fig. 6

Sparsity patterns for the Galerkin FEM system with various p and h

The number of nonzero entries (nnz) in a row of the FEM system corresponding to an un-modified interior spline is less than or equal to (p + 1)^d. Additionally, from [29], the ratio of outer splines to inner splines decreases like O(h). Thus, most rows in the system have nnz ≤ (p + 1)^d, and the number of entries in the system is less than or equal to C_DoFh^−d for a constant C_DoF > 0. For the example problem here, we expect to see order, two growth in the number of nonzero entries for the system as the mesh size is decreased independent of the choice of p. This is verified in Table 29 where we present the number of nonzero entries for various p and h and demonstrate that the experimental order of growth (EOG) for the systems used in this paper is approximately two.

Table 29 Number of nonzero entries for the system (nnz) and the experimental order of growth (EOG) for the WEB-spline basis with various grid width h and spline degree p for the example geometry shown in Fig. 4