Andrew Hassell

Alex H. Barnett, Andrew Hassell

We present a new algorithm for numerical computation of large eigenvalues and associated eigenfunctions of the Dirichlet Laplacian in a smooth, star-shaped domain in $\mathbb{R}^d$, $d\ge 2$. Conventional boundary-based methods require a root-search in eigenfrequency $k$, hence take $O(N^3)$ effort per eigenpair found, using dense linear algebra, where $N=O(k^{d-1})$ is the number of unknowns required to discretize the boundary. Our method is O(N) faster, achieved by linearizing with respect to $k$ the spectrum of a weighted interior Neumann-to-Dirichlet (NtD) operator for the Helmholtz equation. Approximations $\hat{k}_j$ to the square-roots $k_j$ of all O(N) eigenvalues lying in $[k - ε, k]$, where $ε=O(1)$, are found with $O(N^3)$ effort. We prove an error estimate $$ |\hat k_j - k_j| \leq C \Big(\frac{ε^2}{k} + ε^3 \Big), $$ with $C$ independent of $k$. We present a higher-order variant with eigenvalue error scaling empirically as $O(ε^5)$ and eigenfunction error as $O(ε^3)$, the former improving upon the 'scaling method' of Vergini--Saraceno. For planar domains ($d=2$), with an assumption of absence of spectral concentration, we also prove rigorous error bounds that are close to those numerically observed. For $d=2$ we compute robustly the spectrum of the NtD operator via potential theory, Nyström discretization, and the Cayley transform. At high frequencies (400 wavelengths across), with eigenfrequency relative error $10^{-10}$, we show that the method is $10^3$ times faster than standard ones based upon a root-search.

A. H. Barnett, Andrew Hassell

We study eigenfunctions and eigenvalues of the Dirichlet Laplacian on a bounded domain $Ω\subset\RR^n$ with piecewise smooth boundary. We bound the distance between an arbitrary parameter $E > 0$ and the spectrum $\{E_j \}$ in terms of the boundary $L^2$-norm of a normalized trial solution $u$ of the Helmholtz equation $(Δ+ E)u = 0$. We also bound the $L^2$-norm of the error of this trial solution from an eigenfunction. Both of these results are sharp up to constants, hold for all $E$ greater than a small constant, and improve upon the best-known bounds of Moler--Payne by a factor of the wavenumber $\sqrt{E}$. One application is to the solution of eigenvalue problems at high frequency, via, for example, the method of particular solutions. In the case of planar, strictly star-shaped domains we give an inclusion bound where the constant is also sharp. We give explicit constants in the theorems, and show a numerical example where an eigenvalue around the 2500th is computed to 14 digits of relative accuracy. The proof makes use of a new quasi-orthogonality property of the boundary normal derivatives of the eigenmodes, of interest in its own right.

Alex Barnett, Andrew Hassell, Melissa Tacy

For smooth bounded domains in $\mathbb{R}$, we prove upper and lower $L^2$ bounds on the boundary data of Neumann eigenfunctions, and prove quasi-orthogonality of this boundary data in a spectral window. The bounds are tight in the sense that both are independent of eigenvalue; this is achieved by working with an appropriate norm for boundary functions, which includes a `spectral weight', that is, a function of the boundary Laplacian. This spectral weight is chosen to cancel concentration at the boundary that can happen for `whispering gallery' type eigenfunctions. These bounds are closely related to wave equation estimates due to Tataru. Using this, we bound the distance from an arbitrary Helmholtz parameter $E>0$ to the nearest Neumann eigenvalue, in terms of boundary normal-derivative data of a trial function $u$ solving the Helmholtz equation $(Δ-E)u=0$. This `inclusion bound' improves over previously known bounds by a factor of $E^{5/6}$. It is analogous to a recently improved inclusion bound in the Dirichlet case, due to the first two authors. Finally, we apply our theory to present an improved numerical implementation of the method of particular solutions for computation of Neumann eigenpairs on smooth planar domains. We show that the new inclusion bound improves the relative accuracy in a computed Neumann eigenvalue (around the $42000$th) from 9 digits to 14 digits, with little extra effort.

A. H. Barnett, Andrew Hassell

We consider the "Method of particular solutions" for numerically computing eigenvalues and eigenfunctions of the Laplacian $Δ$ on a smooth, bounded domain Omega in RR^n with either Dirichlet or Neumann boundary conditions. This method constructs approximate eigenvalues E, and approximate eigenfunctions u that satisfy $Δu=Eu$ in Omega, but not the exact boundary condition. An inclusion bound is then an estimate on the distance of E from the actual spectrum of the Laplacian, in terms of (boundary data of) u. We prove operator norm estimates on certain operators on $L^2(\partial Ω)$ constructed from the boundary values of the true eigenfunctions, and show that these estimates lead to sharp inclusion bounds in the sense that their scaling with $E$ is optimal. This is advantageous for the accurate computation of large eigenvalues. The Dirichlet case can be treated using elementary arguments and has appeared in SIAM J. Num. Anal. 49 (2011), 1046-1063, while the Neumann case seems to require much more sophisticated technology. We include preliminary numerical examples for the Neumann case.