Nonlinear Bloch Waves in the Gross-Pitaevskii Equation

Bifurcation of Nonlinear Bloch Waves from the Spectrum in the Gross-Pitaevskii Equation

Tomáš Dohnal T. Dohnal Department of Mathematics, Technical University Dortmund D-44221 Dortmund, Germany  and  Hannes Uecker H. Uecker Institute of Mathematics, Carl von Ossietzky University Oldenburg D-26111 Oldenburg, Germany
July 3, 2019

We rigorously analyze the bifurcation of stationary so called nonlinear Bloch waves (NLBs) from the spectrum in the Gross-Pitaevskii (GP) equation with a periodic potential, in arbitrary space dimensions. These are solutions which can be expressed as finite sums of quasi-periodic functions, and which in a formal asymptotic expansion are obtained from solutions of the so called algebraic coupled mode equations. Here we justify this expansion by proving the existence of NLBs and estimating the error of the formal asymptotics. The analysis is illustrated by numerical bifurcation diagrams, mostly in 2D. In addition, we illustrate some relations of NLBs to other classes of solutions of the GP equation, in particular to so called out–of–gap solitons and truncated NLBs, and present some numerical experiments concerning the stability of these solutions.

Key words and phrases:
periodic nonlinear Schrödinger equation, nonlinear Bloch wave, Lyapunov-Schmidt decomposition, asymptotic expansion, bifurcation, delocalization
2000 Mathematics Subject Classification:
Primary: 35Q55, 37K50 ; Secondary: 35J61

1. Introduction

The nonlinear Schrödinger/Gross–Pitaevskii (GP) equation in dimensions,


with a real potential is a canonical model in mathematics and physics. It appears in various contexts, e.g., nonlinear optics [33, 17], and Bose–Einstein condensation [26, 4]. See also, e.g., [34, 40, 19] for mathematical and modeling background. Plugging into (1.1), where is the frequency of time–harmonic waves in nonlinear optics, and where is called the chemical potential in Bose–Einstein condensation, we obtain the stationary problem


Here we consider the case that the potential is real and periodic. For simplicity, we let be periodic in each coordinate direction, i.e.,

where denotes the -th Euclidean unit vector in . In other words, we consider the periodic lattice . We make the basic assumption that for some , where . This smoothness assumption on ensures -smoothness of linear Bloch waves, i.e., solutions of (1.2) with . See §1.1 for a review of spectral properties of

and linear Bloch waves. For suitable the spectrum of shows so called spectral gaps and in recent years a focus has been on the bifurcation of so called gap solitons from the the zero solution at band edges into the gaps. These are localized solutions, which, in the near edge asymptotics have small amplitude and long wave modulated shape. In detail, the asymptotic expansion at with is


where , are Bloch waves at the edge , and the are localized solutions of a system of (spatially homogeneous) nonlinear Schrödinger equations. See, for instance, [32, 12, 15, 23], and the references therein.

Here we seek solutions of (1.2) which can be expressed as finite a sum of quasi-periodic functions and call such solutions nonlinear Bloch waves (NLBs), with quasi–periodicities determined from a selected finite subset of the Bloch waves at . NLBs have been studied in, for instance, [12, 38, 42, 43, 9], where in [38, 42, 43] the approaches are numerical and formal. They have been observed even experimentally, see e.g. [10] for experiments in Bose-Einstein condensates. In [12] the special case of a bifurcation of NLBs into an asymptotically small spectral gap for a separable periodic potential in two dimensions is studied rigorously. In [9] the bifurcation of single component () NLBs in one dimension is proved, including results on secondary bifurcations and exchange of stability. Similarly to Bloch waves in linear lattices NLBs can be understood as the fundamental bounded oscillatory states of the nonlinear system. From the applied point of view one motivation for studying NLBs is the continuation of gap–solitons to “out-of-gap” solitons, i.e., the continuation of localized solutions from one band edge across the gap and into the spectrum on the other side of the gap, where their tails start interacting with the NLBs. For this reason, the study of bifurcation of NLBs from the zero–solution has been mostly restricted to band edges. Here we show that nonlinear Bloch waves bifurcate in from generic points in the spectrum of , and give their asymptotic expansions in terms of solutions of the so called algebraic coupled mode equations (ACME), together with error estimates.

