Natural domain decomposition algori thms for the solution of time-harmonic elastic waves

Natural domain decomposition algori thms for the solution of time-harmonic elastic waves

R. Brunet Department of Mathematics and Statistics, University of Strathclyde, Glasgow, UK, E-mail:    V. Dolean Department of Mathematics and Statistics, University of Strathclyde, Glasgow, UK, and Laboratoire J.A. Dieudonné, CNRS, University Côte d’Azur, Nice, France. E-mail:    M. J. Gander Université de Genève, 2-4 rue du Lièvre, Genève. E-mail:

We study for the first time Schwarz domain decomposition methods for the solution of the Navier equations modeling the propagation of elastic waves. These equations in the time harmonic regime are difficult to solve by iterative methods, even more so than the Helmholtz equation. We first prove that the classical Schwarz method is not convergent when applied to the Navier equations, and can thus not be used as an iterative solver, only as a preconditioner for a Krylov method. We then introduce more natural transmission conditions between the subdomains, and show that if the overlap is not too small, this new Schwarz method is convergent. We illustrate our results with numerical experiments, both for situations covered by our technical two subdomain analysis, and situations that go far beyond, including many subdomains, cross points, heterogeneous materials in a transmission problem, and Krylov acceleration. Our numerical results show that the Schwarz method with adapted transmission conditions leads systematically to a better solver for the Navier equations than the classical Schwarz method.

Key words. Domain decomposition methods, Schwarz preconditioners, time-harmonic elastic waves, Navier equations.

AMS subject classifications. 65N55, 65N35, 65F10

1 Introduction

Time harmonic problems are difficult to solve by iterative methods in the medium to high frequency regime, see [18] for the case of the Helmholtz equation, which is the prototype of such time harmonic problems with oscillatory solutions. Domain decomposition methods are a natural choice as iterative solvers for such problems, since they are by construction parallel and can still locally use direct solvers without convergence problems. To obtain good domain decomposition convergence for time harmonic problems, adapted transmission conditions are however needed between subdomains. Such transmission conditions were first studied for the Helmholtz equation by Desprès in [10, 11], and later optimized variants were introduced and analyzed by Chevalier in his PhD thesis [7], see also Chevalier and Nataf [8], the work by Collino, Delbue, Joly and Piacentini [9], and Gander et al. [24, 23, 25]. Very similar in nature to the Helmholtz equations, high-frequency time-harmonic Maxwell’s equations are also very difficult to solve iteratively, and the design of efficient domain decomposition methods for the intermediate to high frequency regime is even harder. First optimized transmission conditions both for the first and second order formulations of Maxwell’s equations can already be found in the PhD thesis of Chevalier [7, section 4.7] and Collino et al. [9], but were then more systematically developed by Alonso -Rodriguez and Gerardo-Giorda [1], and especially in Dolean et al. [14, 13, 15], see also Peng, Rawat and Lee [29] , and references therein. The Analytic Incomplete LU factorization (AILU) [19], the sweeping preconditioner [16, 17], the source transfer domain decomposition [5, 6], the method based on single layer potentials [32], and the method of polarized traces [35], are all methods in this same class of domain decomposition methods with more effective transmission conditions, which became known under the name optimized Schwarz methods, see [20, 21] for an introduction, and [26] and references therein for a thorough treatment when applied to time harmonic wave propagation problems.

To the best of our knowledge, the use of Schwarz methods for time-harmonic elastic waves modeled by the Navier equations has not been studied so far, and our goal is to investigate classical Schwarz methods, and also a new variant that uses more natural transmission conditions between the subdomains when applied to the Navier equations. To do so, we study the Schwarz methods at the continuous level, for a simplified decomposition as it has become standard with two subdomains, to gain insight into the effect of transmission conditions on the performance of the method. To test the method, we then discretize the problems and implement the Schwarz methods using Restricted Additive Schwarz (RAS) introduced by Cai and Sarkis in [4], which represents a faithful implementation of the continuous parallel Schwarz method of Lions, see [21]. This is especially important when more natural transmission conditions are used, see [31] for Optimized RAS (ORAS).

