Title:

Robust solvers for large indefinite systems in seismic inversion

In this thesis we study and develop iterative methods for solving the linear systems arising from discretisations of the Helmholtz equation. The Helmholtz equation is the simplest model used to describe high frequency wave scattering, and hence arises in many applications. Therefore it is of practical importance for there to be robust numerical methods for the solution of the Helmholtz equation, which is the focus of this thesis. Throughout this thesis we consider as our model problem the Helmholtz equation in a bounded domain subject to an impedance boundary condition which is an approximation of the so called Sommerfeld radiation condition. Once the PDE problem is then discretised with, say, low order finite elements the resulting system matrix is complex, symmetric and nonHermitian. However for large values of the wavenumber k the solution of the Helmholtz equation is highly oscillatory and therefore a number of grid points growing at least as fast as O(k) must be chosen to ensure an accurate solution. The consequence of this is that for large k, or equivalently large domains, the system matrix will be very large. Therefore direct solvers are no longer a viable choice, and hence we turn instead to \textit{cheaper} iterative solvers. In this thesis we consider the Schwarz domain decomposition as our choice of iterative method and preconditioner. Throughout we shall prove theoretical results which show that this algorithm converges with a number of iterations which grows like k with a fractional exponent. This exponent differs depending on the interface condition and overlap. A new theory is developed for the optimised nonoverlapping Schwarz method for the Helmholtz equation with an absorbing term. The optimised method works by first replacing the multiplier in the standard impedance condition with a general complex number determined by minimising the convergence rate of the Schwarz method. The interface condition is known as an optimised interface condition. The consequence of this is that we can improve upon the previous methods with Dirichlet and standard impedance conditions. This is verified theoretically and then numerically by testing the Dirichlet, impedance and optimised interface conditions using the Schwarz method on some model problems. There are however situations where we cannot solve the optimisation problem resulting from the Schwarz method exactly. This situation arises when we consider higher order approximations for the interface condition, and when we use overlap in the optimised method. Instead of solving the problem theoretically we solve it numerically to obtain the multiplier required for the optimised interface condition. What we observe in the numerical experiments is that we can improve the convergence of the Schwarz algorithm further. Finally in this thesis we develop a new preconditioner which combines the ideas of the optimised Schwarz domain decomposition method with the sweeping precondtioner of Engquist and Ying to form a hybrid preconditioner. This hybrid preconditioner is used within GMRES, and tested in numerical computations on large scale geological problems, including some 3D examples. What we find is that this new preconditioner can converge very quickly. However, we also highlight areas which could be improved using other techniques which we outline as further work to be considered.
