Instability of the finite‐difference split‐step method applied to the nonlinear Schrödinger equation. I. standing soliton

We consider numerical instability that can be observed in simulations of solitons of the nonlinear Schrödinger equation (NLS) by a split‐step method (SSM) where the linear part of the evolution is solved by a finite‐difference discretization. The von Neumann analysis predicts that this method is unconditionally stable on the background of a constant‐amplitude plane wave. However, simulations show that the method can become unstable on the background of a soliton. We present an analysis explaining this instability. Both this analysis and the features and threshold of the instability are substantially different from those of the Fourier SSM, which computes the linear part of the NLS by a spectral discretization. For example, the modes responsible for the numerical instability are not similar to plane waves, as for the Fourier SSM or, more generally, in the von Neumann analysis. Instead, they are localized at the sides of the soliton. This also makes them different from “physical” (as opposed to numerical) unstable modes of nonlinear waves, which (the modes) are localized around the “core” of a solitary wave. Moreover, the instability threshold for thefinite‐difference split‐step method is considerably relaxed compared with that for the Fourier split‐step. © 2015 Wiley Periodicals, Inc. Numer Methods Partial Differential Eq 32: 1002–1023, 2016     

Mathematical and numerical treatment of instabilities of non-axisymmetric confined vortices under the Dirichlet boundary conditions

This paper is concerned with an efficient analytical and numerical method for computing instability bounds of non-axisymmetric confined vortices. We propose a new approach based on a recently mathematical method of meshless collocation, when the shifted orthogonal Chebyshev base is adapted to satisfy the boundary conditions. Appropriately modified Tollmien–Schlichting modes used for the linear stability analysis offer the main advantage of dealing with a steady perturbation model of a confined vortex that leads to a variety of information describing the behavior of the full three-dimensional unsteady flow.     

Dark solitons in external potentials

We consider the persistence and stability of dark solitons in the Gross–Pitaevskii (GP) equation with a small decaying potential. We show that families of black solitons with zero speed originate from extremal points of an appropriately defined effective potential and persist for sufficiently small strength of the potential. We prove that families at the maximum points are generally unstable with exactly one real positive eigenvalue, while families at the minimum points are generally unstable with exactly two complex-conjugated eigenvalues with positive real part. This mechanism of destabilization of the black soliton is confirmed in numerical approximations of eigenvalues of the linearized GP equation and full numerical simulations of the nonlinear GP equation. We illustrate the monotonic instability associated with the real eigenvalues and the oscillatory instability associated with the complex eigenvalues and compare the numerical results of evolution of a dark soliton with the predictions of Newton’s particle law for its position.      

Dark solitons in external potentials

We consider the persistence and stability of dark solitons in the Gross–Pitaevskii (GP) equation with a small decaying potential. We show that families of black solitons with zero speed originate from extremal points of an appropriately defined effective potential and persist for sufficiently small strength of the potential. We prove that families at the maximum points are generally unstable with exactly one real positive eigenvalue, while families at the minimum points are generally unstable with exactly two complex-conjugated eigenvalues with positive real part. This mechanism of destabilization of the black soliton is confirmed in numerical approximations of eigenvalues of the linearized GP equation and full numerical simulations of the nonlinear GP equation. We illustrate the monotonic instability associated with the real eigenvalues and the oscillatory instability associated with the complex eigenvalues and compare the numerical results of evolution of a dark soliton with the predictions of Newton’s particle law for its position.     

On the Boyd–Kadomtsev system for the three-wave coupling problem

The three-wave coupling system is widely used in plasma physics, specially for the Brillouin instability simulations. We study here a related system obtained with an infinite speed of light. After showing that it is well posed, we propose a numerical method which is based on an implicit time discretization. This method is illustrated on test cases and an extension to the problem with finite speed of light is proposed. To cite this article: R. Sentis, C. R. Acad. Sci. Paris, Ser. I 347 (2009).     

Instability analysis of the split‐step Fourier method on the background of a soliton of the nonlinear Schrödinger equation

We analyze a numerical instability that occurs in the well‐known split‐step Fourier method on the background of a soliton. This instability is found to be very sensitive to small changes of the parameters of both the numerical grid and the soliton, unlike the instability of most finite‐difference schemes. Moreover, the principle of “frozen coefficients,” in which variable coefficients are treated as “locally constant” for the purpose of stability analysis, is strongly violated for the instability of the split‐step method on the soliton background. Our analysis quantitatively explains all these features. It is enabled by the fact that the period of oscillations of the unstable Fourier modes is much smaller than the width of the soliton. Our analysis is different from the von Neumann analysis in that it requires spatially growing or decaying harmonics (not localized near the boundaries) as opposed to purely oscillatory ones. © 2010 Wiley Periodicals, Inc. Numer Methods Partial Differential Eq 28: 641–669, 2012     

Solution of sparse linear least squares problems using givens rotations

