ON GROUND STATE (IN-)STABILITY IN MULTI-DIMENSIONAL CUBIC-QUINTIC SCHR ¨ODINGER EQUATIONS

. We consider the nonlinear Schr¨odinger equation with a focusing cubic term and a defocusing quintic nonlinearity in dimensions two and three. The main interest of this article is the problem of orbital (in)stability of ground state solitary waves. We recall the notions of energy minimizing versus action minimizing ground states and prove that, in general, the two must be considered as nonequivalent. We numerically investigate the orbital stability of least action ground states in the radially symmetric case, confirming existing conjectures or leading to new ones.


Basic setting
This work is concerned with the time-evolution corresponding to the cubic-quintic nonlinear Schrödinger equation (NLS) in dimensions  = 2, or  = 3, and subject to initial data The quintic modification of the cubic Schrödinger equation is a model which was introduced in the onedimensional case in [27], as an approximate model in the framework of nonlinear optics.Equation (1.1) appeared more recently in the context of Bose-Einstein condensation, with  = 2 or 3: see e.g.[1,12,25], and [24] for a review.In space dimensions  = 2 or 3, the impact of the quintic term on the dynamical properties of the solution  is stronger than in  = 1, as we shall discuss below.
Depending on the space dimension, which we always assume at most three to simplify the discussion, the nonlinearity in this model is seen to be: focusing  2 -subcritical plus defocusing  2 -critical ( = 1), focusing  2critical plus defocusing  1 -subcritical ( = 2), or focusing  2 -supercritical plus defocusing  1 -critical ( = 3).Recall that for the purely focusing cubic NLS, solitons exist in every dimension and finite time blow-up is possible provided  2 (see e.g.[8]).The presence of the quintic nonlinearity prevents finite time blow-up in  = 2 or 3 (see Prop. 1.1 below), and also affects the stability of solitary waves: understanding this latter aspect more precisely is the main motivation for this paper.For a more precise discussion on the role of criticality in combined power nonlinearities see [29].
The NLS (1.1) formally enjoys the following basic conservation laws: (1) Mass: Momentum: ū(, )∇(, ), (3) Energy: As evoked above, one important effect of the defocusing, quintic term is to prevent finite time blow-up which may occur in the purely cubic case.Indeed, the conservation of the energy, combined with Hölder's inequality, shows that the focusing, cubic part cannot be an obstruction to the existence of a global in-time solution.More precisely, we have, in view of [8] for  = 2 and [34] for  = 3: Proposition 1.1 (Global well-posedness).Let  = 2, 3.For any initial data  0 ∈  1 (R  ), the equation (1.1) has a unique solution  ∈ (R;  1 (R  )), such that  |=0 =  0 .This solution obeys the conservation of mass, energy, and momentum.
We note that in [21], numerical simulations are presented, in which the influence of a small defocusing quintic term on the time-evolution of a focusing cubic NLS is studied.In  = 2 and 3, and for initial data consisting of Gaussians, one obtains a time-periodic (multi-focusing) solution, similar to the one depicted in Figure 11.

