Alexander Ostermann

Alexander Ostermann, Katharina Schratz

We introduce low regularity exponential-type integrators for nonlinear Schrödinger equations for which first-order convergence only requires the boundedness of one additional derivative of the solution. More precisely, we will prove first-order convergence in $H^r$ for solutions in $H^{r+1}$ ($r>d/2$) of the derived schemes. This allows us lower regularity assumptions on the data in the energy space than for instance required for classical splitting or exponential integration schemes. For one dimensional quadratic Schrödinger equations we can even prove first-order convergence without any loss of regularity. Numerical experiments underline the favorable error behavior of the newly introduced exponential-type integrators for low regularity solutions compared to classical splitting and exponential integration schemes.

Marco Caliari, Peter Kandolf, Alexander Ostermann et al.

The Leja method is a polynomial interpolation procedure that can be used to compute matrix functions. In particular, computing the action of the matrix exponential on a given vector is a typical application. This quantity is required, e.g., in exponential integrators. The Leja method essentially depends on three parameters: the scaling parameter, the location of the interpolation points, and the degree of interpolation. We present here a backward error analysis that allows us to determine these three parameters as a function of the prescribed accuracy. Additional aspects that are required for an efficient and reliable implementation are discussed. Numerical examples that illustrate the performance of our Matlab code are included.

Alexander Ostermann, Frédéric Rousset, Katharina Schratz

We present a new filtered low-regularity Fourier integrator for the cubic nonlinear Schrödinger equation based on recent time discretization and filtering techniques. For this new scheme, we perform a rigorous error analysis and establish better convergence rates at low regularity than known for classical schemes in the literature so far. In our error estimates, we combine the better local error properties of the new scheme with a stability analysis based on general discrete Strichartz-type estimates. The latter allow us to handle a much rougher class of solutions as the error analysis can be carried out directly at the level of $L^2$ compared to classical results \black in dimension $d$, \black which are limited to higher-order (sufficiently smooth) Sobolev spaces $H^s$ with $s>d/2$. In particular, we are able to establish a global error estimate in $L^2$ for $H^1$ solutions which is roughly of order $τ^{ {1\over 2} + { 5-d \over 12} }$ in dimension $d \leq 3$ ($τ$ denoting the time discretization parameter). This breaks the "natural order barrier" of $τ^{1/2}$ for $H^1$ solutions which holds for classical numerical schemes (even in combination with suitable filter functions).

Lukas Einkemmer, Alexander Ostermann

A rigorous convergence analysis of the Strang splitting algorithm for Vlasov-type equations in the setting of abstract evolution equations is provided. It is shown that under suitable assumptions the convergence is of second order in the time step τ. As an example, it is verified that the Vlasov-Poisson equation in 1+1 dimensions fits into the framework of this analysis. Also, numerical experiments for the latter case are presented.

Lukas Einkemmer, Alexander Ostermann

Splitting methods constitute a well-established class of numerical schemes for the time integration of partial differential equations. Their main advantages over more traditional schemes are computational efficiency and superior geometric properties. In the presence of non-trivial boundary conditions, however, splitting methods usually suffer from order reduction and some additional loss of accuracy. For diffusion-reaction equations with inhomogeneous oblique boundary conditions, a modification of the classic second-order Strang splitting is proposed that successfully resolves the problem of order reduction. The same correction also improves the accuracy of the classic first-order Lie splitting. The proposed modification only depends on the available boundary data and, in the case of non Dirichlet boundary conditions, on the computed numerical solution. Consequently, this modification can be implemented in an efficient way, which makes the modified splitting schemes superior to their classic versions. The framework employed in our error analysis also allows us to explain the fractional orders of convergence that are often encountered for classic Strang splitting. Numerical experiments that illustrate the theory are provided.

Erwan Faou, Alexander Ostermann, Katharina Schratz

We analyze the convergence of the exponential Lie and exponential Strang splitting applied to inhomogeneous second-order parabolic equations with Dirichlet boundary conditions. A recent result on the smoothing properties of these methods allows us to prove sharp convergence results in the case of homogeneous Dirichlet boundary conditions. When no source term is present and natural regularity assumptions are imposed on the initial value, we show full-order convergence of both methods. For inhomogeneous equations, we prove full-order convergence for the exponential Lie splitting, whereas order reduction to 1.25 for the exponential Strang splitting. Furthermore, we give sufficient conditions on the inhomogeneity for full-order convergence of both methods. Moreover our theoretical convergence results explain the severe order reduction to 0.25 of splitting methods applied to problems involving inhomogeneous Dirichlet boundary conditions. We include numerical experiments to underline the sharpness of our theoretical convergence results.

Lukas Einkemmer, Alexander Ostermann

A rigorous convergence analysis of the Strang splitting algorithm with a discontinuous Galerkin approximation in space for the Vlasov--Poisson equations is provided. It is shown that under suitable assumptions the error is of order $\mathcal{O}(τ^2+h^q +h^q / τ)$, where $τ$ is the size of a time step, $h$ is the cell size, and $q$ the order of the discontinuous Galerkin approximation. In order to investigate the recurrence phenomena for approximations of higher order as well as to compare the algorithm with numerical results already available in the literature a number of numerical simulations are performed.