We describe a direct method for solving sparse linear least squares problems. The storage required for the method is no more than that needed for the conventional normal equations approach. However, the normal equations are not computed; orthogonal transformations are applied to the coefficient matrix, thus avoiding the potential numerical instability associated with computing the normal equations. Our approach allows full exploitation of sparsity, and permits the use of a fixed (static) data structure during the numerical computation. Finally, the method processes the coefficient matrix one row at a time, allowing for the convenient use of auxiliary storage and updating operations.     

Numerical Challenges for Resolving Spike Dynamics for Two One-Dimensional Reaction-Diffusion Systems

Asymptotic and numerical methods are used to highlight different types of dynamical behaviors that occur for the motion of a localized spike-type solution to the singularly perturbed Gierer–Meinhardt and Schnakenberg reaction-diffusion models in a one-dimensional spatial domain. Depending on the parameter range in these models, there can either be a slow evolution of a spike toward the midpoint of the domain, a sudden oscillatory instability triggered by a Hopf bifurcation leading to an intricate temporal oscillation in the height of the spike, or a pulse-splitting instability leading to the creation of new spikes in the domain. Criteria for the onset of these oscillatory and pulse-splitting instabilities are obtained through asymptotic and numerical techniques. A moving-mesh numerical method is introduced to compute these different behaviors numerically, and results are compared with corresponding results computed using a method of lines based software package.     

A control type method for solving the Cauchy–Stokes problem

In this work, we propose an approximate optimal control formulation of the Cauchy problem for the Stokes system. Here the problem is converted into an optimization one. In order to handle the instability of the solution of this ill-posed problem, a regularization technique is developed. We add a term in the least square function which happens to vanish while the algorithm converges. The efficiency of the proposed method is illustrated by numerical experiments.     

On the structure and growth rate of unstable modes to the Rossby‐Haurwitz wave

The normal mode instability study of a steady Rossby‐Haurwitz wave is considered both theoretically and numerically. This wave is exact solution of the nonlinear barotropic vorticity equation describing the dynamics of an ideal fluid on a rotating sphere, as well as the large‐scale barotropic dynamics of the atmosphere. In this connection, the stability of the Rossby‐Haurwitz wave is of considerable mathematical and meteorological interest. The structure of the spectrum of the linearized operator in case of an ideal fluid is studied. A conservation law for perturbations to the Rossby‐Haurwitz wave is obtained and used to get a necessary condition for its exponential instability. The maximum growth rate of unstable modes is estimated. The orthogonality of the amplitude of a non‐neutral or non‐stationary mode to the Rossby‐Haurwitz wave is shown in two different inner products. The analytical results obtained are used to test and discuss the accuracy of a numerical spectral method used for the normal mode stability study of arbitrary flow on a sphere. The comparison of the numerical and theoretical results shows that the numerical instability study method works well in case of such smooth solutions as the zonal flows and Rossby‐Haurwitz waves. © 2004 Wiley Periodicals, Inc. Numer Methods Partial Differential Eq, 2005     

Thermal instability of a layered viscoelastic rectangular prism under high-frequency shear loading

Within the framework of a coupled problem of thermoviscoelasticity, using the method of numerical modeling, we have studied the thermal instability of a rectangular prism, composed of copper and polyethylene or polymethylmethacrylate layers, in the course of its dissipative heating. The prism is subjected to high-frequency force or kinematic shear. We have established that, in the case of polyethylene, thermal instability occurs under conditions of force loading and is absent under a kinematic one. However, for polymethylmethacrylate, thermal instability takes place for both cases of loading. This effect is attributable to the existence of temperature intervals where the shear and volume loss compliances for each of the polymers increases with the temperature. We have also revealed that the critical values of thermal instability for a prism with metal layers are substantially higher than those for a homogeneous prism. In addition to thermal instability, thermal resonance instability is also possible. It is caused by the jump of the thermal state from the low- to high-temperature branch of the soft-type resonance characteristic.     

Spike patterns in a reaction–diffusion ODE model with Turing instability

We explore a mechanism of pattern formation arising in processes described by a system of a single reaction–diffusion equation coupled with ordinary differential equations. Such systems of equations arise from the modeling of interactions between cellular processes and diffusing growth factors. We focus on the model of early carcinogenesis proposed by Marciniak‐Czochra and Kimmel, which is an example of a wider class of pattern formation models with an autocatalytic non‐diffusing component. We present a numerical study showing emergence of periodic and irregular spike patterns because of diffusion‐driven instability. To control the accuracy of simulations, we develop a numerical code on the basis of the finite‐element method and adaptive mesh grid. Simulations, supplemented by numerical analysis, indicate a novel pattern formation phenomenon on the basis of the emergence of nonstationary structures tending asymptotically to a sum of Dirac deltas. Copyright © 2013 John Wiley & Sons, Ltd.     

Three dimensional electromagnetic scattering T-matrix computations

