A Fast Direct Solver for a Class of Elliptic Partial Differential Equations

We describe a fast and robust method for solving the large sparse linear systems that arise upon the discretization of elliptic partial differential equations such as Laplace’s equation and the Helmholtz equation at low frequencies. While most existing fast schemes for this task rely on so called “iterative” solvers, the method described here solves the linear system directly (to within an arbitrary predefined accuracy). The method is described for the particular case of an operator defined on a square uniform grid, but can be generalized other geometries. For a grid containing N points, a single solve requires O(Nlog ² N) arithmetic operations and storage. Storing the information required to perform additional solves rapidly requires O(Nlog N) storage. The scheme is particularly efficient in situations involving domains that are loaded on the boundary only and where the solution is sought only on the boundary. In this environment, subsequent solves (after the first) can be performed in operations. The efficiency of the scheme is illustrated with numerical examples. For instance, a system of size 10⁶×10⁶ is directly solved to seven digits accuracy in four minutes on a 2.8 GHz P4 desktop PC.