Our paper is structured as follows: in Section LABEL:sec:classical, we present and analyze the classical Schwarz algorithm applied to the Navier equations. We prove for a simplified two subdomain setting at the continuous level that the Schwarz algorithm is not a convergent iterative method in this case. We then introduce new transmission conditions in Section LABEL:sec:optimal and show first that there exist transmission conditions which make the Schwarz method converge in a finite number of steps. These transmission conditions involve however non local operators, and we thus introduce a local, low frequency approximation for the Navier equations, for which we prove convergence of the new Schwarz method provided the overlap is not too small. In Section LABEL:sec:numerical we study these new Schwarz methods numerically, first for a two subdomain decomposition covered by our analysis, but then also for the case of many subdomains with cross points and material heterogeneities. Our numerical results show that the new Schwarz method performs much better than the classical one when used as a preconditioner for a Krylov method.

2 Classical Schwarz algorithm for the Navier Equations

We are interested in solving the Navier equations in the frequency domain,


where the operator is defined by . To study the basic (non)-convergence properties of the Schwarz algorithm applied to the Navier equations (LABEL:NavierEq), we consider the domain and decompose it into two unbounded overlapping subdomains and , with overlap parameter . The classical parallel Schwarz algorithm then starts with an initial guess on subdomain , , and solves for iteration index


To study the convergence properties of this algorithm, we use a Fourier transform in the direction. We denote by the Fourier parameter and the Fourier transformed solution,

The convergence factor for each Fourier mode of (LABEL:ClassicalSchwarz) is given in

Lemma 2.1 (Convergence factor of classical Schwarz)

For a given initial guess , , the classical Schwarz algorithm (LABEL:ClassicalSchwarz) with overlap multiplies at each iteration the error in each Fourier mode with the convergence factor


where the eigenvalues of the iteration matrix are


and are given by


Proof. The convergence factor can be obtained by a direct computation working on the error equations, as it is shown in the short publication [3].     

We show in Figure LABEL:Fig3a a plot of the modulus of the convergence factor (LABEL:RhoClassical) as function of the Fourier mode for an example of the parameters in the Navier equations.

Figure 2.1: Modulus of the convergence factor of the classical Schwarz method for , , for different values of the overlap .

We see that the classical Schwarz method converges for high frequencies, , diverges for intermediate frequencies, , and stagnates for low frequencies . We prove in the next theorem that this behavior holds for all choices of parameters in the Navier equations, and thus the classical Schwarz method is not an effective iterative solver for these equations.

Theorem 2.2 ((Non-) Convergence of the overlapping classical Schwarz method)

The convergence factor (LABEL:ClassicalSchwarz) of the overlapping classical Schwarz method (LABEL:ClassicalSchwarz) applied to the Navier equations (LABEL:NavierEq) satisfies


where the last two results are shown to hold for overlap small.

Proof. The proof is quite technical. To simplify the notation, we define for the case when the roots in (LABEL:lambda12) are complex the quantities


We have to treat five cases: three intervals for , and two values separating the intervals: in the first interval , , and the eigenvalues (LABEL:r+r-) become

The square of their modulus is given by


where the complex sign is defined as

and we introduced the quantities , , and ,

The terms and appearing in the square root are real and defined by , which gives after some computations

Then we obtain by a direct computation that


We now show that in (LABEL:rcase1) vanishes identically: we get on the one hand


and on the other hand, we have


and we obtain by adding the three terms from (LABEL:part1b) to each other


This leads, by adding (LABEL:part1a) and (LABEL:part1c) indeed to . We next show that also in (LABEL:rcase1) vanishes identically: we get

and for the term involving and

By analyzing the signs of the different terms, we obtain for the complex sign

and after a lengthy computation we obtain

where is a complicated factor depending on . A direct computation for the second factor of shows that independently of the value of , we get . We can thus conclude from (LABEL:rcase1) that and therefore the algorithm stagnates in the first interval , see the first interval in Figure LABEL:Fig3a.

At the boundary between the first and second interval, where , we have that and , and therefore the eigenvalues in (LABEL:r+r-) become

and being positive we have equivalently

and hence the algorithm stagnates also when the first interval is closed on the right, i.e. for .

In the second interval, , we have that and , and hence the eigenvalues in (LABEL:r+r-) become

We compute the modulus of the eigenvalues and expand them for overlap parameter small to find

