Simon J. A. Malham

Stochastic Lie Group Integrators

We present Lie group integrators for nonlinear stochastic differential equations with noncommutative vector fields whose solution evolves on a smooth finite-dimensional manifold. Given a Lie group action that generates transport along the manifold, we pull back the stochastic flow on the manifold to the Lie group via the action, and subsequently pull back the flow to the corresponding Lie algebra via the exponential map. We construct an approximation to the stochastic flow in the Lie algebra via closed operations, and then push back to the Lie group and then to the manifold, thus ensuring that our approximation lies in the manifold. We call such schemes stochastic Munthe-Kaas methods after their deterministic counterparts. We also present stochastic Lie group integration schemes based on Castell-Gaines methods. These involve using an underlying ordinary differential integrator to approximate the flow generated by a truncated stochastic exponential Lie series. They become stochastic Lie group integrator schemes if we use Munthe-Kaas methods as the underlying ordinary differential integrator. Further, we show that some Castell-Gaines methods are uniformly more accurate than the corresponding stochastic Taylor schemes. Lastly we demonstrate our methods by simulating the dynamics of a free rigid body such as a satellite and an autonomous underwater vehicle both perturbed by two independent multiplicative stochastic noise processes.

Unsteady fronts in an autocatalytic system

Travelling waves in a model for autocatalytic reactions have, for some parameter regimes, been suggested to have oscillatory instabilities. This instability is confirmed by various methods, including linear stability analysis (exploiting Evans function) and direct numerical simulations. The front instability sets in when the order of the reaction, m, exceeds some threshold, mc(;τ), that depends on the inverse of the Lewis number, τ. The stability boundary, m =mc(;τ), is found numerically for m order one. In the limit m ≫ 1 (in which the system becomes similar to combustion systems with Arrhenius kinetics), the method of matched asymptotic expansions is employed to find the asymptotic front speed and show that mc ∼ (;τ – 1)–1 as ;τ → 1. Just beyond the stability boundary, the unstable rocking of the front saturates supercritically. If the order is increased still further, period–doubling bifurcations occur, and, for small tau there is a transition to chaos through intermittency after the disappearance of a period–4 orbit.

Length scales in solutions of the Navier-Stokes equations

A set of ladder inequalities for the 2d and 3d forced Navier-Stokes equations on a periodic domain (0, L)d is developed, leading to a natural definition of a set of length scales. The authors discuss what happens to these scales if intermittent fluctuations in the vorticity field occur, and they consider how these scales compare to those derived from the attractor dimension and the number of determining modes. Their methods are based on estimates of ratios of norms which appear to play a natural role and which make many of the calculations comparatively easy. In 3d they cannot preclude length scales which are significantly shorter than the Kolmogorov length. In 2d their estimate for a length scale l turns out to be (l/L)-2<or=cG(1+logG)1/2 where G is the Grashof number. This estimate of l is shorter than that derived from the attractor dimension. The reason for this is discussed in detail.

Grassmannian spectral shooting

We present a new numerical method for computing the pure-point spectrum associated with the linear stability of coherent structures. In the context of the Evans function shooting and matching approach, all the relevant information is carried by the flow projected onto the underlying Grassmann manifold. We show how to numerically construct this projected flow in a stable and robust manner. In particular, the method avoids representation singularities by, in practice, choosing the best coordinate patch representation for the flow as it evolves. The method is analytic in the spectral parameter and of complexity bounded by the order of the spectral problem cubed. For large systems it represents a competitive method to those recently developed that are based on continuous orthogonalization. We demonstrate this by comparing the two methods in three applications: Boussinesq solitary waves, autocatalytic travelling waves and the Ekman boundary layer.

Efficient Strong Integrators for Linear Stochastic Systems

