Perturbation theory for nonlinear Schrödinger equations

Treating the nonlinear term of the Gross–Pitaevskii nonlinear Schrödinger equation as a perturbation of an isolated discrete eigenvalue of the linear problem one obtains a Rayleigh–Schrödinger power series. This power series is proved to be convergent when the parameter representing the intensity of the nonlinear term is less in absolute value than a threshold value, and it gives a stationary solution to the nonlinear Schrödinger equation.


Introduction
Nonlinear Schrödinger equation (hereafter NLS) is a research topic with a large variety of applications [24]: from problems in nonlinear optics to the analysis of quantum dynamics of Bose-Einstein condensates.In particular, the study of its stationary solutions has attracted increasing attention, and, apart from the few cases in which the solution exists in explicit form, the analysis has mainly focused variational methods or on approximation methods based on both semiclassical and perturbative techniques.
Variational methods are widely used in order to construct bound states for NLS with a linear potential, typically by solving a minimization problem; for instance, this is done by [18] where they proved the existence of a small amplitude stationary solution that bifurcates from the zero solution (see also [21,22] where nonlinear scattering is considered for NLS with, respectively, one and two nonlinear bound states and the references there in).Similarly, in [5] variational methods have been applied to prove that for the minimizer of the nonlinear Hartree energy functional a symmetry breaking effect occurs.
For what concerns semiclassical methods, they have been successfully used in this framework where several authors have been able to demonstrate the existence, in the semiclassical limit using variational techniques, of stationary solution concentrated around the critical points of the potential [1,13,16,27].Also the occurrence of bifurcation phenomena has been discussed in the semiclassical limit [19].
On the other side, the perturbative approach takes up the underlying idea of Rayleigh-Schrödinger series expansion, where the solution is written as a formal series of powers whose coefficients are determined recursively and where the convergence of the series is under investigation.Typically in these cases the perturbation is represented by the nonlinear term, and the unperturbed Schrödinger equation, where the nonlinear term is absent, admits isolated eigenvalues.Several applications of this idea have been developed over the years [2,3,8,12,25,26] limited, in general, to a formal analysis of the series without proving its convergence.In fact, we should emphasize that the problem of convergence of the power series has been solved for some kind of nonlinear Schrödinger equations; more precisely, the spinless real Hartree-Fock model and the Thomas-Fermi-Von-Weizsäcker model has been considered by [8] proving, in particular, that in the first model the Rayleigh-Schrödinger perturbation series has a positive convergence radius.
Finally, it should be mentioned that numerical methods based on discrete Galerkin approximations or spectral splitting methods are widely and effectively used for the study of time-dependent NLS (see [4,6,7,20,23] and references therein).
In this paper we aim to give a rigorous basis to the perturbative approach for computing the stationary solution of the NLS by going so far as to demonstrate, under fairly general assumptions, the convergence of the Rayleigh-Schrödinger series when the perturbative parameter, which measures the intensity of the nonlinear perturbation, is less in absolute value than a given threshold.In this way it is shown that the steady states associated with isolated and nondegenerate eigenvalues of the linear operator transform into stationary solutions of the NLS when nonlinearity is switched on, and the latter can be computed very efficiently through the convergent perturbative series.Finally, it is also possible to give a lower estimate of the radius of convergence of the power series.
The paper is organized as follows.In Section 2 we describe the model, we write the formal power series of the stationary solutions and we state the convergence result in Theorem 1.In Section 3 we state and prove some technical preliminary results.In Section 4 we obtain the convergence of the perturbative series proving thus Theorem 1.In Sections 5 and 6 we discuss a couple of one-dimensional examples: namely in Section 5 we consider the case of an infinite well potential, in this case we are also able to compare the perturbative results with the exact ones; in Section 6 we compute the perturbative series in the case where the potential is the harmonic one.The discussion of these two models is, in some sense, "pedagogical"; indeed, by means of numerical experiments it is possible to see that the coefficients of the power series expansion rapidly decreases and then one can guess the convergence radius of the power series.Finally, in Section 7 we draw some closing comments.A small technical appendix closes the paper.

