Jitse Niesen

Veerle Ledoux, Simon J. A. Malham, Jitse Niesen et al.

We present a numerical method for computing the pure-point spectrum associated with the linear stability of multi-dimensional travelling fronts to parabolic nonlinear systems. Our method is based on the Evans function shooting approach. Transverse to the direction of propagation we project the spectral equations onto a finite Fourier basis. This generates a large, linear, one-dimensional system of equations for the longitudinal Fourier coefficients. We construct the stable and unstable solution subspaces associated with the longitudinal far-field zero boundary conditions, retaining only the information required for matching, by integrating the Riccati equations associated with the underlying Grassmannian manifolds. The Evans function is then the matching condition measuring the linear dependence of the stable and unstable subspaces and thus determines eigenvalues. As a model application, we study the stability of two-dimensional wrinkled front solutions to a cubic autocatalysis model system. We compare our shooting approach with the continuous orthogonalization method of Humpherys and Zumbrun. We then also compare these with standard projection methods that directly project the spectral problem onto a finite multi-dimensional basis satisfying the boundary conditions.

Per Christian Moan, Jitse Niesen

We revisit an algorithm by Skeel et al. for computing the modified, or shadow, energy associated with the symplectic discretization of Hamiltonian systems. By rephrasing the algorithm as a Richardson extrapolation scheme arbitrary high order of accuracy is obtained, and provided error estimates show that it does capture the theoretical exponentially small drift associated with such discretizations. Several numerical examples illustrate the theory.

Shami A Alsallami, Jitse Niesen, Frank W Nijhoff

Modified Hamiltonians are used in the field of geometric numerical integration to show that symplectic schemes for Hamiltonian systems are accurate over long times. For nonlinear systems the series defining the modified Hamiltonian usually diverges. In contrast, this paper constructs and analyzes explicit examples of nonlinear systems where the modified Hamiltonian has a closed-form expression and hence converges. These systems arise from the theory of discrete integrable systems. We present cases of one- and two-degrees symplectic mappings arising as reductions of nonlinear integrable lattice equations, for which the modified Hamiltonians can be computed in closed form. These modified Hamiltonians are also given as power series in the time step by Yoshida's method based on the Baker-Campbell-Hausdorff series. Another example displays an implicit dependence on the time step which could be of relevance to certain implicit schemes in numerical analysis. In the light of these examples, the potential importance of integrable mappings to the field of geometric numerical integration is discussed.

Jitse Niesen, Will M. Wright

We develop an algorithm for computing the solution of a large system of linear ordinary differential equations (ODEs) with polynomial inhomogeneity. This is equivalent to computing the action of a certain matrix function on the vector representing the initial condition. The matrix function is a linear combination of the matrix exponential and other functions related to the exponential (the so-called phi-functions). Such computations are the major computational burden in the implementation of exponential integrators, which can solve general ODEs. Our approach is to compute the action of the matrix function by constructing a Krylov subspace using Arnoldi or Lanczos iteration and projecting the function on this subspace. This is combined with time-stepping to prevent the Krylov subspace from growing too large. The algorithm is fully adaptive: it varies both the size of the time steps and the dimension of the Krylov subspace to reach the required accuracy. We implement this algorithm in the Matlab function phipm and we give instructions on how to obtain and use this function. Various numerical experiments show that the phipm function is often significantly more efficient than the state-of-the-art.

Per Christian Moan, Jitse Niesen

The Magnus series is an infinite series which arises in the study of linear ordinary differential equations. If the series converges, then the matrix exponential of the sum equals the fundamental solution of the differential equation. The question considered in this paper is: When does the series converge? The main result establishes a sufficient condition for convergence, which improves on several earlier results.

Simon Malham, Jitse Niesen

We consider the numerical evaluation of the Evans function, a Wronskian-like determinant that arises in the study of the stability of travelling waves. Constructing the Evans function involves matching the solutions of a linear ordinary differential equation depending on the spectral parameter. The problem becomes stiff as the spectral parameter grows. Consequently, the Gauss--Legendre method has previously been used for such problems; however more recently, methods based on the Magnus expansion have been proposed. Here we extensively examine the stiff regime for a general scalar Schrödinger operator. We show that although the fourth-order Magnus method suffers from order reduction, a fortunate cancellation when computing the Evans matching function means that fourth-order convergence in the end result is preserved. The Gauss--Legendre method does not suffer from order reduction, but it does not experience the cancellation either, and thus it has the same order of convergence in the end result. Finally we discuss the relative merits of both methods as spectral tools.

Jitse Niesen

The Evans function is a tool for assessing the stability of travelling wave solutions for partial differential equations. A recent paper (math.NA/0605581) analyzes the order reduction experienced when evaluating the Evans function numerically. The details of some lengthy calculations were excluded from that paper for clarity. The purpose of this technical report is to make these details publicly available. The report is not intended to be read on its own; the reader is referred to math.NA/0605581 for background and references.