Marvin Knöller, Alexander Ostermann, Katharina Schratz

Standard numerical integrators suffer from an order reduction when applied to nonlinear Schrödinger equations with low-regularity initial data. For example, standard Strang splitting requires the boundedness of the solution in $H^{r+4}$ in order to be second-order convergent in $H^r$, i.e., it requires the boundedness of four additional derivatives of the solution. We present a new type of integrator that is based on the variation-of-constants formula and makes use of certain resonance based approximations in Fourier space. The latter can be efficiently evaluated by fast Fourier methods. For second-order convergence, the new integrator requires two additional derivatives of the solution in one space dimension, and three derivatives in higher space dimensions. Numerical examples illustrating our convergence results are included. These examples demonstrate the clear advantage of the Fourier integrator over standard Strang splitting for initial data with low regularity.

Hermann Mena, Alexander Ostermann, Lena-Maria Pfurtscheller et al.

The efficient numerical integration of large-scale matrix differential equations is a topical problem in numerical analysis and of great importance in many applications. Standard numerical methods applied to such problems require an unduly amount of computing time and memory, in general. Based on a dynamical low-rank approximation of the solution, a new splitting integrator is proposed for a quite general class of stiff matrix differential equations. This class comprises differential Lyapunov and differential Riccati equations that arise from spatial discretizations of partial differential equations. The proposed integrator handles stiffness in an efficient way, and it preserves the symmetry and positive semidefiniteness of solutions of differential Lyapunov equations. Numerical examples that illustrate the benefits of this new method are given. In particular, numerical results for the efficient simulation of the weather phenomenon El Niño are presented.

Alexander Ostermann, Chiara Piazzola, Hanna Walach

We propose a numerical integrator for determining low-rank approximations to solutions of large-scale matrix differential equations. The considered differential equations are semilinear and stiff. Our method consists of first splitting the differential equation into a stiff and a non-stiff part, respectively, and then following a dynamical low-rank approach. We conduct an error analysis of the proposed procedure, which is independent of the stiffness and robust with respect to possibly small singular values in the approximation matrix. Following the proposed method, we show how to obtain low-rank approximations for differential Lyapunov and for differential Riccati equations. Our theory is illustrated by numerical experiments.

Marlis Hochbruck, Alexander Ostermann

Since their introduction in 1967, Lawson methods have achieved constant interest in the time discretization of evolution equations. The methods were originally devised for the numerical solution of stiff differential equations. Meanwhile, they constitute a well-established class of exponential integrators. The popularity of Lawson methods is in some contrast to the fact that they may have a bad convergence behaviour, since they do not satisfy any of the stiff order conditions. The aim of this paper is to explain this discrepancy. It is shown that non-stiff order conditions together with appropriate regularity assumptions imply high-order convergence of Lawson methods. Note, however, that the term regularity here includes the behaviour of the solution at the boundary. For instance, Lawson methods will behave well in the case of periodic boundary conditions, but they will show a dramatic order reduction for, e.g., Dirichlet boundary conditions. The precise regularity assumptions required for high-order convergence are worked out in this paper and related to the corresponding assumptions for splitting schemes. In contrast to previous work, the analysis is based on expansions of the exact and the numerical solution along the flow of the homogeneous problem. Numerical examples for the Schrödinger equation are included.

Alexander Ostermann, Chunmei Su

We introduce two exponential-type integrators for the "good" Bousinessq equation. They are of orders one and two, respectively, and they require lower regularity of the solution compared to the classical exponential integrators. More precisely, we will prove first-order convergence in Hrfor solutions in H^{r+1} with r > 1/2 for the derived first-order scheme. For the second integrator, we prove second-order convergence in Hrfor solutions in H^{r+3} with r > 1/2 and convergence in L2for solutions in H^3. Numerical results are reported to illustrate the established error estimates. The experiments clearly demonstrate that the new exponential-type integrators are favorable over classical exponential integrators for initial data with low regularity.

Alexander Ostermann, Chunmei Su

We propose an explicit numerical method for the periodic Korteweg-de Vries equation. Our method is based on a Lawson-type exponential integrator for time integration and the Rusanov scheme for Burgers' nonlinearity. We prove first-order convergence in both space and time under a mild Courant-Friedrichs-Lewy condition $τ=O(h)$, where $τ$ and $h$ represent the time step and mesh size, respectively. Numerical examples illustrating our convergence result are given.

Marco Caliari, Alexander Ostermann, Chiara Piazzola

The Schrödinger equation in the presence of an external electromagnetic field is an important problem in computational quantum mechanics. It also provides a nice example of a differential equation whose flow can be split with benefit into three parts. After presenting a splitting approach for three operators with two of them being unbounded, we exemplarily prove first-order convergence of Lie splitting in this framework. The result is then applied to the magnetic Schrödinger equation, which is split into its potential, kinetic and advective parts. The latter requires special treatment in order not to lose the conservation properties of the scheme. We discuss several options. Numerical examples in one, two and three space dimensions show that the method of characteristics coupled with a nonequispaced fast Fourier transform (NFFT) provides a fast and reliable technique for achieving mass conservation at the discrete level.