Orbital stability of action minimizing ground states
A particular class of global solutions are time-periodic solitary waves of the form (, ) =   (), with  ∈ R and  satisfying For  3, solitary waves  exist provided that the frequency  satisfies the (necessary and sufficient) condition: 0 <  < 3  16 , see [7].Given an admissible  ∈ (0, 3  16 ), we may then look for least action ground states, i.e. solutions   () which minimize the action   () = () +  () among all nontrivial stationary solutions  ∈  1 (R  ).Indeed, it is known from [6,10] that every minimizer of the action   is of the form for some constant  ∈ R,  0 ∈ R  , and with   the unique positive, radial solution to (1.3).In the following, we are mainly interested in the orbital stability of these specific solutions.To this end, we note that, as in the case of more standard, homogeneous nonlinearities (e.g. the cubic case), the NLS (1.1) enjoys three important invariances: (i) Spatial translation: if (, ) solves (1.1), then so does (,  −  0 ), for any given  0 ∈ R  .
The first two invariants are seen to be present in formula (1.4).In combination with the third one, these invariants motivate the following standard notion of stability (see e.g.[8]): Let  be a solution of (1.3).The standing wave   () is orbitally stable in  1 (R  ), if for all  > 0, there exists  > 0 such that if then the solution to (1.1) with Otherwise, the standing wave is said to be unstable.
Following the breakthrough due to M. Weinstein, Grillakis, Shatah and Strauss introduced a general stability/instability criterion in [14] (see also [11]).Assuming certain spectral properties of the linearization of (1.3) about   (which are satisfied in the present cubic-quintic case, see e.g.[7,17,22]), one has the following dichotomy: (i) If    (  ) > 0, then     () is orbitally stable, (ii) If    (  ) < 0, then     () is unstable.This criterion has proven extremely useful in the case of homogeneous nonlinearities, as well as in the case of mixed nonlinearities in  = 1 thanks to an explicit formula, cf.[15,26].In particular, when  = 1, all ground states solitary waves for (1.1) are orbitally stable.However, in the case  = 2 or 3, only partial results are currently available by using the criterion above, see Sections 2.2.2 and 2.2.3, respectively.In our numerical simulations, we will only consider radial perturbation of ground states, and thus remain in the radial framework.The notion of orbital stability then coincides with asymptotic stability up to a phase.
Remark 1.3.The dependence of  upon  in Definition 1.2 is unknown, in general.The proof of stability via the above criterion provides a rather explicit dependence of (a possible)  as a function of , when    (  ) > 0 is known.On the other hand, the stability of the set of constrained energy minimizers (cf.Definition 2.1 below) is obtained by a non-constructive argument.In numerical simulations, tuning initial perturbations of a solitary wave which are sufficiently large to be visible, but not too large (to still adhere to the notion of stability), requires a subtle balance.

Constrained energy minimizers
Since solitary waves may be obtained by other means than minimizing the action   , one may want to look for alternative approaches to orbital stability.An important such alternative is obtained if for  > 0, we denote and assume that the constrained minimization problem has a solution.We call such minimizers energy ground states, in order to make the distinction with least action ground states as clear as possible.Denote by ℰ() the set of such solutions, i.e.

ℰ(𝜌)
Now, let  ∈ ℰ(): Then there exists a Lagrange multiplier Λ such that and thus,  solves the stationary Schrödinger equation (1.3) for some (unknown)  ∈ (0, 3  16 ).Observe that if When the nonlinearity is homogeneous (and  2 -subcritical), this inclusion becomes an equality, see [8,9].However, for non-homogeneous nonlinearities like in (1.1), relating these two constructions of solitary waves (i.e., action minimizing ground states versus constrained energy minimizers) is not obvious at all, and the issue is possibly more complex than it may appear at a glance.First, a-priori nothing guarantees that an element of ℰ() minimizes the action.Second, and this is more subtle: consider a least action ground state   , and let  =  (  ).It is not obvious, and not necessarily true, that   ∈ ℰ().In particular, the map  ↦ →  may not be one-to-one.We will prove that, unlike in the case of homogeneous nonlinearities, we may indeed have   ̸ ∈ ℰ(), cf.Theorem 2.5 below.This fact should be compared to the recent results of [16]: In there, the authors establish for a large class of nonlinearities (including the cubic-quintic one in space dimension  3), that all energy minimizing ground states are least action ground states.In addition, they show that if  is obtained as the Lagrange multiplier associated to the mass constrained  () = , then any least action solution of (1.3) at this value of  is a constrained energy minimizer with the same mass .This paper is now organized as follows: In Section 2, we shall review several known mathematical results on the (in-)stability of ground state solitary waves in  = 1, 2, 3.In particular, we recall the fact that the set of energy minimizers is orbitally stable.We then prove that the dynamics of this set can be distinguished from dispersive behavior and that in  = 3 the sets of action and energy minimizers are not equivalent.In Section 3 we numerically construct action ground states and also collect several of their qualitative properties.Numerical evidence for the orbital stability of these action ground states in  = 2 is then given in Section 4, where we will also describe the numerical algorithm used to simulate the time evolution of (1.1).Finally, we shall turn to the question of orbital stability and instability of 3D action ground states in Section 5, where we will provide numerical evidence for several conjectures on the particular nature of the instability.