The infinite T-matrix method is a powerful tool for electromagnetic scattering simulations, particularly when one is interested in changes in orientation of the scatterer with respect to the incident wave or changes of configuration of multiple scatterers and random particles, because it avoids the need to solve the fully reconfigured systems. The truncated T-matrix (for each scatterer in an ensemble) is often computed using the null-field method. The main disadvantage of the null-field based T-matrix computation is its numerical instability for particles that deviate from a sphere. For large and/or highly non-spherical particles, the null-field method based truncated T-matrix computations can become slowly convergent or even divergent. In this work, we describe an electromagnetic scattering surface integral formulation for T-matrix computations that avoids the numerical instability. The new method is based on a recently developed high-order surface integral equation algorithm for far field computations using basis functions that are tangential on a chosen non-spherical obstacle. The main focus of this work is on the mathematical details required to apply the high-order algorithm to compute a truncated T-matrix that describes the scattering properties of a chosen perfect conductor in a homogeneous medium. We numerically demonstrate the stability and convergence of the T-matrix computations for various perfect conductors using plane wave incident radiation at several low to medium frequencies and simulation of the associated radar cross of the obstacles.     

On the effect of electron-mass variation in numerical simulations of the Farley-Buneman instability

The Farley-Buneman instability is the plasma instability in the E region of the Earth’s ionosphere. Many studies on instability simulations use algorithms based on the particle method, which has many shortcomings. In particular, the solutions obtained with this method involve numerical noises owing to a finite number of the particles, since the electron mass is increased in order to reduce the computational complexity of the algorithm. In this study, the effect of an increase in the electron mass on the qualitative and quantitative characteristics of the instability development process are considered. For this purpose, a software package developed on the basis of a mathematical model consisting of the Poisson equation for the electric potential, the fluid equation for electrons, and the kinetic equation for ions is used. The equations are solved numerically, which makes it possible to avoid the problems inherent in the particle method. An important result obtained in this paper is that even a slight increase in the electron mass leads to a considerable variation of the instability development parameters. This variation manifests itself as an increase in the wavelength of plasma fluctuations and a decrease in the strength of the turbulent electric field.     

Stability of periodic traveling flexural‐gravity waves in two dimensions

In this work, we solve the Euler's equations for periodic waves traveling under a sheet of ice. These waves are referred to as flexural‐gravity waves. We compare and contrast two models for the effect of the ice: a linear model and a nonlinear model. The benefit of this reformulation is that it facilitates the asymptotic analysis. We use it to derive the nonlinear Schrödinger equation that describes the modulational instability of periodic traveling waves. We compare this asymptotic result with the numerical computation of stability using the Fourier–Floquet–Hill method to show they agree qualitatively. We show that different models have different stability regimes for large values of the flexural rigidity parameter. Numerical computations are also used to analyze high‐frequency instabilities in addition to the modulational instability. In the regions examined, these are shown to be the same regardless of the model representing ice.     

Barotropic instability of shear flows

We consider barotropic instability of shear flows for incompressible fluids with Coriolis effects. For a class of shear flows, we develop a new method to find the sharp stability conditions. We study the flow with Sinus profile in details and obtain the sharp stability boundary in the whole parameter space, which corrects previous results in the fluid literature. Our new results are confirmed by more accurate numerical computation. The addition of the Coriolis force is found to bring fundamental changes to the stability of shear flows. Moreover, we study dynamical behaviors near the shear flows, including the bifurcation of nontrivial traveling wave solutions and the linear inviscid damping. The first ingredient of our proof is a careful classification of the neutral modes. The second one is to write the linearized fluid equation in a Hamiltonian form and then use an instability index theory for general Hamiltonian partial differential equations. The last one is to study the singular and nonresonant neutral modes using Sturm-Liouville theory and hypergeometric functions.     

Instability of solitary waves for a generalized Benney–Luke equation

We study linear instability of solitary wave solutions of a one-dimensional generalized Benney–Luke equation, which is a formally valid approximation for describing two-way water wave propagation in the presence of surface tension. Further, we implement a finite difference numerical scheme which combines an explicit predictor and an implicit corrector step to compute solutions of the model equation which is used to validate the theory presented.     

Recovering a tsunami source and designing an observational system based on an <Emphasis Type="Italic">r</Emphasis>-solution method

A theoretical approach to recovering the initial tsunami waveform proposed by the author is applied to a real tsunami event. The approach is based on least squares and truncated singular value decomposition techniques. The problem of recovering a tsunami source from available measurements of the incomingwave is considered as an inverse problem ofmathematical physics. To avoid potential instability of the numerical solution, an r-solution method is used. An observation system for retrieving the initial tsunami waveform is optimally configured by using this approach. Several numerical simulations with real tsunami data are presented.     

一个拟极小化左共轭梯度算法及其数值表现

左共轭梯度法是求解大型稀疏线性方程组的一种新兴的Krylov子空间方法.为克服该算法数值表现不稳定、迭代中断的缺点,本文对原方法进行等价变形,得到左共轭梯度方向的另一迭代格式,给出一个拟极小化左共轭梯度算法.数值结果证实了该变形算法与原算法的相关性.