We present numerical schemes for the strong solution of linear stochastic differential equations driven by an arbitrary number of Wiener processes. These schemes are based on the Neumann (stochastic Taylor) and Magnus expansions. First, we consider the case when the governing linear diffusion vector fields commute with each other, but not with the linear drift vector field. We prove that numerical methods based on the Magnus expansion are more accurate in the mean-square sense than corresponding stochastic Taylor integration schemes. Second, we derive the maximal rate of convergence for arbitrary multidimensional stochastic integrals approximated by their conditional expectations. Consequently, for general nonlinear stochastic differential equations with noncommuting vector fields, we deduce explicit formulae for the relation between error and computational costs for methods of arbitrary order. Third, we consider the consequences in two numerical studies, one of which is an application arising in stochastic linear-quadratic optimal control.

Computing the Maslov index for large systems

We address the problem of computing the Maslov index for large linear symplectic systems on the real line. The Maslov index measures the signed intersections (with a given reference plane) of a path of Lagrangian planes. The natural chart parameterization for the Grassmannian of La- grangian planes is the space of real symmetric matrices. Linear system evolu- tion induces a Riccati evolution in the chart. For large order systems this is a practical approach as the computational complexity is quadratic in the order. The Riccati solutions, however, also exhibit singularites (which are traversed by changing charts). Our new results involve characterizing these Riccati sin- gularities and two trace formulae for the Maslov index as follows. First, we show that the number of singular eigenvalues of the symmetric chart represen- tation equals the dimension of intersection with the reference plane. Second, the Cayley map is a diffeomorphism from the space of real symmetric matrices to the manifold of unitary symmetric matrices. We show the logarithm of the Cayley map equals the arctan map (modulo 2i) and its trace measures the angle of the Langrangian plane to the reference plane. Third, the Riccati flow under the Cayley map induces a flow in the manifold of unitary symmetric matrices. Using the natural unitary action on this manifold, we pullback the flow to the unitary Lie algebra and monitor its trace. This avoids singularities, and is a natural robust procedure. We demonstrate the effectiveness of these approaches by applying them to a large eigenvalue problem. We also discuss the extension of the Maslov index to the infinite dimensional case.

Computing stability of multidimensional traveling waves

We present a numerical method for computing the pure-point spectrum associated with the linear stability of multidimensional traveling 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 ...

EVALUATING THE EVANS FUNCTION: ORDER REDUCTION IN NUMERICAL METHODS

We consider the numerical evaluation of the Evans function, a Wronskian-like determinant that arises in the study of the stability of travel- ling 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 Schrodinger 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 conver- gence in the end result. Finally we discuss the relative merits of both methods as spectral tools.

Algebraic structure of stochastic expansions and efficient simulation

We investigate the algebraic structure underlying the stochastic Taylor solution expansion for stochastic differential systems. Our motivation is to construct efficient integrators. These are approximations that generate strong numerical integration schemes that are more accurate than the corresponding stochastic Taylor approximation, independent of the governing vector fields and to all orders. The sinhlog integrator introduced by Malham & Wiese (2009) is one example. Herein, we show that the natural context to study stochastic integrators and their properties is the convolution shuffle algebra of endomorphisms; establish a new whole class of efficient integrators; and then prove that, within this class, the sinhlog integrator generates the optimal efficient stochastic integrator at all orders.

Stochastic expansions and Hopf algebras

We study solutions to nonlinear stochastic differential systems driven by a multi-dimensional Wiener process. A useful algorithm for strongly simulating such stochastic systems is the Castell–Gaines method, which is based on the exponential Lie series. When the diffusion vector fields commute, it has been proved that, at low orders, this method is more accurate in the mean-square error than corresponding stochastic Taylor methods. However, it has also been shown that when the diffusion vector fields do not commute, this is not true for strong order one methods. Here, we prove that when there is no drift, and the diffusion vector fields do not commute, the exponential Lie series is usurped by the sinh-log series. In other words, the mean-square error associated with a numerical method based on the sinh-log series is always smaller than the corresponding stochastic Taylor error, in fact to all orders. Our proof uses the underlying Hopf algebra structure of these series, and a two-alphabet associative algebra of shuffle and concatenation operations. We illustrate the benefits of the proposed series in numerical studies.