Stability for energy ground states
As a preliminary step, we shall recall the Pohozaev identities for the cubic-quintic case (for a derivation, see e.g.[7]): if  ∈  1 (R  ) solves the stationary Schrsdinger equation (1.3), then as well as 2) The aforementioned admissible range for  ∈ (0, 3  16 ) is one of the consequences of these identities.In addition, if  = 2, and after multiplying (2.1) by 2 and subtracting (2.2), we find Therefore, any solitary wave in 2D has negative energy.In the 3D case, this is not necessarily so, as we will see below.
Next, we recall the notion of orbital stability for the set of energy minimizers, as introduced in [9]: Definition 2.1.We say that solitary waves are ℰ()-orbitally stable, if for all  > 0, there exists  > 0 such that if then the solution to (1.1) with This notion is weaker than the one given in Definition 1.2, in the sense that ℰ() may be a large set.In particular, we cannot infer orbital stability of individual members  ∈ ℰ(), as required in Definition 1.2.In the case of homogeneous nonlinearities, however, one can prove that the set ℰ() consists of only a single element  (up to translation and phase conjugation), and thus one recovers Definition 1.2 from the one above.To prove ℰ()orbital stability via the concentration-compactness method, the main step consists in showing that the minimal energy  min () < 0. In our case this yields: Theorem 2.2 (From [7], Section 5).Let  = 2 or 3.
In particular, it seems reasonable to expect that no dispersion is possible near elements of ℰ(), in the sense that a solution (, •) within Definition 2.1 cannot satisfy We will use this criterion as a guiding principle for interpreting several of our numerical findings below.However, it is not clear a priori that for  > 0, such that ℰ() ̸ = ∅, we have The property  min () < 0 makes it possible to rule out this scenario.
Proposition 2.3 (Non-dispersion of energy ground states).Let  = 2 or 3, and  > 0. If  min () < 0, then In particular, there exists  0 () such that for 0 <   0 (), any solution  provided by Definition 2.1 satisfies The statement of this proposition does not involve the above parameter  0 ().For practical application, we want to emphasize that if for a given mass , a stationary solution has negative energy, then solutions around ℰ() cannot disperse.

Further mathematical results
In the following we review some of the known results on orbital (in-)stability of ground states in  = 1, 2, 3.Moreover, we shall prove that 3D action ground states are not necessarily energy minimizers.

The purely cubic case
For the cubic Schrödinger equation with  3, the Cauchy problem is globally well-posed for  = 1 (both in  2 (R) and  1 (R), since the nonlinearity is  2 -subcritical), while finite time blow-up is possible if  = 2 or 3, see e.g.[8].Regarding the solitary waves, the analogue of (1.3) is The corresponding Pohozaev identities imply that such a non-trivial solution exists only if  > 0. Conversely, for any given  > 0, (2.5) has a unique positive radial solution.As a matter of fact, since the nonlinearity is homogeneous, the role of  is explicit: consider  cubic the unique positive radial solution in the case  = 1, Then for any  > 0, is a positive radial solution to (2.5).By uniqueness of such solutions,   is an action ground state, minimizing In particular, we readily compute Recalling the Grillakis-Shatah-Strauss stability criterion, this directly implies that cubic ground states are orbitally stable in  = 1, and unstable in  = 3 (in fact, we have strong instability by blow-up, see e.g.[8]).The 2D case is  2 -critical and instability stems from the fact that the cubic ground state has exactly zero energy (as seen from [32]).Arbitrarily small perturbations can therefore make the energy negative, which consequently leads to finite time blow-up of the associated solution  by a standard virial argument.

Cubic-quintic case in 2D
For  = 2, it follows from Theorem 1.9 of [7] that any solitary wave has a mass larger than that of the cubic ground state, i.e. for any  ∈ (0, 3  16 ), and any solution to (1.3) where  cubic is the radial, positive solution to (2.6).
Recall that   denotes the action ground state.Having in mind the Grillakis-Shatah-Strauss theory, the following asymptotic results have been proved in [7,23]: for  ≈ 0 or  ≈ 3  16 , the map  ↦ →  (  ) is increasing.This implies orbital stability in the sense of Definition 1.2, at least for some range of the frequency  close to the critical values.The numerical plots of  (  ) given in [23] (see also Sect. 3 below) suggest that  ↦ →  (  ) is indeed increasing on the whole range  ∈ (0, 3  16 ), and hence: Conjecture 2.4.In  = 2, all cubic-quintic action ground state solutions are orbitally stable.
In Section 4 we shall give further numerical evidence for this conjecture to be true, by performing several simulations of the time-evolution of perturbed action ground states in 2D.