We thus obtain that is bigger than one for small and the method diverges, see the middle interval in Figure LABEL:Fig3a111Numerically we observe that also for a large overlap, the algorithm diverges, see Figure LABEL:Fig3a, but this seems to be difficult to prove..

Between the second and third interval, where , we have that and , and hence the eigenvalues in (LABEL:r+r-) become

We thus obtain

and the algorithm stagnates for .

In the last interval, , and by expanding from (LABEL:r+r-) for small, we get

since . We can thus conclude that

see the last interval in Figure LABEL:Fig3a, where we also see that , since all the real exponentials involved in the expressions of are decreasing to as increases.     

We see from Theorem LABEL:ClassicalSchwarzTheorem that the classical Schwarz method with overlap can not be used as an iterative solver to solve the Navier equations, since the method stagnates for low frequencies and even diverges for intermediate frequencies; only high frequencies are converging. A precise estimate for how fast the classical Schwarz method diverges depending on the overlap is given in the short publication [3].

3 New Transmission Conditions for the Schwarz algorithm

A remedy for the divergence problems of the classical Schwarz method is to introduce different transmission conditions, and to consider the new Schwarz method


where the traction operators , , are defined by , and the operators are two by two matrix valued operators one can choose to obtain better convergence. The traction operators play for the Navier equations the role the Neumann condition plays for the Poisson equation. Like we obtained the convergence factor of the classical Schwarz algorithm using a Fourier transform in Lemma LABEL:th:convclas, we can obtain the convergence factor in the case where more general transmission operators with Fourier symbols are used.

Lemma 3.1

For a given initial guess , , the general Schwarz algorithm with overlap (LABEL:GeneralOptimizedSchwarz) has for each Fourier mode the convergence factor






and are given in (LABEL:lambda12).

Proof. This result is obtained by a direct calculation, replacing the solutions in Fourier space into the transmission conditions of the general Schwarz algorithm (LABEL:GeneralOptimizedSchwarz), for details, see the PhD thesis [2, Lemma 2.3].     

3.1 An Optimal Schwarz Method

The new transmission conditions in (LABEL:GeneralOptimizedSchwarz) are a very powerful tool to fix convergence problems of the classical Schwarz method, and are used in many modern domain decomposition methods for time harmonic wave propagation, like the sweeping preconditioner, source transfer and the method of polarized traces, which are all variants of the so called optimized Schwarz methods [20, 21]; for a review, see [26]. To see how powerful this idea is, we start by introducing the best possible choice, namely transparent boundary conditions (TBC) as transmission conditions in (LABEL:GeneralOptimizedSchwarz), which leads to what is called an optimal Schwarz method222Optimal here is not used in the sense of scalability, but really means faster convergence is not possible!:

Theorem 3.2 (Convergence of the optimal Schwarz algorithm.)

If one chooses in the new Schwarz algorithm (LABEL:GeneralOptimizedSchwarz) the operators with the Fourier symbols


where and are given in (LABEL:lambda12), the resulting algorithm converges in two iterations, and this for all values of the overlap , even without overlap, .

Proof. If we replace ) defined in (LABEL:OptimalChoice) into (LABEL:OptiB1), the convergence factor obtained vanishes identically and the algorithm thus converges in two iterations, independently of any initial guess and the overlap .     

To use the optimal choice (LABEL:OptimalChoice) as transmission operators in practice, one needs to back transform the associated TBC into the physical domain, and the corresponding are non local operators, because of the inverse transform with square root terms at the interfaces, like it is the case for many TBCs. It is therefore of interest to design local approximations for the optimal transmissions conditions, like in the development of absorbing boundary conditions (ABCs), which will lead to a new class of practical, so called optimized Schwarz algorithms. We approximate the symbols of the optimal transmission conditions in (LABEL:OptimalChoice) by polynomial symbols in which correspond to derivatives after the Fourier backtransform, and are thus local operators.

3.2 Optimized Schwarz Methods

We have seen in Section LABEL:sec:classical that the classical Schwarz method converges well for high frequency error components, large, but stagnates for low frequency components and even diverges for intermediate range frequencies, see Figure LABEL:Fig3a. It is therefore natural to approximate the operators in the transmission conditions using a low frequency expansion in the Fourier variable of the optimal choice given in Theorem LABEL:OptimalSchwarzThm. This leads to the so called Taylor transmission conditions (TTC), which have the symbols