Computing eigenvalues of ordinary differential equations

Discretisations of differential eigenvalue problems have a sensitivity to perturbations which is asymptotically least as h ? 0 when the differential equation is in first order system form. Both second and fourth order accurate discretisations of the first order system are straightforward to derive and lead to generalised eigenvalue problems of the form ?? where both A and B are narrow-banded, block bidiagonal (hence unsymmetric) matrices, and typically B is singular. Solutions of the differential equation associated with eigenvalues of small magnitude are best determined by the discretisations. Thus Krylov subspace methods (for example) require A to be invertible and seek large solutions of ?? This already requires rational methods in principle. It follows that rapidly convergent methods based on inverse iteration applied to the original formulation as a nonstandard generalised eigenvalue problem prove attractive for the narrow-banded systems considered here. Also they have the advantage that they are applicable under the weaker condition A ? ?B ? =? . We have had extensive experience with a method combining aspects of Newton's method and inverse iteration and having a convergence rate of 3.56 . Our implementation combines this basic algorithm with a limiting form of Weilandt deflation to find a sequence of eigenvalues. It has proved extremely satisfactory in a range of applications. This formulation has the further advantage that it is easy to insert the eigenvalue calculation inside an outer loop to satisfy a constraint on an auxiliary parameter. Examples to illustrate both the robustness of the deflation and the flexibility of the approach are provided