Cubic-quintic case in 3D
As established in [17,23], when  = 3, it holds: (i) On the one hand, as  → 0, it holds: where  cubic is the positive, radial solution to (2.6).(ii) On the other hand: According to the Grillakis-Shatah-Strauss theory, this implies that cubic-quintic action ground states in 3D are unstable near  min = 0, and orbitally stable near  max = 3  16 .Numerical plots in [23] show a U-shaped curve for  ↦ →  (  ).This suggests the existence of a unique unstable branch and a unique stable branch.We shall numerically investigate the nature of instability in this case in Section 5.
Recalling the fact that the set of (constrained) energy minimizers with negative energy is indeed orbitally stable, cf.Theorem 2.2, we shall now show that solutions to (1.3) may have positive energy when  = 3.Indeed, recalling the main ideas from [17, Section 2.4] (see [17], p. 493 for details), we can rescale   via The new unknown   then solves ∆  +   −  3  +  5  = 0, and, as established in [17], This implies where the last simplification stems from Pohozaev identities for  cubic (discard the  6 norms from (2.1) and (2.2)).Therefore, there exists  0 > 0 such that Recalling that  (  ) → ∞ as  → 0, this shows that there exist unstable action ground states with positive energy and arbitrarily large mass.On the other hand, Theorem 2.2 shows that there exists  0 > 0 such that for all  >  0 , the minimization problem (1.5) has a solution, and  min () < 0. In summary this yields: Theorem 2.5.Not all action ground states in  = 3 are energy ground states.
To our knowledge, this is the first rigorous statement which shows that the two notions of action versus energy ground states, in general, need to be considered as independent (this is also a consequence of [33], Appendix E, as pointed out in [23] -even though it is not stated this way -, since at least for some values of  close to zero,  (  )/ < 0).Note that our theorem is consistent with the results of [16] which establishes that the converse is true, i.e. all energy ground states are action ground states.From the point of view of dynamics, the results above leave open the possibility of having orbitally stable least action ground states, which are not members of ℰ().

Previous results in the 1D case
In the case  = 1, it is fairly natural to generalize the nonlinearity in (1.1) to recover features similar to those of (1.1) when  = 2 or 3.More precisely, consider with 1 <  < .In dimension one, all algebraic nonlinearities are energy-subcritical and explicit solution formulas for   are available in some cases (see in particular [15]).In addition, the focusing part is  2 -critical for  = 4.We therefore expect (2.7) to behave similarly to (1.1) in  = 2, if we choose  = 4, and similarly to (1.1) in  = 3, if  > 4. Indeed, using [14,15], it is proved in [26] that for  = 4, all standing waves (not necessarily ground states) are orbitally stable, while for  > 4, some are orbitally stable (for  ≈  max , computed analogously to the value 3  16 for (1.1)), and some are unstable (for  ≈ 0).Numerical simulations have addressed the case  > 4, see [5,13,28].In particular, [13] reports simulations for perturbations of unstable action ground states, showing two possible dynamics: full dispersion, or convergence to another (stable) soliton.
Remark 2.6.In the case  = 1 and  > 4, the conclusion of Theorem 2.5 remains true, using the same proof.

