A system of ODEs for a Perturbation of a Minimal Mass Soliton

A system of ODEs for a Perturbation of a Minimal Mass Soliton

Jeremy Marzuola marzuola@math.uni-bonn.de Sarah Raynor raynorsg@math.wfu.edu  and  Gideon Simpson simpson@math.toronto.edu Mathematics Institute, Bonn University
Endenicher Allee 60, D-53115 Bonn, Germany
Mathematics Department, Wake Forest University
P.O. Box 7388, 127 Manchester Hall, Winston-Salem, NC, 27109 USA
Mathematics Department, University of Toronto
40 St. George St., Toronto, Ontario, Canada M5S 2E4

We study soliton solutions to a nonlinear Schrödinger equation with a saturated nonlinearity. Such nonlinearities are known to possess minimal mass soliton solutions. We consider a small perturbation of a minimal mass soliton, and identify a system of ODEs similar to those from [16], which model the behavior of the perturbation for short times. We then provide numerical evidence that under this system of ODEs there are two possible dynamical outcomes, which is in accord with the conclusions of [31]. For initial data which supports a soliton structure, a generic initial perturbation oscillates around the stable family of solitons. For initial data which is expected to disperse, the finite dimensional dynamics follow the unstable portion of the soliton curve.

1. Introduction

We consider the initial value problem for the nonlinear Schrödinger equation (NLS) in :


where the nonlinearity is a saturated nonlinearity of the form


where for and for . For large, (1.1) behaves as though it were subcritical while for small, it behaves as though it were supercritical. This guarantees both existence of soliton solutions and global well-posedness in .

For our purposes, must be chosen substantially larger than the critical exponent, , in order to allow sufficient regularity when linearizing the equation. For our numerical analysis, we work in one spatial dimension, with the specific nonlinearity


The equation (1.1) is globally well-posed in with the usual norm

where is the usual Sobolev space with norm

and is the weighted Sobolev space with norm

This is commonly referred to as the space with finite variance. The global well-posedness initial data in follows from the standard well-posedness theory for semilinear Schrödinger equations. Additionally, we assume that is spherically symmetric, which implies is also spherically symmetric for all . Proofs can be found in numerous references including [15] and [43].

A soliton solution of (1.1) is a function of the form


where and is a positive, spherically symmetric, exponentially decaying solution of the equation:


For our particular nonlinearity, for any there is a unique solitary wave solution to (1.5), see [8] and [29].

For large the solitons are stable, while for small they are unstable. A precise stability criterion identifying stable and unstable regions is provided in [22] and [36], generalizing earlier work on stability in [44], [45]. This amounts to examining the relation , defining a soliton curve. Where it is increasing(decreasing) as a function of , the solitons are stable(unstable). Several such curves appear in Figure 1.

Figure 1. Plots of the soliton curves ( with respect to ) for a subcritical nonlinearity, critical nonlinearity, supercritical nonlinearity, and the saturated nonlinearity (1.3). The curves for the monomial nonlinearities are found analytically, while the curves for the saturated nonlinearities are found numerically using the method discussed in Section 4.1. All are in . Here , as the supercritical and saturated cases diverge as .

As can be seen numerically in Figure 1, the nonlinearity spawns a soliton of minimal mass. Though certain asymptotic methods can be used to describe the increasing nature of the curve as (multiscale methods) and (variational methods), we forego an analytic description of the soliton curve and focus on the minimal mass soliton shown to exist in the numerical plot. In [16], Comech and Pelinovsky demonstrated that the minimal mass soliton possesses a fundamentally nonlinear instability. They accomplished this by finding a small perturbation that forces the solution a fixed distance away from the minimal mass soliton in finite time. Their technique reduces to studying an ODE modeling the perturbation for short times. For appropriate data, the ODE is unstable.

We conjecture that though the minimal mass soliton may be unstable on short time scales, on a longer time scale the solution will ultimately relax to the stable branch of the soliton curve. This conjecture is part of a larger conjecture that solutions which do not disperse as must eventually converge towards the stable portion of the soliton curve. For nonlinearities with a specific two power structure, dynamics of this type were observed by Pelinovsky, Afanasjev, and Kivshar, who modelled the behavior of solutions near a minimal mass soliton by a second order ODE via adiabatic expansion in [31]. By contrast, our method uses the full dynamical system of modulation parameters to find a -dimensional system of ODEs which is structured to allow for eventual recoupling to the continuous spectrum. This conjecture has also been explored numerically by Buslaev and Grikurov for two power nonlinearities in [12], where they found that a solution which is initially a perturbation of an unstable soliton tends to approach and then oscillate around a stable soliton.