In addition to the rigorous analysis we illustrate our results numerically. For this we focus on 2D, as this is much richer than 1D, and use the same potential as in [15], i.e.



This represents a square geometry with smoothed-out edges. The function in (1.4) is extended periodically to to obtain . The numerical band structure of over the Brillouin zone , and also along the boundary of the irreducible Brillouin zone, is plotted in Fig. 1(a),(b), respectively. We denote the so called high symmetry points in for by

(a) (b)
Figure 1. (a): Band structure of over the Brillouin zone for the periodic potential (1.4); (b): along the boundary of the, so called, irreducible Brillouin zone.
Example 1.

Figure 2 shows a numerical bifurcation diagram of single component () NLBs for , calculated with the package pde2path [36, 13], together with example plots on the bifurcating branches. A branch of NLBs bifurcates from the zero solution at for any attained by one of the band functions at , i.e. at the coordinate of any of the points in Fig. 1 (b). See §1.1 for the definition of band functions.

Figure 2. Example 1: Bifurcation diagram of the first four bifurcating branches for , i.e. branches bifurcating from points a-d in Fig. 1 (b). Spectral bands are indicated by the black dashed line. The sign in the branch labels stands for . Small panels: example solution plots of NLBs from the bifurcation diagram, over the fundamental cell . At bifurcation we choose a real Bloch wave. Then the imaginary parts are small near bifurcation, and we only plot them for a. Roughly horizontal axis corresponds to in all plots.

In §7 we explain the method behind Fig. 2, and study in detail the bifurcations of NLBs at the points marked A,B,C in Fig. 1(b), relating the numerical calculations to our analysis.

As already said, one motivation for studying NLBs are the intriguing properties of their interaction with localized solutions, which we illustrate numerically in §8. For instance, when a gap soliton is continued from the gap into the spectrum, we get a so called “out–of–gap” soliton (OGS) with oscillating (delocalized) tails, see also [41, 24]. In 1D, numerically these OGS can be seen to be homoclinic orbits approaching NLBs, and essentially the same happens in 2D. Moreover, the NLB can form building blocks of so called truncated NLBs (tNLBs), see also [4, 38]. These are localized solutions for in a gap which are close to a NLB on some finite interval but approach as . Then, continuing a branch of tNLBs from the gap into the spectrum, the same interaction scenario as for GS happens, i.e., the tails of the tNLBs pick up NLBs bifurcating from the gap edge into the spectrum, and the tNLBs become delocalized, for which we use the acronym dtNLB. Note that both gap solitons and tNLBs have been observed experimentally, see e.g. [5] for experiments in optical lattices. However, even in 1D at present it is unclear how to analyze OGS, tNLBs and dtNLBs rigorously, i.e., so far there only exist heuristic asymptotics, see §8 for further comments.

Stability of most of these solutions is an open problem. Thus, at the end of the paper we also give a numerical outlook on this, and obtain stability of some NLBs in 1D, and, consequently, stability of some tNLBs and some OGS and dtNLBs. In 2D, we did not find stable NLBs for the potential from (1.4), and we close with summarizing the open questions. The broad spectrum of applications of NLBs clearly motivates our rigorous bifurcation analysis.

In the remainder of this introduction we explain the linear band structure, a simple analytical bifurcation result, formulate the main theorem, and describe the structure of the paper in more detail.

1.1. Linear Bloch waves

For in the Brillouin zone, , consider the Bloch eigenvalue problem


where is the th Euclidean unit vector in . The spectrum of is continuous and is given by the union of the bands defined by the band structure , i.e.

see Theorem 6.5.1 in [16]. The functions are called band functions. The Bloch waves have the form with for all and all . Clearly, both and are periodic in each component of . We assume the normalization

For a given point in the band structure, i.e. with for some , also the point lies in the band structure, which follows from the symmetry


This symmetry is due to the equivalence of complex conjugation and replacing in the eigenvalue problem (1.5). Hence, we also have the conjugation symmetry of the Bloch waves, namely