Numerical algorithm
In this section, we shall discuss a numerical approach for constructing least action ground state solutions to (1.3) in dimensions  = 2 and 3. To this end, we first note that since   is real and radially symmetric, it solves 1 2 where  = ||.In order to get an equation with regular coefficients (which consequently allows for a more efficient numerical approximation), we introduce the new independent variable in which (3.1) reads Since it is known that cubic-quintic ground states are exponentially decreasing, we choose an  0 ≫ 1 such that   ( 0 ) vanishes within numerical precision (which is of the order of 10 −16 here since we work in double precision).Below  0 = 10 3 , while in the next section we shall also consider examples with  0 = 10 4 .The numerical task is thus to find a non-trivial solution to (3.3) for given  ∈ (0, 3  16 ), such that   (numerically) satisfies the homogenous Dirichlet condition   ( 0 ) = 0.
Similarly, the (radial) derivative of  is approximated via the derivative of the Lagrange polynomial, i.e.
At the collocation points, this implies   (ℓ) ≈ Q, since the interpolation polynomial is obviously linear in the ℓ  ,  = 0, . . .,  .Here, ℓ is the vector with components ℓ  ,  is the Chebyshev differentiation matrix [30,31] and Q is the vector with components ((ℓ  )),  = 0, . . .,  .With the above discretization, equation (3.3) is approximated by a system of nonlinear equations for the vector Q which can be formally written in the form F(Q) = 0.The homogenous Dirichlet condition for  =  0 is thereby implemented by eliminating the column and the line corresponding to  0 , cf. [30] for more details.This shows that F() = 0 is an  -dimensional system of nonlinear equations for the  components ((ℓ  )),  = 1, . . .,  .This system will be solved via a Newton iteration, where Jac F denotes the Jacobian of F, and where Q () ,  = 0, 1, . . .denotes the -th iterate.Note, however, that Q = 0 is always a trivial solution to (3.3), which needs to be avoided during the iteration process.In order for this algorithm to converge to our desired, non-trivial solution , one needs to identify a suitable initial iterate Q (0) .To do so, we shall apply a tracing or continuation technique as follows: We introduce in (3.3) a deformation parameter  ∈ [0, 1], such that for  = 0 we have only the focusing cubic nonlinearity, while for  = 1 we obtain the full cubic-quintic equation, i.e. we effectively solve  instead of only (3.3).The cubic solitons  ,=0 are numerically known, see, e.g.[2,20] (no explicit ground state formula exists in dimensions  2).We can thus solve the discretized equation (3.5) for  = 0.1 and, say,  = 0.1 via the Newton iteration described above.The resulting solution  =0.1,=0.1 is then used as an initial iterate for the same equation for  = 0.2, and so on, until the cubic-quintic case  = 1 is reached.In a second step, we use the ground states obtained for  = 0.1 as an initial iterate for slightly smaller or larger 's within the admissible range 0 <  < 3 16 = 0.1875.In this way all examples presented below can be treated.As an example we show in Figure 1 the ground state solutions   () at  = 0.1 for the cubic NLS in blue and the cubic-quintic NLS in red.
It can be seen that the situation is qualitatively different depending on the spatial dimension.Whereas in 2D, the cubic-quintic ground state has a slightly greater maximum and is slightly faster decaying than  cubic , in 3D the cubic-quintic ground state has a much smaller maximum and a considerably larger support.

Numerical ground states in 2D
We first consider the case  = 2 with  = 400 collocation points: In Figure 2    In Figure 3 we depict the ground-state mass  (  ) and energy (  ) as a function of .These plots are based on a total library of roughly 100 numerical ground state solutions  on the shown range of .The corresponding mass-and energy-integrals are thereby computed with the Clenshaw-Curtis algorithm in  =  2 , a spectral integration method based on the same Chebyshev collocation points as before, see [30].Both  (  ) and (  ) appear to be monotonic in .In particular, the monotonicity of  () indicates orbital stability in the sense of Definition 1.2, in view of the Grillakis-Shatah-Strauss theory.

Numerical ground states in 3D
In the case  = 3, we use the same numerical parameters as before: In Figure 4, we show on the left the ground states for several values of .It can again be seen that the maximum of   increases with , at least up to some value  * ≈ 0.1.For larger values of , however, the  ∞ (R 3 )-norm of   is seen to be decreasing again.
Analogously to the 2D case, the solutions become more localized with increasing .Note, however, that despite its exponential decay, the 3D soliton is less localized than in the case of the purely focusing, cubic NLS, see Figure 1 on the right.The 3D ground states in Figure 4 on the left are also found to be less peaked than the corresponding solutions in dimension 2, see Figure 2.
In contrast to the 2D case, the ground state mass  (  ) is no longer monotonically increasing as a function of .Looking at Figure 5, we see that, instead,  (  ) has a minimum at  crit ≈ 0.026.We consequently expect orbital instability of ground states   for  <  crit , a phenomenon we shall study in more detail in Section 5.In Figure 6, the corresponding ground state energy (  ) is seen to have a maximum at the same  crit ≈ 0.026.
The appearance of an unstable branch is clearly visible when the energy (  ) is plotted as a function of the mass  (  ), see Figure 7.
Our figure is in good agreement with Figure 2 of [17] (where the constants are different because the factors are different from (1.1)).Clearly, ground states   corresponding to the upper branch cannot correspond to constrained energy minimizers with mass  =  (  ).

