Search

In this paper, we consider the task of efficiently computing the numerical solution of evolutionary complex Ginzburg--Landau equations on Cartesian product domains with homogeneous Dirichlet/Neumann or periodic boundary conditions. To this aim, we employ for the time integration high-order exponential methods of splitting and Lawson type with constant time step size. These schemes enjoy favorable stability properties and, in particular, do not show restrictions on the time step size due to the underlying stiffness of the models. The needed actions of matrix exponentials are efficiently realized by using a tensor-oriented approach that suitably employs the so-called $\mu$-mode product (when the semidiscretization in space is performed with finite differences) or with pointwise operations in Fourier space (when the model is considered with periodic boundary conditions). The overall effectiveness of the approach is demonstrated by running simulations on a variety of two- and three-dimensional (systems of) complex Ginzburg--Landau equations with cubic or cubic-quintic nonlinearities, which are widely considered in literature to model relevant physical phenomena. In fact, we show that high-order exponential-type schemes may outperform standard techniques to integrate in time the models under consideration, i.e., the well-known second-order split-step method and the explicit fourth-order Runge--Kutta integrator, for stringent accuracies.

We propose a second order exponential scheme suitable for two-component coupled systems of stiff evolutionary advection--diffusion--reaction equations in two and three space dimensions. It is based on a directional splitting of the involved matrix functions, which allows for a simple yet efficient implementation through the computation of small-sized exponential-like functions and tensor-matrix products. The procedure straightforwardly extends to the case of an arbitrary number of components and to any space dimension. Several numerical examples in 2D and 3D with physically relevant (advective) Schnakenberg, FitzHugh--Nagumo, DIB, and advective Brusselator models clearly show the advantage of the approach against state-of-the-art techniques.

In this manuscript, we propose an efficient, practical and easy-to-implement way to approximate actions of $\varphi$-functions for matrices with $d$-dimensional Kronecker sum structure in the context of exponential integrators up to second order. The method is based on a direction splitting of the involved matrix functions, which lets us exploit the highly efficient level 3 BLAS for the actual computation of the required actions in a $\mu$-mode fashion. The approach has been successfully tested on two- and three-dimensional problems with various exponential integrators, resulting in a consistent speedup with respect to a technique designed to compute actions of $\varphi$-functions for Kronecker sums.

In this paper we consider an approach to improve the performance of exponential Runge--Kutta integrators and Lawson schemes} in cases where the solution of a related, but usually much simpler, problem can be computed efficiently. While for implicit methods such an approach is common (e.g. by using preconditioners), for exponential integrators this has proven more challenging. Here we propose to extract a constant coefficient differential operator from the semilinear advection-diffusion-reaction equation for which, in many situations, efficient methods are known to compute the required matrix functions. Both a linear stability analysis and {\color{black} extensive} numerical experiments show that the resulting schemes can be unconditionally stable. In fact, we find that exponential integrators of Runge--Kutta type and Lawson schemes can have better stability properties than similarly constructed implicit-explicit schemes. We also derive two new Lawson type integrators that further improve on these stability properties. The overall effectiveness of the approach is highlighted by a number of performance comparisons on examples in two and three space dimensions.

In this paper we consider mean-field optimal control problems with selective action of the control, where the constraint is a continuity equation involving a non-local term and diffusion. First order optimality conditions are formally derived in a general framework, accounting for boundary conditions. Hence, the optimality system is used to construct a reduced gradient method, where we introduce a novel algorithm for the numerical realization of the forward and the backward equations, based on exponential integrators. We illustrate extensive numerical experiments on different control problems for collective motion in the context of opinion formation and pedestrian dynamics.

We present a method for computing actions of the exponential-like $\varphi$-functions for a Kronecker sum $K$ of $d$ arbitrary matrices $A_\mu$. It is based on the approximation of the integral representation of the $\varphi$-functions by Gaussian quadrature formulas combined with a scaling and squaring technique. The resulting algorithm, which we call PHIKS, evaluates the required actions by means of $\mu$-mode products involving exponentials of the small sized matrices $A_\mu$, without forming the large sized matrix $K$ itself. PHIKS, which profits from the highly efficient level 3 BLAS, is designed to compute different $\varphi$-functions applied on the same vector or a linear combination of actions of $\varphi$-functions applied on different vectors. In addition, thanks to the underlying scaling and squaring techniques, the desired quantities are available simultaneously at suitable time scales. All these features allow the effective usage of PHIKS in the exponential integration context. In fact, our newly designed method has been tested on popular exponential Runge--Kutta integrators of stiff order from one to four, in comparison with state-of-the-art algorithms for computing actions of $\varphi$-functions. The numerical experiments with discretized semilinear evolutionary 2D or 3D advection--diffusion--reaction, Allen--Cahn, and Brusselator equations show the superiority of the proposed $\mu$-mode approach.

In this manuscript, we present a common tensor framework which can be used to generalize one-dimensional numerical tasks to arbitrary dimension $d$ by means of tensor product formulas. This is useful, for example, in the context of multivariate interpolation, multidimensional function approximation using pseudospectral expansions and solution of stiff differential equations on tensor product domains. The key point to obtain an efficient-to-implement BLAS formulation consists in the suitable usage of the $\mu$-mode product (also known as tensor-matrix product or mode- $n$ product) and related operations, such as the Tucker operator. Their MathWorks MATLAB/GNU Octave implementations are discussed in the paper, and collected in the package KronPACK. We present numerical results on experiments up to dimension six from different fields of numerical analysis, which show the effectiveness of the approach.