For we have and the point must be understood as the -periodic image within . When is one of the so called high symmetry points, i.e. for all , then and are identified via this periodicity. Equation (1.7) then implies that is real. This can be seen directly from the eigenvalue problem (1.5), where for all implies that the boundary condition is real such that a real eigenfunction must exist. Note that the above symmetries rely only on the realness of .

1.2. The Bifurcation Problem

Remark 1.

In the simplest scenario we can look for real solutions of (1.2) with the quasi-periodic boundary conditions given by a single vector , i.e.


In this case the realness condition on requires


such that the sought solution is -periodic or -antiperiodic in each coordinate direction. We study bifurcations in the parameter . Classical theory for bifurcations at simple eigenvalues, e.g., Theorem 3.2.2 in [29], shows that if for exactly one , i.e. is a simple eigenvalue of under the boundary conditions (1.8), then is a bifurcation point. To this end define and study on under the boundary conditions (1.8). We have for all and . As is a simple eigenvalue, we have the one dimensional kernel

Because with (1.8) and (1.9) is self adjoint, we have . The transversality condition of Theorem 3.2.2 in [29] thus holds because . As a result, the theorem guarantees the existence of a unique non-trivial branch of solutions bifurcating from .

Remark 2.

Without the restriction to real solutions the eigenvalue is never simple due to invariances. In the real variables , where , the problem becomes

Since (1.2) possesses the phase invariance and the complex conjugation invariance, we get that is invariant, i.e.

Bifurcations can now be studied using the equivariant branching lemma, see e.g. [28, §5], by restricting to a fixed point subspace of a subgroup of . The only nontrivial subgroup is with the fixed point subspace being the vectors with corresponding to real solutions of (1.2). Therefore, this leads again to real solutions. Nevertheless, more complicated solutions than the single component ones in Remark 1 can be studied. The most general real ansatz is


with , with for all , such that for all , and with for . Here means equality modulo 1 in each coordinate. While the use of the equivariant branching lemma should describe the bifurcation problem and produce the effective Lyapunov-Schmidt reduction, we choose to carry out a detailed analysis without this tool in order to obtain more explicit results. This will allow us to provide estimates of the asymptotic approximation error.

Our aim is to prove a general bifurcation theorem for NLBs, and, moreover, to derive and justify an effective asymptotic model related to the Lyapunov-Schmidt reduction of the bifurcation problem including an estimate on the asymptotic error. In our approach we select a frequency in the spectrum and choose points in the level set of the band structure at , such that for each we have for some . Our method requires that each of the points is either one of the high symmetry points or belongs to a pair with . See (H1)–(H6) on page 3.2 for a summary of our assumptions. We seek NLBs bifurcating from and having the asymptotic form


at . The coefficients , i.e. the (complex) amplitudes of the waves, are given by solving the ACME as an effective algebraic system of equations. Generally a sum of quasiperiodic functions with the quasiperiodicity of each given by one of the vectors cannot be an exact solution of (1.2) as the nonlinearity generates functions with other quasiperiodicities. Our ansatz for the exact solution is thus

with and for . Importantly, this set (defined in (3.2)) can be a proper subset of the level set. The subset may be finite even if the level set is, for instance, uncountable. In fact, our assumption (H4) ensures the finiteness. Besides, the subset can be much smaller than and hence the effective ACME-system can be rather small. The subset has to satisfy only (H2-H6).

The major assumptions of our analysis are rationality (assumption (H4)) and certain non-resonance conditions (H5) on the -vectors . In addition, the solutions of the coupled mode equations need to satisfy certain symmetry (“reversibility”) and non-degeneracy conditions, see Definitions 2 and 3, in order for us to guarantee that (1.11) approximates a solution of (1.2).The main result is the following

Theorem 1.

Assume (H1)-(H6). There exist and such that for all the following holds. If the ACMEs (2.3) have a reversible non-degenerate solution , then (1.2) with has a nonlinear Bloch wave solution of the form (3.4), and