Orbital stability of action ground states in 2D
4.1.Numerical method for the time-evolution In this section we will numerically study the time-evolution of (1.1) resulting from initial data  0 given by perturbations of action ground states.We will only consider perturbations which conserve the radial symmetry.To consider more general perturbations, a full 3D code would be necessary which is beyond the scope of the present paper.This allows us to use the change of variables (3.2) and effectively solve (1.1) in the (non-singular) form We thereby use the same discretization for  ∈ [0, ∞) in terms of Chebyshev collocation points as detailed in the previous section.
After the spatial discretization in , equation (4.1) is then approximated via a system of ordinary differential equations.These equations are then integrated in time using a time-splitting method in which the linear step is solved numerically via an implicit fourth-order Runge-Kutta method, see [19] for more details.The accuracy of this time-integration algorithm is henceforth controlled via the analytically conserved quantity (), which in our case nevertheless depends on time due to unavoidable numerical errors.As discussed in [18], the numerical conservation of the relative mass tends to overestimate the numerical error by one to two orders of magnitude.We shall always aim at a numerical error below (10 −3 ), i.e. below plotting accuracy.This means that we ensure a relative energy-conservation We shall use a single computational domain Ω = [0,  0 ] for which we impose a homogeneous Dirichlet condition (,  0 ) = 0, for all  0. We mostly choose  0 = 10 3 , but in some unstable situations we shall also take  0 = 10 4 .As a basic test case, we first propagate the three-dimensional ground state   , numerically found at  = 0.1.We thereby use   = 10 4 time-steps until a final time  f = 10.We find that the hereby obtained numerical solution , at  =  f , satisfies max i.e. the same order of accuracy as reached in [19].
Remark 4.1.In general, it is not unproblematic to work with a homogeneous Dirichlet boundary condition on a finite numerical domain, since this could lead to unwanted reflections of the emitted radiation at the boundary, see, e.g., the discussion in [3].In our case, however, only small, rapidly decreasing perturbations of ground states are considered.It is thus possible to work on sufficiently large computational domains Ω, on which the radiation can separate from the bulk before spurious reflections from the boundary lead to noticeable effects.

Time-evolution of perturbed 2D ground states
In this subsection, we shall study the time-evolution of perturbed ground states to (1.1) in dimension  = 2.To this end, we first consider the case where  0 () =   (),  > 0. (4.2) Here  > 0 is a perturbation parameter and   is a numerically obtained action ground state, at a certain admissible frequency  ∈ (0, 3  16 ).We first study the case where  = 0.1 and  = 0.99, and use   = 10 4 time steps to reach the indicated final time  f = 20.As expected, the solution  to (1.1), effectively given by (4.1), is found to be close to the exact time-periodic state   (, ) =     ().
To this end, we show on the left of Figure 8 the  ∞ (R 2 )-norm of the solution as a function of time.It can be seen that it approaches a final state as  → 20.The latter is found to be very close (in absolute value) to the unperturbed ground state  =0.1 .Note that the  ∞ -difference is of the order of 10 −4 and thus, much smaller than the initial perturbation.As a second case, we consider the same 2D initial data (4.2), but with  = 1.001.In Figure 9 we again show the  ∞ -norm of the solution as a function of time.Similarly as before, a final state is reached and its maximum is again found to be very close to the unperturbed ground state  =0.1 .In both cases, we find that the difference is largest for  close to the origin.
In order to illustrate that the qualitative picture found before is not due to our specific choice of perturbations, we shall also consider ground states perturbed by a small Gaussian-like perturbation, i.e.  and  = 0.05.On the left the case with the "+" sign.On the right the case with the "−" sign.
Note that we only consider smooth perturbations in this paper in order to allow for spectral accuracy in the radial coordinate, i.e., an exponential decrease of the numerical error with the number of collocation points.In Figure 10 we show the behavior in time of the respective  ∞ -norms for the two choices ±.In both situations the difference between || at  f = 20 and  =0.1 is found to be of the order (10 −4 ).Moreover, the error (not depicted here for the sake of readability) is again found to be largest close to the origin.
If similar perturbations are applied to other ground states   , the resulting solution  behaves qualitatively similarly.Our numerical tests therefore support Conjecture 2.4.However, we also find that the smaller the choice of  ∈ (0, 3  16 ), the longer it takes for the solution  to reach its final state.In fact, for small enough , damping effects within the time-oscillations of || become almost invisible, even if one computes up to much larger times  f = 400, see Figure 11.