In this paper, we propose a $\mu$-mode integrator for computing the solution of stiff evolution equations. The integrator is based on a $d$-dimensional splitting approach and uses exact (usually precomputed) one-dimensional matrix exponentials. We show that the action of the exponentials, i.e. the corresponding batched matrix-vector products, can be implemented efficiently on modern computer systems. We further explain how $\mu$-mode products can be used to compute spectral transforms efficiently even if no fast transform is available. We illustrate the performance of the new integrator by solving, among the others, three-dimensional linear and nonlinear Schr\"odinger equations, and we show that the $\mu$-mode integrator can significantly outperform numerical methods well established in the field. We also discuss how to efficiently implement this integrator on both multi-core CPUs and GPUs. Finally, the numerical experiments show that using GPUs results in performance improvements between a factor of $10$ and $20$, depending on the problem.

Rational exponential integrators (REXI) are a class of numerical methods that are well suited for the time integration of linear partial differential equations with imaginary eigenvalues. Since these methods can be parallelized in time (in addition to the spatial parallelization that is commonly performed) they are well suited to exploit modern high performance computing systems. In this paper, we propose a novel REXI scheme that drastically improves accuracy and efficiency. The chosen approach will also allow us to easily determine how many terms are required in the approximation in order to obtain accurate results. We provide comparative numerical simulations for a shallow water equation that highlight the efficiency of our approach and demonstrate that REXI schemes can be efficiently implemented on graphic processing units., Comment: 25 pages, 5 figures, submitted to Journal of Computational Physics

We present an anisotropic extension of the isotropic osmosis model that has been introduced by Weickert et al.~(Weickert, 2013) for visual computing applications, and we adapt it specifically to shadow removal applications. We show that in the integrable setting, linear anisotropic osmosis minimises an energy that involves a suitable quadratic form which models local directional structures. In our shadow removal applications we estimate the local structure via a modified tensor voting approach (Moreno, 2012) and use this information within an anisotropic diffusion inpainting that resembles edge-enhancing anisotropic diffusion inpainting (Weickert, 2006, Gali\'c, 2008). Our numerical scheme combines the nonnegativity preserving stencil of Fehrenbach and Mirebeau (Fehrenbach, 2014) with an exact time stepping based on highly accurate polynomial approximations of the matrix exponential. The resulting anisotropic model is tested on several synthetic and natural images corrupted by constant shadows. We show that it outperforms isotropic osmosis, since it does not suffer from blurring artefacts at the shadow boundaries.

Although Fourier series approximation is ubiquitous in computational physics owing to the Fast Fourier Transform (FFT) algorithm, efficient techniques for the fast evaluation of a three-dimensional truncated Fourier series at a set of \emph{arbitrary} points are quite rare, especially in MATLAB language. Here we employ the Nonequispaced Fast Fourier Transform (NFFT, by J. Keiner, S. Kunis, and D. Potts), a C library designed for this purpose, and provide a Matlab and GNU Octave interface that makes NFFT easily available to the Numerical Analysis community. We test the effectiveness of our package in the framework of quantum vortex reconnections, where pseudospectral Fourier methods are commonly used and local high resolution is required in the post-processing stage. We show that the efficient evaluation of a truncated Fourier series at arbitrary points provides excellent results at a computational cost much smaller than carrying out a numerical simulation of the problem on a sufficiently fine regular grid that can reproduce comparable details of the reconnecting vortices.

The Schr\"odinger 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\"odinger 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.

We extensively study the numerical accuracy of the well-known time splitting Fourier spectral method for the approximation of singular solutions of the Gross-Pitaevskii equation. In particular, we explore its capability of preserving a steady-state vortex solution, whose density profile is approximated by a very accurate diagonal Pad\'e expansion of order 8, here explicitly derived for the first time. Although the Fourier spectral method turns out to be only slightly more accurate than a time splitting finite difference scheme, the former is reliable and efficient. Moreover, at a post-processing stage, it allows an accurate evaluation of the solution outside grid points, thus becoming particularly appealing when high resolution is needed, such as in the study of quantum vortex interactions.

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.

We compute the first eigenpair for variable exponent eigenvalue problems. We compare the homogeneous definition of first eigenvalue with previous nonhomogeneous notions in the literature. We highlight the symmetry breaking phenomena

By virtue of numerical arguments we study a bifurcation phenomenon occurring for a class of minimization problems associated with the quasi-linear Schrodinger equation., Comment: 9 pages

We provide some numerical computations for the soliton dynamics of the nonlinear Schr\"odinger equation with an external potential. After computing the ground state solution $r$ of a related elliptic equation we show that, in the semi-classical regime, the center of mass of the solution with initial datum modelled on $r$ is driven by the solution of a Newtonian type law. Finally, we provide some examples and analyze the numerical errors in the two dimensional case when $V$ is an harmonic potential., Comment: 13 pages, 4 figures

We consider a system of Gross-Pitaevskii equations in R^2 modelling a mixture of two Bose-Einstein condensates with repulsive interaction. We aim to study the qualitative behaviour of ground and excited state solutions. We allow two different harmonic and off-centered trapping potentials and study the spatial patterns of the solutions within the Thomas-Fermi approximation as well as phase segregation phenomena within the large-interaction regime., Comment: 21 pages, to appear in Dyn. Partial Diff. Equat