There are three relatively straightforward generalizations of the result. Firstly, the periodic lattice can be replaced by any lattice with linearly independent vectors . Of course, the periodicity cell and the Brillouin zone have to be redefined accordingly. Except for the examples in §6 the results, in particular Theorem 1, hold for a general lattice. Secondly, the nonlinearity can be replaced by other locally Lipschitz nonlinearities which are phase invariant and satisfy for . This may, however, change the powers of in the expansion and the error estimate. Also, the linear operator can be generalized to self adjoint second order differential operators with periodic coefficients such that the asymptotic distribution of eigenvalues remains that in (3.6).

1.3. The Structure of the Paper

In §2 we present a formal asymptotic approximation of nonlinear Bloch waves and a derivation of the ACMEs as effective amplitude equations. In §3 we pose conditions on the solution ansatz and the band structure which are necessary for our analysis, and apply the Lyapunov-Schmidt decomposition to the bifurcation problem. The invertible part of the decomposition is estimated in §4. The singular part and its relation to the ACMEs is described in §5, where also the proof of the main theorem is completed. In §6 we present the ACMEs and their solutions in the scalar case () and in the case of two equations (). Section 7 presents numerical computations of nonlinear Bloch waves in two dimensions for and . The convergence rate of the approximation error is confirmed by numerical tests. Finally, in §8 we give a numerical outlook on the interaction of localized solutions with NLBs, first for some 1D and 2D GS, and second for tNLBs, and we report numerical experiments on stability of NLBs and other solutions.

2. Formal Asymptotics

Let and choose vectors in the level set of the band structure at . For the asymptotics of nonlinear Bloch waves near we make an analogous ansatz to that used in [12, §3] for nonlinear Bloch waves near band edges in (1.2) with a separable periodic potential. Formally we write


where the amplitudes are to be determined and where satisfies the quasiperiodicity given by the vector .

Substituting (2.1) in (1.2) we get at for each



The condition in the sum ensures that the nonlinear terms have the same quasi-periodicity as . Nonlinear terms generated by the ansatz (2.1) and having other quasi-periodicity than one of those defined by have been ignored in this formal calculation.

Imposing the solvability condition, i.e. making the right hand side -orthogonal to , we get the algebraic coupled mode equations (ACMEs)


To make the approximation (2.1) rigorous, we must account for the nonlinear terms left out above and provide an estimate on the correction .

3. Solution Ansatz, Assumptions, Lyapunov-Schmidt Decomposition

As mentioned above, one of the difficulties of the analysis is that for a sum of functions with distinct quasi-periodic conditions the nonlinearity can generate functions with a new quasi-periodicity. If the -points defining these new quasi-periodic boundary conditions lie in the -level set of the band structure, then a resonance with the kernel of the linear operator occurs. Also, if the points generated by a repeated iteration of the nonlinearity merely converge to the level set, our techniques fail because a lower bound on the inverse of the linear operator cannot be obtained. These obstacles are avoided if for a selected assumptions (H4) and (H5) below hold.

We select points in the -level set of the band structure. Suppose we seek solutions of (1.2) with given by the sum of quasiperiodic functions. The ansatz with quasiperiodic such that for all can be a solution of (1.2) only if each term generated by the nonlinearity applied to this sum has quasiperiodicity defined by one of the vectors in , i.e. if the consistency condition



is satisfied. In other words the consistency condition (3.1) says that all combinations for must lie in , with from (2.2).

An example of a consistent ansatz for is with , like e.g. for in [15]. On the other hand, for , where with , see Sec. in [15], the ansatz is inconsistent. It is also inconsistent for typical in the interior of with generic in the level set. Therefore, we drop the consistency condition and pursue the more general case where the nonlinearity generates quasiperiodic functions with quasi-periodicity vectors not necessarily contained in .

Hence, we define the set of -points generated by iterations of the nonlinear operator


and write, with ,


At this point is possible but as explained below, our assumption (H4) ensures , i.e. only finitely many new vectors are generated. Thus we can search for a solution in the form of the sum of finitely many quasiperiodic functions


with . The choice of the function space for is made clear below.