Lukas Einkemmer, Martina Moccaldi, Alexander Ostermann

Strang splitting is a well established tool for the numerical integration of evolution equations. It allows the application of tailored integrators for different parts of the vector field. However, it is also prone to order reduction in the case of non-trivial boundary conditions. This order reduction can be remedied by correcting the boundary values of the intermediate splitting step. In this paper, three different approaches for constructing such a correction in the case of inhomogeneous Dirichlet, Neumann, and mixed boundary conditions are presented. Numerical examples that illustrate the effectivity and benefits of these corrections are included.

Lukas Einkemmer, Alexander Ostermann

In this paper we consider splitting methods in the presence of non-homogeneous boundary conditions. In particular, we consider the corrections that have been described and analyzed in Einkemmer, Ostermann 2015 and Alonso-Mallo, Cano, Reguera 2016. The latter method is extended to the non-linear case, and a rigorous convergence analysis is provided. We perform numerical simulations for diffusion-reaction, advection-reaction, and dispersion-reaction equations in order to evaluate the relative performance of these two corrections. Furthermore, we introduce an extension of both methods to obtain order three locally and evaluate under what circumstances this is beneficial.

Lukas Einkemmer, Alexander Ostermann

In this paper we propose a method to solve the Kadomtsev--Petviashvili equation based on splitting the linear part of the equation from the nonlinear part. The linear part is treated using FFTs, while the nonlinear part is approximated using a semi-Lagrangian discontinuous Galerkin approach of arbitrary order. We demonstrate the efficiency and accuracy of the numerical method by providing a range of numerical simulations. In particular, we find that our approach can outperform the numerical methods considered in the literature by up to a factor of five. Although we focus on the Kadomtsev--Petviashvili equation in this paper, the proposed numerical scheme can be extended to a range of related models as well.

Andreas Kofler, Tijana Levajković, Hermann Mena et al.

In this paper, we combine deterministic splitting methods with a polynomial chaos expansion method for solving stochastic parabolic evolution problems. The stochastic differential equation is reduced to a system of deterministic equations that we solve explicitly by splitting methods. The method can be applied to a wide class of problems where the related stochastic processes are given uniquely in terms of stochastic polynomials. A comprehensive convergence analysis is provided and numerical experiments validate our approach.

Robert Altmann, Alexander Ostermann

We consider Lie and Strang splitting for the time integration of constrained partial differential equations with a nonlinear reaction term. Since such systems are known to be sensitive with respect to perturbations, the splitting procedure seems promising as we can treat the nonlinearity separately. This has some computational advantages, since we only have to solve a linear constrained system and a nonlinear ODE. However, Strang splitting suffers from order reduction which limits its efficiency. This is caused by the fact that the nonlinear subsystem produces inconsistent initial values for the constrained subsystem. The incorporation of an additional correction term resolves this problem without increasing the computational cost. Numerical examples including a coupled mechanical system illustrate the proven convergence results.

Christian Knapp, Alexander Kendl, Antti Koskela et al.

The equations of motion of a single particle subject to an arbitrary electric and a static magnetic field form a Poisson system. We present a second-order time integration method which preserves well the Poisson structure and compare it to commonly used algorithms, such as the Boris scheme. All the methods are represented in a general framework of splitting methods. We use the so-called $ϕ$ functions, which give efficient ways for both analyzing and implementing the algorithms. Numerical experiments show an excellent long term stability for the new method considered.

Gregor Ehrensperger, Alexander Ostermann, Felix Schwitzer

This paper presents innovative algorithms to efficiently compute erosions and dilations of run-length encoded (RLE) binary images with arbitrary shaped structuring elements. An RLE image is given by a set of runs, where a run is a horizontal concatenation of foreground pixels. The proposed algorithms extract the skeleton of the structuring element and build distance tables of the input image, which are storing the distance to the next background pixel on the left and right hand sides. This information is then used to speed up the calculations of the erosion and dilation operator by enabling the use of techniques which allow to skip the analysis of certain pixels whenever a hit or miss occurs. Additionally the input image gets trimmed during the preprocessing steps on the base of two primitive criteria. Experimental results show the advantages over other algorithms. The source code of our algorithms is available in C++.

Lukas Einkemmer, Alexander Ostermann

For diffusion-reaction equations employing a splitting procedure is attractive as it reduces the computational demand and facilitates a parallel implementation. Moreover, it opens up the possibility to construct second-order integrators that preserve positivity. However, for boundary conditions that are neither periodic nor of homogeneous Dirichlet type order reduction limits its usefulness. In the situation described the Strang splitting procedure is not more accurate than Lie splitting. In this paper, we propose a splitting procedure that, while retaining all the favorable properties of the original method, does not suffer from order reduction. We demonstrate our results by conducting numerical simulations in one and two space dimensions with inhomogeneous and time dependent Dirichlet boundary conditions. In addition, a mathematical rigorous convergence analysis is conducted that confirms the results observed in the numerical simulations.