The purpose of this work is to numerically explore this conjecture. Following [16], we break the perturbation into the discrete and continuous parts relative to the linearization of the Schrödinger operator around the soliton. The discrete portion yields a four dimensional system of nonlinear ODEs. We further simplify the system expanding the equations in powers of the dependent variables and dropping cubic and higher terms.

An obstacle in studying these ODEs is that the signs and magnitudes of the coefficients are not self-evident, necessitating numerical methods. We compute these numbers, which are intimately related to the minimal mass soliton, using the spectral method. The use of the function for numerically solving differential equations dates to Stenger [39]. It has been successfully used in a wide variety of linear and nonlinear, time dependent and independent, differential equations, [4, 7, 10, 14, 19, 25, 32, 42]. In this work, we first numerically solve (1.5) for the soliton as a nonlinear collocation problem. We then use this information to compute the generalized kernel of the operator after linearization about the soliton.

With these coefficients in hand, we numerically integrate the ODE system, plotting the results. We find that there are two different types of behavior for the finite dimensional system, depending on the initial data. If the initial data represents a solution which our nonlinear solver indicates can support a soliton, then we find that the solution is oscillatory. It is initially attracted to the stable side of the curve, and, over intermediate time scales, proceeds to oscillate around the minimal mass soliton. If we initialize with this type of data but with the unstable conditions found in [16], the ODEs initially move in the unstable direction but quickly reverse, before commencing oscillation. On the other hand, if we begin with initial conditions which are expected to disperse as , our data indicate that the finite dimensional dynamics push the solution along the unstable soliton curve towards the value rather quickly. This solution matches well to the solution for (1.1) with corresponding initial data for as long as the mass conservation of the solution allows, after which our model continues to follow the unstable soliton curve but the actual solution disperses. In [31], the authors observed similar dynamics, with both oscillatory and dispersive regimes.

These ODEs are an approximation valid on a short time interval. This study is the beginning of an analysis to show that perturbations of the minimal mass soliton are attracted to the stable side of the soliton curve. In a forthcoming work we hope to show how the continuous-spectrum part of the perturbation interacts with the discrete-spectrum perturbation. Based on the work of Soffer and Weinstein,[38] we expect coupling to the continuous spectrum to cause radiation damping, which will ultimately cause the solution to have damped oscillations and select a soliton on the stable side of the curve.

This paper is organized as follows. In section 2, we introduce preliminaries and necessary definitions. In section 3, we derive the system of ODEs. In section 4, we explain our numerical methods for finding the coefficients of the ODEs. In section 5, we show the numerical solutions of the ODEs and explain our results. Finally, in section 6 we present our conclusions and plans for future work. An appendix contains details of our numerical method for computation of the soliton and related coefficients.

Acknowledgments This project began out of a conversation with Catherine Sulem and JM. JM was partially funded by an NSF Postdoc at Columbia University and a Hausdorff Center Postdoc at the University of Bonn. In addition, JM would like to thank the University of North Carolina, Chapel Hill for graciously hosting him during part of this work. SR would like to thank the University of Chicago for their hospitality while some of this work was completed. GS was funded in part by NSERC. In addition, the authors wish to thank Dmitry Pelinovksy, Mary Pugh, Catherine Sulem, and Michael Weinstein for helpful comments and suggestions.

2. Definitions and Setup

For data , there are several conserved quantities. Particularly important invariants are:

Conservation of Mass (or Charge):
Conservation of Energy:


Detailed proofs of these conservation laws can be easily arrived at by using energy estimates or Noether’s Theorem, which relates conservation laws to symmetries of an equation. See [43] for details.

With this type of nonlinearity, it is known that soliton solutions to NLS exist and are unique. Existence of solitary waves for nonlinearities of the type (1.2) is proved by in [8] in using ODE techniques and in higher dimensions by minimizing the functional

with respect to the functional