Main results
2.1.Assumptions.We consider the time-independent nonlinear Schrödinger equation where H = −∆ + V is a linear operator formally defined on L 2 (R d ).The nonlinear term plays the role of perturbation and its strength ν ∈ C is a small perturbative parameter.
Hypothesis 1.The potential V is assumed to be a real-valued piecewise continuous function bounded from below: for some Γ ∈ R.
Remark 1.We assume that the potential V (x) is a piecewise continuous function bounded from below for the sake of simplicity.In fact, we must remark that one could extend our treatment to the case where some milder conditions on V (x) are assumed; however, we don't dwell on those details here.On the other side, it might be interesting to consider the case in which V (x) is given by means of an attractive Dirac's δ [11]; this case does not fall under the Hypothesis 1.
Hence, H admits a self-adjoint extension, still denoted by H, on a self-adjointness domain D(H) ⊂ L 2 (R d ).
Remark 2. Since ϕ 0 ∈ D(H) and the potential V is bounded from below then it follows that ϕ 0 ∈ H 1 because follows from this fact and from the Gagliardo-Nirenberg inequality [9] for some positive constant C p,d and where .
2.2.Formal solutions.We look for a formal stationary solution to (1) close to the solution to the linear problem (3) by means of a formal power series where ν n e n and ψ N (x, ν) and where e n and ϕ n are defined by induction as follows.In fact, E and ψ depend on the perturbative parameter ν; sometimes, for simplicity, we will omit this dependence when this fact does not cause misunderstanding.
Remark 3. We should underline that the following formulas make sense provided that the vectors u n and v n below belongs to L 2 (R d ) and ϕ n ∈ D(H) ∩ L 6 (R d ); we'll discuss this point in Section 3.
Let e ℓ and ϕ ℓ be defined for any ℓ = 0, 1, . . ., n − 1, where ⟨ϕ 0 , ϕ ℓ ⟩ L 2 = 0 for any ℓ = 1, 2, . . ., n − 1, and let We define and By construction it follows that , where Π ⊥ = 1 − Π and Π is the projection operator on the space spanned by ϕ 0 .Hence, the resolvent operator [H − e 0 ] −1 is bounded on Π ⊥ L 2 and we can define Lemma 1.Let e n and ϕ n ∈ Π ⊥ L 2 be defined by induction for any n ≥ 1 as in ( 8) and (9).Let E N and ψ N be defined as in (7).Let Then r N is a power series in ν with finitely many terms where all the coefficients of the powers ν n , with n ≤ N , are exactly zero.
Remark 4. Since e 0 is a simple and isolated eigenvalue of the selfadjoint operator H and since φ n ⊥ ΠL 2 for any n ≥ 1 then: Proof.By formally substituting ( 7) and ( 6) in (1) we then have to check that This equation can be written as where u n and v n are defined above.By equating the term with the same power of the perturbative parameter ν we have that which is satisfied by assumption, and If we multiply both side by ϕ 0 then then φ n ⊥ ϕ 0 and thus we get Here we state our main result.
Remark 5.It is worth noting that the stationary solution ψ given by ( 14) is not normalized to one, that is, to the value of the norm of the unperturbed eigenvector ϕ 0 , which is assumed, for convenience of argument, to be equal to 1.In fact, a simple calculation gives that In particular and where g(ν) is the analytic function obtained by means of the perturbative procedure for ν in a neighborhood of ν = 0 and such that g(0) > 0. If one looks for a normalized solution may act as follows.Let be the normalized stationary solution to the equation where E is still given by ( 14) and where Such a relation is invertible with inverse function In conclusion, if one look for the normalized solution to the equation for a given value of the parameter ν let ν ⋆ be such that ν(ν ⋆ ) = ν, let ψ and E be the perturbative solutions given by ( 14) corresponding to such a value of ν ⋆ ; then ψ/∥ψ∥ L 2 and E are the normalized solution to (16).
In addition, by means of the scaling ψ = ν −2 ω then (1) takes the form of the ν-normalized equation where we have just seen that . Thus, for ν in a neighborhood of 0, we can find a continuous curve (E(ν), ∥ω∥ L 2 ), near the point (e 0 , 0), for solution to (17).Recall that the analysis of the slope of this curve is important in the stability analysis of the stationary state (see, e.g., the "slope condition" in [14]).