We make the following assumptions:

  • for some , where ;

  • and are points in the -level set of the band structure, i.e., there are such that

  • each point is repeated according to the multiplicity of at . In detail, if band functions touch at , then points in equal and ;

  • the points have rational coordinates, i.e.

  • the intersection of the set with the level set of the band structure at is exactly the set , i.e.


  • for each the reflection w.r.t. the origin lies in the set too, i.e.

    where and denotes congruence with respect to the -periodicity in each component.

With (H3) the bifurcation from multiple Bloch eigenvalues is allowed. In one dimension () multiplicity is at most two, which occurs for so called finite band potentials, see e.g. [27, 11], and only at or . In higher dimensions () crossing or touching of band functions is abundant in generic geometries (potentials ). Our results thus apply also to Dirac points in two dimensions studied, e.g., in [18].

Due to the rationality condition (H4) the sought solution (3.4) is, in fact, periodic. (H4) also ensures that the set is finite (). Indeed, iterating the operator on a set of points with rational coordinates on a -dimensional torus generates a periodic orbit, i.e. only finitely many distinct points are generated, and the number depends solely on . Condition (H4) is satisfied, e.g. if is a subset of the high symmetry points of , i.e. for all . This is frequently the case for the locations of extrema defining a spectral edge. In general, (H4) is, however, a serious limitation, and removing this assumption would be a major improvement.

The non-resonance condition (H5) is satisfied, for instance, if , i.e. for at one of the band edges, and are all the extremal points of the band structure at which the edge is attained.

The symmetry condition (H6) is needed in the persistence step of the proof, see §5. Note that if , then also by (1.6) and the periodicity in . For , clearly, . For is (e.g. for with ) and then is the -periodic image of within (for the example ). Moreover, also (H6) is automatically satisfied if is a subset of the high symmetry points of .

Note again that as well as may be proper subsets of . This is, for instance, the case in 2D if , where we could choose (if is simple), and or , which yields two decoupled scalar bifurcation problems; see §6.2 for further discussion.

The remaining two assumptions in Theorem 1 are non-degeneracy and reversibility of , defined as follows.

Definition 2.

is a non-degenerate solution of (2.3) if the Jacobian111Strictly speaking, the problem should first be rewritten in real variables to define a Jacobian, see the discussion above Lemma 10, but for brevity we use this compact symbolic notation here. , where , has a simple zero eigenvalue.

Note that due to the phase invariance of (2.3) the Jacobian is singular.

Definition 3.

is reversible if

where is given by .

Reversibility is a symmetry of the solution. The motivation for restricting to reversible non-degenerate solutions is to ensure the invertibility of in the fixed point iteration for the singular part of the Lyapunov-Schmidt decomposition in §5. Within the phase invariance is, indeed, no longer present. The choice of in Definition 3 is natural and based on the intrinsic symmetry (1.7) of the Bloch eigenfunctions which ensures the complex conjugation symmetry among the coefficients in (2.3) and, hence, the possibility of reversible solutions. Note that (1.7) follows directly from .

Next, we assume (H1-H6) and use the Lyapunov-Schmidt decomposition in Bloch variables together with the Banach fixed point theorem to prove the main result, i.e., Theorem 1, which justifies the formal asymptotics for solutions at .

3.1. Lyapunov-Schmidt Decomposition

Due to the completeness of the Bloch waves in we can expand


As the following lemma shows, working with in the space is equivalent to working with , where

Lemma 4.

For the following norm equivalence holds. There exist constants such that

where is related to by (3.5).

The proof is analogous to that of Lemma 3.3 in [8], see also [15, §4.1]. The main ingredients are firstly the fact that for large enough (such that for all and , e.g. ) the squared norm is equivalent to

Secondly, one uses the asymptotic distribution of bands in dimensions, see [21, p.55]: there are constants such that


For the subsequent analysis we define for each the set of indices producing through the nonlinearity analogously to the definition of , i.e.

For the ansatz (3.4), (3.5) equation (1.2) is equivalent to the algebraic system



Due to the kernel of the linear multiplication operator at in (3.7) we use a Lyapunov-Schmidt decomposition in order to characterize the bifurcation from (i.e. from ). For we let

and write