Stable branch
In this section, we shall study the question of (in-)stability of cubic-quintic ground states in dimension  = 3.In view of Figure 5, we expect ground states   with  >  crit ≈ 0.026 to be orbitally stable.That this is indeed the case, is strongly suggested by our numerical results below.
To this end, we first consider multiplicative perturbations of   on the stable branch: In Figures 12 and 13 we study the time-evolution of (1.1) with initial data of the form (4.2).On the left of Figure 12 we show the  In Figure 13 we study the analogous situation with  = 1.001:Again, the (absolute value of the) solution  seems to settle around  f = 15 on the stable unperturbed ground state  0.1 .In both cases, the error between || and  =0.1 is again found to be of the order (10 −4 ).

Unstable branch
The situation dramatically changes if we consider perturbations of ground state solutions on the unstable branch, i.e. perturbations of   with  <  crit ≈ 0.026: In Figure 14 we show the solution  to (1.1) obtained from initial data (4.2) with  = 0.01 and  = 0.999.Note that this implies  ( 0 ) <  (  ).The solution is seen to be purely dispersive which is also confirmed by the  ∞ -norm of the solution as a function of time (depicted on the right of the same figure).In fact, we did not discover any stable structure within the time-evolution even if we let the numerical code run for longer times.
If we consider the same ground state as before, but instead choose  = 1.001, we find a different kind of instability.Now the solution  shows oscillations of high amplitude, see Figure 15.These oscillations are even more visible in the  ∞ -norm of the solution, as depicted on the left of Figure 16.One can see that early on the norm is growing strongly but then it appears to show damped oscillations around some final state.We conjecture that the latter corresponds to another ground state on the stable branch.To this end, we compare the maximum of ||, obtained at the final time  f = 500, with the  ∞ -norms in our library of previously computed action ground states   , cf. Figure 4. Indeed we find good agreement of ||, when compared to   with  = 0.047 >  crit , see the right of Figure 16.Thus perturbations of unstable ground states where  ( 0 ) >  (  ), seem to result in solutions which eventually settle on another, stable ground state as  → +∞.Note, however, that  ( =0.01 ) ≈ 79.44 while  ( =0.047 ) ≈ 77.05.If the final state had the same mass as the unperturbed initial state, this would correspond to an  ≈ 0.0495.This shows that a certain amount of mass is lost through radiation.
(5.1)For both | 0 | = 0 and | 0 | ̸ = 0, the solution in the case with the "+" sign looks very similar to the one depicted in Figure 16.This fact becomes particularly clear when one compares the time-evolution of the  ∞ -norm of  depicted in Figure 17, with the one from Figure 16.
By comparing the maximum of || found at the final time  f = 600 with the  ∞ -norm of a stable ground state, we find good agreement with  =0.048 .The latter has mass  ( =0.048 ) ≈ 77.95.Unfortunately, we are unable to clearly decide whether the final state is closer to  =0.048than to  =0.047 .
In the case of initial data (5.1) with the "−" sign, we again find that the solution is completely dispersed, see Figure 18.This is consistent with our earlier findings above which indicate that perturbation with  ( 0 ) <  (  ) lead to purely dispersive solutions.
We finally note that the situation is qualitatively similar for other values of  on the unstable branch.In our last example, we choose  = 0.007 <  crit within the initial  0 given by (5.1).We again find that perturbations with an initial mass smaller than  ( =0.007 ) are purely dispersive.Perturbations with mass larger than  ( =0.007 ) lead to damped oscillations around some final state, see Figure 19.
The asymptotic final state appears to be close to  =0.044 on the stable branch.The mass of the unperturbed initial data  ( =0.007 ) ≈ 90.57 is seen to be bigger than  ( =0.044 ) ≈ 75.37, showing again that a nonnegligible part of the initial mass has been radiated away.
Remark 5.2.The same instability scenario was found (numerically) for perturbed solitary wave solutions to the generalized BBM equation in [4], and for a version of NLS with derivative nonlinearity in [2].In particular, analogously to our situation, a perturbation which lowered the mass of the initial data below the one of the (unstable) solitary wave always resulted in purely dispersive solutions.In all of these cases, it remains an interesting open question to find a possible selection criterion for the specific value  which describes the (stable) asymptotic state   .