L p estimates
As anticipated in Remark 3 it turns out that formulas (8) and ( 9) make sense provided that v n and φ n belongs to L 2 .Hence, we have to prove that ϕ n belongs to L 2 ∩ L 6 for any n.In order to obtain a L p -norm estimate of the vectors ϕ n we make use of the Gagliardo-Nirenberg inequality (4).Lemma 2. Let V (x) be a potential bounded from below (2); let p and C p,d as given in ( 4) and (5).Concerning the H 1 and L p norms of ϕ n we have that for some constants independent of n.

4.
Are the formal series (7) convergent as N goes to infinity?Proof of Theorem 1 In order to prove the convergence of the perturbation series we give the following results.
Lemma 4. Let Proof.In order to prove the result above we remark that and from Lemma 2, where Hence, the above result follows since and □ Lemma 5. Let us assume that for some α > 0 and where Then for any n ≥ 1 and some Remark 8.By construction, (22) holds true for j = 0.
Proof.Indeed, from (22) it turns out that where we set A simple estimate proves that where J(n) has been defined and estimated in Appendix A. Therefore, (n + 1) 2 from which the statement follows.□ Remark 9. From Lemma 5 and from (20) it follows that and for some β ≥ 4b 1 = 4|e 1 | and where and where α > 0 has been introduced in Lemma 5. Then for some C 2 ≤ 2.7.
Proof.The proof immediately follows since where J(n) ≤ 2.7 (n+1) 2 (see Appendix A). □ Remark 11.In fact, the estimate of the constants C 1 and C 2 are far to be optimal.Numerical analysis suggests that a sharp estimate for the term I defined in (24) has the form Concerning C 2 from Appendix A numerical analysis proves that C 2 ≤ 1.52 .Remark 12. From Lemma 6 and from (21) it follows that Remark 13.From Remarks 9 and 12 and from (19) it follows that 2 , where Collecting Lemma 4, Lemma 5 and Lemma 6, we have that if ( 22), ( 25) and ( 26) hold true.
In particular, if we choose and α > 0 large enough such that then we have that ( 22), ( 25) and ( 26) hold true for j = n, too.
In conclusion, we have proved that Lemma 7.There exists four positive constants α > 0 large enough, β > 0, γ > 0 and δ > 0 independent of n such that the following estimates 2 , hold true for any n = 1, 2, . ... Remark 14.From Remark 13 and from Lemma 2 it follows that ψ N is norm convergent in H 1 .

Finally:
Theorem 2. Let d = 1, 2, 3 and let ν be such that |ν| < e −α where α > 0 is large enough as given in Lemma 7. Then the power series E N is absolutely convergent, and the power series ψ N is norm convergent in L 2 and L 6 , and the power series Proof.Convergence of E N and ψ N directly comes from Lemma 7. Concerning the convergence of N n=0 ν n Hϕ n we simply remark that for some C 4 > 0, and thus the formal power series □ So far we have proved that there exists vectors u , w , φ ∈ L 2 , v ∈ L 6 and z ∈ H 1 such that is a bounded operator on the eigenspace orthogonal to ϕ 0 , and where the convergence of the infinite sum has to be intended in the space L 2 .Hence u ∈ D(H) .
Furthermore, we immediately have that In conclusion, there exists a vector ψ ∈ D(H) such that and Thus we have proved the following result.Finally, it's not hard to see that ψ is a stationary solution associated to the energy E to (1).Indeed: Proof.A simple straightforward calculation gives that where the two power series where From the above results immediately follows that as N goes to infinity.From these facts and since (32) then Theorem 1 is proved.