Then, using a minimizing sequence and Schwarz symmetrization, one infers the existence of the nonnegative, spherically symmetric, decreasing soliton solution. Once we know that minimizers are radially symmetric, uniqueness can be established via a shooting method, showing that the desired soliton occurs at only one initial value, [29].

Of great importance is the fact that and are differentiable with respect to . This can be determined from the works of Shatah, namely [34], [35]. Differentiating (1.5), and all with respect to , we have the relation

Numerics show that if we plot with respect to for the saturated nonlinearity, the soliton curve goes to as goes to or and has a global minimum at some ; see Figure 1. This will be explored in detail in a subsequent numerical work by Marzuola [27].

We are interested in the stability of these explicit solutions under perturbations of the initial data.

Definition 2.1.

The soliton is said to be orbitally stable if, , such that, for any initial data such that , for any , there is some such that .

Definition 2.2.

The soliton is said to be asymptotically stable, if, such that if , then for large , such that disperses as a solution to the corresponding linear problem would.

Variational techniques developed in [44], [45] and generalized to an abstract setting in [22] and [36] tell us that when is convex, or , we are guaranteed stability under small perturbations, while for we are guaranteed that the soliton is unstable under small perturbations. For brief reference on this subject, see Chapter 4 of [43]. For nonlinearities that are twice differentiable at the origin and of monomial type at infinity (which would include our saturated nonlinearities), asymptotic stability has been studied for a finite collection of strongly orbitally stable solitons by Buslaev and Perelman[13], Cuccagna[17], and Rodnianski, Soffer and Schlag[33].

At a minimum of , soliton instability is more subtle, because it is due solely to nonlinear effects. See [16], where this purely nonlinear instability is proved to occur by reducing the behavior of the discrete part of the spectrum to an ODE that is unstable for certain initial conditions.

2.1. Linearization about a Soliton

Throughout this section, we use vector notation, , to represent complex functions. Any function written without vector notation is assumed to be real. For example, the complex valued scalar function will be written . In this notation, multiplication by is represented by the matrix . We denote by the complex vector , where is the real profile of the soliton with parameter . For simplicity, we suppress the subscript, writing in place of .

For later reference, we now explicitly characterize the linearization of NLS about a soliton solution. First consider the linear evolution of the perturbation of a soliton via the ansatz:


with . For the purposes of finding the linearized hamiltonian at we do not need to allow the parameters and to modulate, but when we develop our full system of equations in Section 3 parameter modulation will be taken into account. Inserting (2.1) into the equation we know that since is a soliton solution we have


(This calculation is explained in more detailed at the start of Section 3.) Here we have used the following calculation of the nonlinear terms of the perturbation equation:


The linear terms will be absorbed into the linearized operator , while the quadratic terms are handled explicitly; the terms in the expansion of the equation around will be denoted by in the sequel. In this work, after expansion in powers of , we drop all terms of order greater than two.

We are interested in linearizing this equation:



Definition 2.3.

A Hamiltonian, is called admissible if the following hold:
1) There are no embedded eigenvalues in the essential spectrum,
2) The only real eigenvalue in is ,
3) The values are non-resonant .

Definition 2.4.

Let (NLS) be taken with nonlinearity . We call admissible if there exists a minimal mass soliton, , for (NLS) and the Hamiltonian, , resulting from linearization about is admissible in terms of Definition 2.3.

The spectral properties we need for the linearized Hamiltonian equation in order to prove stability results are precisely those from Definition 2.3. However note that it is sometimes possible to numerically solve this sort of problem even if Definition 2.3 does not hold; see, for example, [12]. Notationally, we refer to as the projection onto the finite dimensional discrete spectral subspace of relative to . Similarly, represents projection onto the continuous spectral subspace for .

In this work, we must simply assume that is an admissible nonlinearity. However, this assumption is justified by the observed dynamics. Great care must be taken in studying the spectral properties of a linearized operator; although admissibility is expected to hold generically, certain algebraic conditions on the soliton structure itself must factor into the analysis, often requiring careful numerical computations. See [18] as an introduction to such methods and the difficulties therein. To this end, in the forthcoming work [28], two of the authors will look at analytic and computational methods for verifying these spectral conditions.

2.2. The Discrete Spectral Subspace

We approximate perturbations of the minimal mass soliton by projecting onto the discrete spectral subspace of the linearized operator. We now describe, in detail, the discrete spectral subspace at the minimal mass.