Figure 1 .
Figure 1.Ground state solutions  =0.1 to the cubic NLS in blue and the cubic-quintic NLS in red: on the left for  = 2 and on the right for  = 3.

Figure 2 .
Figure 2. Left: Ground state solutions to (3.1) in dimension  = 2 for several values of .Right: The  ∞ -norm of these states as a function of .
we show on the left a plot of the ground state function   () for various values of .It is seen that the maximum of the ground states increases with .The solutions also become more localized with increasing .On the right of the same figure, we show the  ∞ (R 2 )-norm of the ground states as function of .For convenience, we only consider values of  ∈ [0.005, 0.16].

Figure 4 .
Figure 4. Left: Ground state solutions to (3.1) in dimension  = 3 for several values of .Right: The  ∞ -norm of these states as a function of .

Figure 5 .
Figure 5. Left:  (  ) as a functions of , for cubic-quintic ground states   in dimension  = 3. Right: a close-up of the same curve near  crit .

Figure 6 .
Figure 6.Left: (  ) as a functions of , for cubic-quintic ground states   in dimension  = 3. Right: a close-up of the same curve near  crit .

Figure 8 .
Figure 8.The 2D solution  to (1.1) for initial data (4.2) with  = 0.1 and  = 0.99.On the left the  ∞ -norm as a function of time.On the right the difference between || and  =0.1 at the final time  f = 20.

Figure 9 .
Figure 9.The 2D solution  to (1.1) for initial data (4.2) with  = 0.1 and  = 1.001.On the left the  ∞ -norm as a function of time.On the right the difference between || and  =0.1 at the final time  f = 20.

Figure 10 .
Figure 10. ∞ -norm as a function of time for the 2D solution  to (1.1) with initial data (4.3)and = 0.1.On the left the case with "+" sign.On the right the case with "−" sign.

Figure 11 .
Figure 11. ∞ -norm as a function of time for the 2D solution  to (1.1) with initial data (4.3)and = 0.05.On the left the case with the "+" sign.On the right the case with the "−" sign.

Figure 12 .
Figure 12.The 3D solution  to (1.1) for initial data (4.2) with  = 0.1 and  = 0.99.On the left the  ∞ -norm as a function of time.On the right the difference between || and  =0.1 at the final time  f = 15.

Figure 13 .
Figure 13.The 3D solution  to (1.1) for initial data (4.2) with  = 0.1 and  = 1.001.On the left the  ∞ -norm as a function of time.On the right the difference between || and  =0.1 at the final time  f = 15.

Figure 14 .
Figure 14.Left: The absolute value of the solution  to (1.1) in dimension  = 3 obtained from initial data (4.2) with  = 0.01 and  = 0.999.Right: The  ∞ -norm of  as a function of time.

Figure 16 .
Figure 16.Solution to (1.1) in dimension  = 3 for initial data (4.2) with  = 0.01 and  = 1.001:On the left the  ∞ -norm of  as a function of time.On the right || at the final time (in blue) together with the ground state  =0.047(in green).

Figure 17 .
Figure 17.Solution to (1.1) in dimension  = 3 for initial data of the form (5.1) with | 0 | = 1, and  = 0.01: On the left the  ∞ -norm of  as a function of time.On the right || at the final time (in blue) together with the ground state  =0.048(in green).

Figure 18 .
Figure 18.The  ∞ -norm of the solution  in dimension  = 3 obtained from initial data (5.1) with the "−" sign.On the left the case with  0 = 0 and on the right the one with  0 = 1.