A Toy model -infinite well potential
Let us consider, in dimension one, the infinite well potential of the form: 5.1.Linear time-independent Schrödinger Equation.The linear operator H is formally defined as follows: with Dirichlet boundary conditions By means of a straightforward calculation it follows that the spectrum of H is purely discrete and it is given by means of simple eigenvalues with associated normalized eigenvectors The resolvent operator is given by 5.2.Perturbation theory.By making use of the perturbation formula we compute now the coefficients of the formal power series (7) where is the first unperturbed eigenvalue with associated unperturbed eigenvector Remark 16.Here, we have considered, for argument's sake, the formal power series (7) associated to the first eigenvalue λ 1 .Similarly, the same method may be applied to the unperturbed eigenvalues λ j for any j > 1.
The perturbation theory exploited in Lemma 1 gives that and In Tables 1-2 we compute the values of E N , ∥ψ N ∥ L 2 and ∥r N ∥ L 2 , where r N is the remainder term defined by (10), for different values of N and for ν = ±0.1 and ν = ±1.It turns out that the formal power series (7) rapidly converges and that the norm of the remainder term r N rapidly decreases when N increases.
with Dirichlet boundary conditions (33).If we restrict our attention to the case of ν ∈ R and E ∈ R then we known that the stationary solution is, up to a constant Table 2. Infinite well potential -Table of values corresponding to the case of focusing nonlinearities when ν = −0.1 and ν = −1.
phase factor, a real-valued function.The proof of this result is quite similar to the one of Lemma 3.7 given by [17].Indeed, if we multiply both sides of (35) by ψ, we obtain that and similarly from wich follows that for some constant C. Recalling that ψ(±π) = 0 then C = 0 and thus θ = arg(ψ) is a constant term.Therefore, stationary solutions ψ to (35) may be assumed to be real-valued and they satisfy to the equation The general solution to such an equation has the form [10] where x 0 and ζ are arbitrary constants and where The Dirichlet boundary conditions imply that x 0 = −π and that 2ζπ is a zero of the Jacobian Elliptic function sn(x, k), i.e.: where K(k) and E(k) are the complete elliptic integral of first and second kind.The norm of the wavefunction ψ is given by Hence In order to find the stationary solutions to (36) the quantization conditions read

5.3.1.
Defocusing nonlinearity: ν > 0. When ν > 0 then stationary solutions there exist provided that E − ζ 2 ≥ 0 and k is a real-valued solution to the equation If we remark that the function = +∞ then the equation above (37) has a unique solution k m ∈ (0, 1), for any m = 1, 2, . . .fixed, and then there exists a family of values of the parameter E: 5.3.2.Focusing nonlinearity: ν < 0. On the other hands if ν < 0 then stationary solutions there exist provided that E − ζ 2 ≤ 0 and k = iκ, κ ∈ R, is a purely imaginary complex number; in such a case we recall that sn and where Hence, equation (37) becomes 5.4.Comparison between the perturbative result and the exact one.From Table 1 the perturbative result gives that the stationary solution to (1) for ν = 0.1 and N = 6 has energy with associated wavefunction with norm The value of the solution k to (37), where m = 1 and where the value of ∥ψ∥ L 2 is the one in (42), is given by k = 0.2474031338 and the associated energy E is given by ( 38) in full agreement with (41).Similarly, For ν = 1 and N = 6 then Table 1 gives that  2 gives that with associated wavefunction with norm ∥ψ∥ L 2 ≈ ∥ψ 6 ∥ L 2 = 1.001022099 .

Harmonic oscillator
Let us consider, in dimension one, the harmonic oscillator with potential That is, the linear operator H is defined as follows: It is well known that the spectrum of H is purely discrete and it is given by simple eigenvalues λ j = 2j − 1 , j = 1, 2, . . ., with associated normalized eigenvectors e −x 2 /2 H j−1 (x) , where H j (x) = (−1) j e x 2 d j dx j e −x 2 are the Hermite's polynomials.
By making use of the perturbation formula we compute now the coefficients of the formal power series (7)  In such a case the perturbative procedure gives that Furthermore, where φ 1 = e 1 ϕ 0 − ϕ 3 0 and where the resolvent operator is given by a with infinitely many terms.In numerical calculation we truncate the series for j up to a some large enough positive integer N 2 ; in numerical experiments we observe that N 2 = 60 is a suitable value.Iterating such a procedure we can obtain in Tables 3-4 the numerical values of E N , ∥ψ N ∥ L 2 and ∥r N ∥ L 2 for N = 1, 2, . . ., 6, where r N is the remainder term defined by (10), for ν = ±0.1 and ν = ±1.As in the toy model discussed in Section 5 it turns out that the formal power series seems to rapidly converges for |ν| ≤ 1.Furthermore, a closed expression for J(n) could be given by means of Polygamma functions; however, we don't dwell here on this detail.

Table 1 .
Infinite well potential -Table of values corresponding to the case of defocusing nonlinearities when ν = 0.1 and ν = +1.