Let be the value of the soliton parameter at which the minimal mass soliton occurs. It is proved in [16] (Lemma 3.8) that the discrete specral subspace of at has real dimension . The functions and are in the generalized kernel of at every . Clearly, is purely imaginary and is real. In addition to and , at there are two more linearly independent elements of , the purely imaginary and the purely real .

Applying [16] (Lemma 3.9), and can be extended as continuous functions of in such a way that is purely imaginary, is purely real. We write


with and real-valued functions. The linearized operator, restricted to this subspace, is


where is a differentiable function that is equal to at .

Before proceeding to the derivation of the ODEs, it is helpful to make a minor change of basis. Our goal is that, in the new basis, , , which will make it easier for us to compute the dual basis. Replace by

To preserve the relationship , we need to replace by

To preserve the relationship , we replace by

To preserve , we get that must be replaced by

With these substitutions, the matrix on remains the same and we obtain the relationship . From here on we will assume that we are working with this modified basis and simply take for .

We will define to be the dual basis to the revised within . That is, the are defined by and

If we make the change of basis described above, then we can compute the as follows. Define . Then:

As with the ’s,

for and

for to distinguish between vectors and their scalar components.

3. Derivation of the ODEs

To derive the ODEs we start with a small spherically symmetric perturbation of the minimal mass soliton, then project onto the discrete spectral subspace. Here, we closely follow [16].

We begin with the following ansatz, which allows and to modulate:


Recall, we have assumed to be spherically symmetric, so no other modulation parameters occur. Unlike in [16] we do not assume that the rotation variable is identically zero, so we need to include modulation in our full ansatz. Note that the derivation which follows applies for all nonlinearities in any dimension; specialization is required only to get explicit numerical results. This model includes the full dynamical system for spherically symmetric data and is designed in such a way that coupling to the continuous spectral subspace could easily be reintroduced. The authors plan to analyze the effect of that coupling, which is expected to be dissipative, in a future work.

Differentiating (3.1) with respect to , we get

where we represent differentiation with respect to by and differentiation with respect to the soliton parameter by . Plugging the above ansatz into the equation and cancelling the phase term yields


Recall that since is a soliton solution, , yielding


We multiply by , solve for , and simplify. At the same time, we collect the and terms with the linear portion of , which yields as defined in (2.5). The remaining terms of the nonlinearity are at least quadratic in ; recall that the quadratic terms are described in (2.3) and denoted .

Defining as the coefficient of in , we have

Then, the above calculations give us


Taking the inner product of (3.4) with each of the as defined in Section 2.2, and applying (2.8) yields the following system:


From this point forward in our approximation we drop the component as a higher order error term. Using the product rule, we solve the left hand side for and put the extra terms from the derivative of the operator that projects onto the discrete spectral subspace onto the right hand side.

We have, as in [16], that

where we have implicity defined

and used that

This gives


There is also coupling to the continuous spectrum through terms such as which we omit. This can be included in the error term and is not analyzed in our finite dimensional system.

To make the system well-determined, we must introduce two orthogonality conditions. The first is , and the second is . These represent the choice of and respectively that minimize the size of . These yields , and , respectively.

The reduced system is then:


In [16], the authors further reduce this system to prove there is an initial nonlinear instability. (Note that they have a slightly different system because they have assumed that .) We are interested in the dynamics on an intermediate time scale; thus, we retain quadratically nonlinear terms in our equations.

Our notation is as follows. First, we have

since . Denote

which is the highest order term and the only one that will figure into our quadratic expansion. Notice that this is a real inner product of functions that normally do not appear in the same component of the complex vectors, because of the in the equation.

Similarly, we have

since . Denote

which is again the highest order term.

Then we have

since . Denote

Finally, we have

since . Denote by

We also write for the term at .

Next, consider the terms . These terms are the components of . We have:

Therefore, the relevant nonzero terms are

We denote


Note that some cancellation will occur with the terms that appear separately in the system of ODEs, leaving only these terms in the finally system.

Finally, the terms must be computed. We are only interested in the quadratic terms, which, according to (2.3) are:


Recall that, since and are , the projection onto the discrete-spectrum of is just and the projection onto the discrete-spectrum of is just . We now have to compute the lowest-order terms of

The multiplier of in is

Similarly, we define