Wave Breaking in Dispersive Fluid Dynamics of the Bose-Einstein Condensate

TL;DR

This paper analyzes wave breaking and dispersive shock wave formation in Bose-Einstein condensates using the Gross-Pitaevskii equation, providing explicit solutions for specific initial profiles and shock velocities.

Contribution

It derives closed-form solutions for the Whitham modulation equations for initial profiles with power-law forms, advancing understanding of wave breaking in BECs.

Findings

Closed-form solutions for n=2,3 cases.

Edge velocities of shock waves for arbitrary n>1.

Application to experimental BEC wave phenomena.

Abstract

The problem of wave breaking during its propagation in the Bose-Einstein condensate to a stationary medium is considered for the case when the initial profile at the breaking instant can be approximated by a power function of the form $(-x)^{1/n}$. The evolution of the wave is described by the Gross-Pitaevskii equation so that a dispersive shock wave is formed as a result of breaking; this wave can be represented using the Gurevich-Pitaevskii approach as a modulated periodic solution to the Gross-Pitaevskii equation, and the evolution of the modulation parameters is described by the Whitham equations obtained by averaging the conservation laws over fast oscillations in the wave. The solution to the Whitham modulation equations is obtained in closed form for $n = 2,3$, and the velocities of the dispersion shock wave edges for asymptotically long evolution times are determined for arbitrary integer values $n > 1$. The problem considered here can be applied for describing the generation of dispersion shock waves observed in experiments with the Bose-Einstein condensate.

Figures

Equations

$$i{\psi }_t+\frac{1}{2}{\psi }_{xx}-|\psi {|}^2\psi =0,$$

$$\psi (x,t)=\sqrt{\rho (x,t)}\exp \left(i{\int }^{x}u(x^{\prime },t)d{x}^{\prime }\right),$$

$$\begin{array}{rl}& {\rho }_t+(\rho u{)}_x=0,\\ & u_t+u{u}_x+{\rho }_x+{\left[\frac{{\rho }_x^2}{8{\rho }^2}-\frac{{\rho }_{xx}}{4\rho }\right]}_x=0.\end{array}$$

$$\begin{array}{rl}\rho =& \frac{1}{4}({\lambda }_4-{\lambda }_3-{\lambda }_2+{\lambda }_1{)}^2+({\lambda }_4-{\lambda }_3)\times \\ & \times ({\lambda }_2-{\lambda }_1){\mathrm{sn}}^2\left(\sqrt{({\lambda }_4-{\lambda }_2)({\lambda }_3-{\lambda }_1)}\theta ,m\right),\end{array}$$

$$u=V-\frac{j}{\rho },$$

$$\begin{array}{rl}j& =\frac{1}{8}(-{\lambda }_1-{\lambda }_2+{\lambda }_3+{\lambda }_4)\times \\ & \times (-{\lambda }_1+{\lambda }_2-{\lambda }_3+{\lambda }_4)({\lambda }_1-{\lambda }_2-{\lambda }_3+{\lambda }_4),\\ \theta & =x-Vt,V=\frac{1}{2}\sum _{i=1}^4{\lambda }_i,\\ m& =\frac{({\lambda }_2-{\lambda }_1)({\lambda }_4-{\lambda }_3)}{({\lambda }_4-{\lambda }_2)({\lambda }_3-{\lambda }_1)}.\end{array}$$

$$L=\frac{2K(m)}{\sqrt{({\lambda }_4-{\lambda }_2)({\lambda }_3-{\lambda }_1)}},$$

$$\rho ={\rho }_0-\frac{a_s}{{\cosh }^2(\sqrt{a_s}(x-V_{s}t))},$$

$$\begin{array}{rl}{\rho }_0& =\frac{1}{4}({\lambda }_4-{\lambda }_1{)}^2,a_s=({\lambda }_4-{\lambda }_2)({\lambda }_2-{\lambda }_1),\\ V_s& =\frac{1}{2}({\lambda }_1+2{\lambda }_2+{\lambda }_4).\end{array}$$

$$\frac{\partial {\lambda }_i}{\partial t}+v_i(\lambda )\frac{\partial {\lambda }_i}{\partial x}=0,i=1,2,3,4,$$

$$\begin{array}{rl}{v}_i(\lambda )& =\left(1-\frac{L}{{\partial }_{i}L}{\partial }_i\right)V,i=1,2,3,4,\\ {\partial }_i& \equiv \partial /\partial {\lambda }_i.\end{array}$$

$$\begin{array}{rl}{v}_1& =\frac{1}{2}\sum {\lambda }_i-\frac{({\lambda }_4-{\lambda }_1)({\lambda }_2-{\lambda }_1)K}{({\lambda }_4-{\lambda }_1)K-({\lambda }_4-{\lambda }_2)E},\\ v_2& =\frac{1}{2}\sum {\lambda }_i+\frac{({\lambda }_3-{\lambda }_2)({\lambda }_2-{\lambda }_1)K}{({\lambda }_3-{\lambda }_2)K-({\lambda }_3-{\lambda }_1)E},\\ v_3& =\frac{1}{2}\sum {\lambda }_i-\frac{({\lambda }_4-{\lambda }_3)({\lambda }_3-{\lambda }_2)K}{({\lambda }_3-{\lambda }_2)K-({\lambda }_4-{\lambda }_2)E},\\ v_4& =\frac{1}{2}\sum {\lambda }_i+\frac{({\lambda }_4-{\lambda }_3)({\lambda }_4-{\lambda }_1)K}{({\lambda }_4-{\lambda }_1)K-({\lambda }_3-{\lambda }_1)E},\end{array}$$

$$\begin{array}{rl}& v_1=\frac{3}{2}{\lambda }_1+\frac{1}{2}{\lambda }_4,v_4=\frac{3}{2}{\lambda }_4+\frac{1}{2}{\lambda }_1,\\ & v_2=v_3=\frac{1}{2}({\lambda }_1+2{\lambda }_2+{\lambda }_4),\end{array}$$

$$\begin{array}{rl}& v_1=\frac{3}{2}{\lambda }_1+\frac{1}{2}{\lambda }_2,v_2=\frac{3}{2}{\lambda }_2+\frac{1}{2}{\lambda }_1,\\ & v_3=v_4=2{\lambda }_4+\frac{({\lambda }_2-{\lambda }_1{)}^2}{2({\lambda }_1+{\lambda }_2-2{\lambda }_4)},\end{array}$$

$$x-v_i(\lambda )t=w_i(\lambda ),i=1,2,3,4,$$

$$\frac{1}{w_i-w_j}\frac{\partial w_i}{\partial {\lambda }_j}=\frac{1}{v_i-v_j}\frac{\partial v_i}{\partial {\lambda }_j},i\ne j.$$

$$w_i(\lambda )=\left(1-\frac{L}{{\partial }_{i}L}{\partial }_i\right)W,i=1,2,3,4.$$

$$\frac{{\partial }^{2}W}{\partial {\lambda }_i\partial {\lambda }_j}-\frac{1}{2({\lambda }_i-{\lambda }_j)}\left(\frac{\partial W}{\partial {\lambda }_i}-\frac{\partial W}{\partial {\lambda }_j}\right)=0.$$

$$W=\frac{{\lambda }^2}{\sqrt{{\prod }_{i=1}^4(\lambda -{\lambda }_i)}},$$

$$W=\sum _{k=0}^{\infty }\frac{W^{(k)}}{{\lambda }^k}=1+W^{(1)}\frac{1}{\lambda }+W^{(2)}\frac{1}{{\lambda }^2}+\dots ,$$

$$w_i^{(k)}=W^{(k)}+(2v_i-s_1){\partial }_i{W}^{(k)},$$

$$x-v_i(\lambda )t=\sum _{k=0}^n{A}_k{w}_i^{(k)}(\lambda ),i=1,2,3,4,$$

$${\rho }_t+(\rho u{)}_x=0,u_t+u{u}_x+{\rho }_x=0.$$

$${\lambda }_{\pm }=\frac{u}{2}\pm \sqrt{\rho },$$

$$\begin{array}{rl}& \frac{\partial {\lambda }_{+}}{\partial t}+v_{+}({\lambda }_{+},{\lambda }_{-})\frac{\partial {\lambda }_{+}}{\partial x}=0,\\ & \frac{\partial {\lambda }_{-}}{\partial t}+v_{-}({\lambda }_{+},{\lambda }_{-})\frac{\partial {\lambda }_{-}}{\partial x}=0,\end{array}$$

$$v_{+}=\frac{3}{2}{\lambda }_{+}+\frac{1}{2}{\lambda }_{-},v_{-}=\frac{1}{2}{\lambda }_{+}+\frac{3}{2}{\lambda }_{-},$$

$$\rho =\frac{1}{4}({\lambda }_{+}-{\lambda }_{-}{)}^2,u={\lambda }_{+}+{\lambda }_{-}.$$

$${\lambda }_{-}=-\sqrt{{\rho }_0}.$$

$$\frac{\partial {\lambda }_{+}}{\partial t}+\left(\frac{3}{2}{\lambda }_{+}-\frac{1}{2}\sqrt{{\rho }_0}\right)\frac{\partial {\lambda }_{+}}{\partial x}=0$$

$$x-\left(\frac{3}{2}{\lambda }_{+}-\frac{1}{2}\sqrt{{\rho }_0}\right)t=w({\lambda }_{+}).$$

Full text

Contribution for the JETP special issue in honor of L. P. Pitaevskii’s 85th birthday

I Introduction

It is well known that with the disregard of viscosity and dispersion effects, nonlinear waves experience “breaking,” i.e., after a certain critical instant, the formal solution to corresponding evolution equations becomes multi-valued (see, for example, LL-6 ). In classical gas dynamics, this problem is eliminated by taking into account weak dissipation effects so that instead of the multi-valuedness domain, a shock wave (i.e., a narrow region of transition from the flow with some values of parameters characterizing the flow to a flow with other values of parameters) appears in the solution. The width of this transition region is proportional to coefficients characterizing dissipative processes; in real conditions, this width is usually of the order of the mean free path of molecules in the gas. For this reason, it can be assumed in the macroscopic theory that this region is a discontinuity in the parameters of the flow of the medium, and when the medium passes through the discontinuity, the mass, momentum, and energy conservation laws must hold. The theory of shock waves formulated on this basis has been profoundly developed and has found numerous applications (see, for example, LL-6 ; CF-50 ).

In modern physics, however, flows of the medium in which dissipation processes can be disregarded in the first approximation are often considered, and then nonlinear wave breaking is eliminated by taking into account dispersion effects which lead to formation of dispersive shock waves (DSWs) (i.e., the evolving regions of the nonlinear flow of the medium) instead of the multi-valued domains. Such effects were studied for the first time in the theory of undular bores in a shallow water flow (see, e.g., bl-54 ), and the general nature of this phenomenon was realized by Sagdeev sagdeev , who indicated that wave breaking in dispersive wave systems leads to the formation of an extended wave structures connecting different states of the flow like a transition in a shock wave connects different states of the medium flow with predominant dissipation. In typical cases, a dispersive shock wave occupies a spatial region expanding with time so that this wave is a sequence of solitons at one of its edge and degenerates into a small-amplitude harmonic wave propagating with the corresponding group velocity at the other edge. The main theoretical approach to the description of DSWs was proposed in the classical work by Gurevich and Pitaevskii gp-73 based on the Whitham theory of modulation of nonlinear waves whitham . In this approach, a DSW is represented in the form of a modulated periodic solution to the corresponding nonlinear wave equation, and the slow evolution of the modulation parameters obeys the Whitham equations obtained by averaging of the conservation laws over fast oscillations of physical variables in the wave. Gurevich and Pitaevskii considered two typical problems of wave breaking, the evolution of the wave is described by the Korteweg-de Vries (KdV) equation. First, a complete analytic solution of the discontinuity decay was obtained in the case when the initial distribution at the breaking instant has a sharp jump. Second, they found the main characteristics of the DSW in vicinity of the breaking point when the initial distribution is described by a cubic parabola. Later, Potemin potemin obtained a full analytic solution to this problem (see also kamch ). The Gurevich and Pitaevskii approach to the DSW theory was developed further and was extended to other equations (see, for example, review eh-16 ).

One of important applications of the DSW theory is the dynamics of the Bose-Einstein condensate, which is described by the Gross-Pitaevskii equation gross ; pit ; for simplicity, we write here this equation in the standard dimensionless form for a 1D flow of the condensate:

[TABLE]

where $\psi$ is the “wave function” of the condensate flowing along the $x$ coordinate; we assume that the interaction between atoms is repulsive, which ensures the stability of its homogeneous state. The theory of Eq. (1) was considered in a huge number of publications. In particular, its solution in the form of a dark soliton was obtained in tsuzuki-1971 , and periodic solutions were obtained, for example, in gk-87 . The integrability of Eq. (1) by the inverse scattering transform method was established in zs-73 and this approach was used in fl-86 ; pavlov-87 for deriving the modulation equations. Finally, the problem of the initial discontinuity decay was analyzed in gk-87 ; eggk-95 , and typical wave breaking was studied in kku-02 . The apparatus developed in these works was applied to the description of the DSW dynamics in the Bose-Einstein condensate.

The DSWs in the condensate were observed experimentally for the first time in cornell ; hoefer , where a shock wave was formed under the action of laser radiation repelling the condensate. The interpretation of these observations as DSWs was reported in kgk-04 , and this experiment was described in hoefer using the theory formulated in gk-87 ; eggk-95 , where it was assumed that a DSW is formed as a result of emergence of a discontinuity. Although such a discontinuity can be formed when the flow of the condensate is induced by a piston moving at a constant velocity hae-08 , such a case is nevertheless quite specific, and the wave propagating to the bulk of a stationary medium in typical situations breaks from a profile with a certain root singularity rather than a sharp discontinuity. For instance, a singularity in the form of a square root appears in the case of the uniformly accelerated motion of the piston kk-10 , and it is clear even from this example that the actual flow of the condensate may have a quite arbitrary singularity at the instant of wave breaking. Here, we consider the DSW formation during wave breaking with an initial singularity of the type $(-x)^{1/n}$. A detailed theory will be developed for $n=2$ and $n=3$, and important DSW characteristics (such as the laws of motion of its edges) will be obtained for an arbitrary integer $n>1$. The theory developed here forms the basis for describing quite general forms of condensate flow with wave breaking.

II Gurevich-Pitaevskii method

Let us first write the basic relations of the Gurevich-Pitaevskii method in the DSW theory as applied to the dynamics of the Bose-Einstein condensate, which obeys the Gross-Pitaevskii equation (1). It is convenient to represent the periodic solutions to this equation in terms of more transparent physical variables by performing the substitution

[TABLE]

so that after the separation of the real and imaginary parts, we obtain the system

[TABLE]

Here, $\rho(x,t)=|\psi(x,t)|^{2}$ is the condensate density and $u(x,t)$ is the condensate flow velocity. The periodic solution can be written in the form

[TABLE]

[TABLE]

where

[TABLE]

This solution depends on four parameters $\lambda_{1}\leq\lambda_{2}\leq\lambda_{3}\leq\lambda_{4}$, in terms of which the main characteristics of the wave can be expressed. In particular, the wavelength is given by

[TABLE]

where $K(m)$ is the complete elliptic integral of the first kind. In the limit $\lambda_{3}\to\lambda_{2}$, when $m\to 1$ and $L\to\infty$, the periodic wave is transformed into the soliton solution

[TABLE]

where background density $\rho_{0}$, along which the dark soliton propagates, its amplitude $a_{s}$ and velocity $V_{s}$ are given by

[TABLE]

In the opposite limit $\lambda_{3}\to\lambda_{4}$, when $m\to 0$, the wave amplitude tends to zero, and it is transformed into a linear harmonic wave propagating over the background with constant density $\rho_{0}$.

DSW parameters $\lambda_{i}$ become slow functions of $x$ and $t$, which change little in one wavelength $L$. Therefore, we can average the conservation laws for Eq. (1) over fast oscillations in the wave and obtain as a result the Whitham equations for modulation parameters $\lambda_{i}$. These equations can be written in the form

[TABLE]

where the characteristic velocities are given by

[TABLE]

The substitution of expressions (6) and (7) into these formulas gives the following expressions for velocities:

[TABLE]

where $E=E(m)$ is the elliptic integral of the second kind. Variables $\lambda_{i}$ are known as Riemann invariants of the system of the Whitham modulation equations. We will also need the limiting expressions for these velocities at the DSW edges. At the soliton edge, where $\lambda_{3}=\lambda_{2}$ and $m=1$, we have

[TABLE]

while at the small-amplitude edge for $\lambda_{3}=\lambda_{4}$ and $m=0$ we have

[TABLE]

(we will not need the expressions for the analogous limit $\lambda_{1}=\lambda_{2}$).

In the generalized hodograph method tsarev the solutions to Eqs. (10) are sought in the form

[TABLE]

where $v_{i}(\lambda)$ are velocities (12) and $w_{i}(\lambda)$ are the sought functions. If these functions have been determined, $x=x(\lambda)$ and $t=t(\lambda)$ turn out to be the functions of parameters $\lambda_{i}$. Since these functions must be inverted and the modulation parameters must become functions $\lambda_{i}=\lambda_{i}(x,t)$, the functions $w_{i}$ cannot be independent of one another. Differentiating Eq. (15) with respect to $\lambda_{j}$, $j\neq i$, and eliminating $t$ from all pair combinations of the resultant relations, we arrive at the system of the Tsarev equations

[TABLE]

In view of their symmetry in $v_{i}$ and $w_{j}$, it is natural to seek their solution in the form analogous to (11) (see gke-92 ; wright ; tian ):

[TABLE]

Then Eqs. (16) are transformed into the system of Euler-Poisson equations $(i\neq j)$

[TABLE]

For our purposes, it is sufficient to know the set of solutions obtained from the generating function

[TABLE]

that satisfies Eqs. (18) for any $\lambda$. The expansion of function (19) in inverse powers of $\lambda$ gives

[TABLE]

where $W^{(k)}=W^{(k)}(\lambda_{1},\lambda_{2},\lambda_{3},\lambda_{4})$ are the required particular solutions to system (18). As a result, using expression (17) we obtain the set of functions $w_{i}^{(k)}$, which give the solutions to Whitham equations (10):

[TABLE]

so that $w_{i}^{(0)}=1,$ $w_{i}^{(1)}=v_{i}$. The Euler-Poisson equation is linear in $W$, like expressions (21) that are linear in $W^{(k)}$; therefore, any of their linear combinations also gives the solution

[TABLE]

where the number of terms $n$ and constant coefficients $A_{k}$ are chosen in accordance with the conditions of the problem.

Let us now show that the above expression of the DSW theory in the Gurevich and Pitaevskii approach make it possible to solve the problem of the DSW formation during wave breaking in the Bose-Einstein condensate.

III Dispersionless limit

Until the instant of breaking, the distributions of density $\rho(x,t)$ and flow velocity $u(x,t)$ are smooth functions of spatial coordinate $x$. Moreover, even after breaking, the DSW occupies a finite spatial region and its edges at the matching points with the smooth distributions must be determined as a part of the solution of the wave breaking problem. In the case of quite smooth functions $\rho$ and $u$ the terms with a large number of derivatives in system (3) can be omitted, which means the disregard of the dispersion effects; therefore, the evolution of smooth distributions can be described by the dispersionless limit equations

[TABLE]

These equations coincide with the “shallow water” equations (see LL-6 ), equivalent to the gas dynamic equations with adiabatic exponent $\gamma=2$. Therefore, their solutions can be obtained using well-known classical methods.

Equations (23) can be transformed to diagonal form by introducing the Riemann invariants

[TABLE]

so that

[TABLE]

where

[TABLE]

here, $\rho$ and $u$ can be expressed in terms of $\lambda_{\pm}$ by the formulas

[TABLE]

Generally, both Riemann invariants $\lambda_{\pm}$ are functions of $x$ and $t$. We are interested, however, in the problem in which a wave propagates to the bulk of the condensate at rest with constant density $\rho_{0}$. It is known LL-6 , that only a flow in the form of a simple wave, in which one of the Riemann invariants is constant, can border on such a state of the gas. Assuming for definiteness that the wave propagates to the right, we can conclude that Riemann invariant $\lambda_{-}$ must be constant and, hence, must have the same value as in the stationary medium bordering the wave:

[TABLE]

Then the second equation in (25) is satisfied automat- ically, while the first equation is transformed into the Hopf equation

[TABLE]

with the well-known general solution

[TABLE]

If function $w(\lambda_{+})$ is known, this solution is defined by the dependence $\lambda_{+}=\lambda_{+}(x,t)$, and this function must be joined at the boundary with the stationary condensate with value $\lambda_{+}=\sqrt{\rho_{0}}$ in this spatial region.

We are interested in the situation when the smooth solution for $\lambda_{+}$ at the breaking instant tends to its boundary value $\sqrt{\rho_{0}}$ as a root function of $x$. Choosing the coordinate system and its origin so that breaking of Riemann invariant $\lambda_{+}$ occurs at instant $t=0$ at the origin $x=0$, we obtain the dependence

[TABLE]

where, to simplify calculations, the units of measurements of length and time are chosen so that the coefficient on the right-hand side be equal to unity. Therefore, the dependence of $\lambda_{+}$ on $x$ for $t<0$ has no singularities, while, at $t=0$ the root singularity appears,

[TABLE]

and the dependence of $\lambda_{+}$ on $x$ for $t>0$ becomes multi-valued (see Fig. 1).

IV Dispersive shock wave

At instant $t$ the DSW occupies the spatial region

[TABLE]

matching to the smooth solution (30) at its boundary point. Comparison of velocity $v_{+}$ in expression (26) with limiting expression (13) shows that the DSW is transformed into expression (30) for $\lambda_{4}=\lambda_{+}$ at boundary $x=x_{-}(t)$; coefficients $A_{k}$ in this case must be chosen so that the right-hand side of relation (22) with $i=4$ be equal to the right-hand side of relation (30). Further, the solution to the Whitham equations is transformed into a harmonic wave at the small-amplitude edge $x=x_{+}(t)$, if $\lambda_{2}=\lambda_{+}=\sqrt{\rho_{0}}$ along the entire DSW, and we have $\lambda_{3}=\lambda_{4}$ at point $x_{+}(t)$. Since $\lambda_{-}=-\sqrt{\rho_{0}}$ at both edges of the SDW, the condition of matching of $\lambda_{1}$ to $\lambda_{-}$ at the DSW edges can be satisfied by setting $\lambda_{1}=-\sqrt{\rho_{0}}$ along the DSW. Therefore, Whitham equations (10) with $i=1,2$ are satisfied by constant solutions

[TABLE]

and only two Riemann invariants $\lambda_{3}$ and $\lambda_{4}$, which satisfy the boundary conditions

[TABLE]

and

[TABLE]

vary along the DSW. These conditions define the solution completely. As a result, the dependence of the Riemann invariants on coordinate $x$ for a fixed value of $t$ has the form shown in Fig. 2. It should be noted that waves with two variable Riemann invariants were called quasi-simple and were studied for the first time in gkm-89 in the theory of the KdV equation. Taking relations (33) into account, we can find the first three coefficients in Eq. (20) in the form

[TABLE]

the knowledge of these coefficients is sufficient for analyzing typical cases with $n=2$ and $n=3$.

IV.1 The case with $n=2$

Since formulas (36) are polynomial, quadratic function $-(\lambda_{4}-\lambda_{2})^{2}$ can be written in the form of a linear combination of the first three expressions (36) with coefficients

[TABLE]

Then formulas (22) with $i=3,4$ and with these values of the coefficients define implicitly the dependencies of $\lambda_{3}$ and $\lambda_{4}$ on $x$ and $t$, which solves in principle the problem in this particular case. At the soliton edge, these formulas are transformed into

[TABLE]

which gives the relation between $t$ and $\lambda_{4}$ at this boundary:

[TABLE]

Substituting the resultant value of $\lambda_{4}$ into any relation from (38) we obtain the law of motion of the soliton edge:

[TABLE]

At the small-amplitude edge for $\lambda_{3}=\lambda_{4}$ both formulas(22) with $i=3,4$ are transformed into the same relation

[TABLE]

This edge moves with the group velocity corresponding to wavenumber $k=2\pi/L$ with $L=2\sqrt{\lambda_{4}^{2}-\lambda_{2}^{2}}$, which gives for the Bogoliubov dispersion law $\omega=k\sqrt{\lambda_{2}^{2}+k^{2}/4}$ the expression $d\omega/dk=2\lambda_{4}^{2}-\lambda_{2}^{2}/\lambda_{4}$. Therefore, the differentiation of expression (41) with respect to $t$ with account of $dx_{+}/dt=d\omega/dk$ determines the dependence of $t$ on the value of $\lambda_{4}$ at this boundary. Introducing parameter $y=\lambda_{4}/\lambda_{2}$ instead of $\lambda_{4}$, we can write this dependence in the form

[TABLE]

and its substitution into relation (41) gives

[TABLE]

Formulas (42) and (43) define parametrically the law of motion $x_{+}=x_{+}(t)$ of the small-amplitude edge. For $t\sim\lambda_{0}y\gg\lambda_{2}^{2}=\rho_{0}$ this law of motion asymptotically takes the simple form

[TABLE]

It should be noted that analogous expressions obtained by solving the problem of motion of the condensate under the action of a uniformly accelerated piston kk-10 can be transformed to the expression obtained above after the transfer of the breaking point to the origin and subtracting the breaking time from $t$.

IV.2 Case with $n=3$

In this case, the calculations are performed analogously. The right-hand sides of formulas (22) with $i=3,4$ now contain function $W_{3}$ and condition $(\ref{eq34})$ gives the values of the coefficients

[TABLE]

At the soliton edge, solution (22) is transformed into

[TABLE]

which gives

[TABLE]

The substitution of this relation into (46) gives the law of motion of the soliton edge:

[TABLE]

At the small-amplitude edge, formulas (22) with $i=3.4$ are transformed into

[TABLE]

The matching condition for the law of motion with the group velocity of a linear wave gives

[TABLE]

and the substitution into expression (49) gives

[TABLE]

For asymptotically long times $t\sim\lambda_{0}^{2}y^{2}\gg\lambda_{2}^{2}=\rho_{0}^{2}$ we obtain

[TABLE]

The laws of motion of the DSW edges as functions of time for $n=2,3$ is shown in Fig. 3. For short times, the motion with the velocity of sound $\sqrt{\rho_{0}}$, prevails, while, for long times, a transition to asymptotic laws occurs.

General formulas (22) together with the specific expressions for functions $w_{i}^{(k)}$, obtained above with the help of formulas (21), (36) and coefficients (37) ($n=2$) and (45) ($n=3$), make it possible to calculate $\lambda_{3}$ and $\lambda_{4}$ as functions of $x$ and $t$, so that their substitution into expressions (4) and (5) gives the density and velocity distribution profiles in a DSW (Fig. 4). The envelopes of the condensate density in the DSW are calculated by the formulas

[TABLE]

Although the formulas are complicated with increasing exponent $n$, important DSW characteristics such as the laws of motion of the edges can be determined without detailed analysis of the complete solution (at least, in the asymptotic limit $t\gg\rho_{0}^{n}$). In the next sections, we will consider this problem.

V Law of motion of the soliton edge

Combining relations (34) with the limiting expression of formula (21), we can write the boundary condition at the soliton edge in the form of a differential equation for the function $W=W(-\lambda_{2},\lambda_{2},\lambda_{2},\lambda_{4})$, which depends only on $\lambda_{4}$:

[TABLE]

the solution to this equation is

[TABLE]

where the integration constant is chosen so that time $t$ in subsequent formulas tends to zero for $\lambda_{4}\to\lambda_{2}$. Since $w_{3}=W$ at this boundary, formulas (22) give

[TABLE]

which is in conformity with relations (38) ($n=2$) and (46) ($n=3$). This gives

[TABLE]

and

[TABLE]

These expressions generalize the formulas obtained above to arbitrary integer values of $n>1$.

It should be noted that, in fact, we can find the law of motion of the soliton edge for an arbitrary monotonic dependence of the initial distribution of Riemann invariant $\lambda_{+}$ of the form $w=w(\lambda_{+}-\sqrt{\rho_{0}})$, by resorting to the considerations used in gkm-89 for deriving the law of motion of the small-amplitude edge for wave breaking in the theory of the KdV equation. Indeed, Whitham equations (10) with $i=3,4$ in the classical hodograph method are transformed into linear differential equations for functions $x=x(\lambda_{3},\lambda_{4})$ and $t=t(\lambda_{3},\lambda_{4})$, one of which for $\lambda_{3}\to\lambda_{2}$ becomes

[TABLE]

On the other hand, the solution at this boundary must match to the smooth solution, which gives

[TABLE]

The differentiation of this relation with respect to $\lambda_{4}$ leads to one more differential equation

[TABLE]

Eliminating $\partial x/\partial\lambda_{4}$ from Eqs. (58) and (59), we obtain the differential equation for $t$, the integration of which gives $(z=\lambda_{4}-\lambda_{2})$

[TABLE]

and, hence,

[TABLE]

These formulas specify the parametric dependence $x_{-}=x_{-}(t)$.

VI Law of motion of the small-amplitude edge

Formulas (41) $(n=2)$ (49) and $(n=3)$ have a simple structure leading to the assumption that for $m\to 0$ $(\lambda_{3}\to\lambda_{4})$ and integer $n$ functions $w_{3}$ and $w_{4}$ must pass to the right-hand side of the relation

[TABLE]

Let us prove this formula in the asymptotic limit $t\to\infty$ $(\lambda_{4}\to\infty)$, when the terms with $\lambda_{2}/\lambda_{4}$ can be neglected. We note that in the limit $\lambda_{2}\to 0$ generating function (19) can be reduced to the generating function of the Legendre polynomials (see, for example, ww )

[TABLE]

i.e.,

[TABLE]

Using the recurrent formula for the derivative of the Legendre polynomial (see ww ) we can easily prove that

[TABLE]

in this approximation. To evaluate function $w_{4}^{(n)}$ in the limit $m\to 0$ we will prove that the following relation holds for $\lambda_{3}\to\lambda_{4}$:

[TABLE]

For this purpose, we note that the argument of the Legendre polynomial in expression (63), which is the ratio of the arithmetic mean to the geometric mean, attains its maximal value for $\lambda_{3}=\lambda_{4}$ and, hence, is quadratic in the small difference $\lambda_{4}-\lambda_{3}=\varepsilon$, so that $P_{n}(1)=1$, i.e., $W_{n}(\lambda_{4},\lambda_{4})=\lambda_{4}^{n}$. Consequently, in the first order in $\varepsilon$ we obtain

[TABLE]

In addition, we note that for $\lambda_{2}\to 0$ only the highest term with $k=n$ is left in the sum in expression (22) since coefficients $A_{k}$ for $k<n$ contain powers of $\lambda_{2}$ as factors. Therefore, with account of relation. $2(V-v_{4})\approx-2\lambda_{4}$ we obtain

[TABLE]

For determining $w_{4}$ completely for $m=0$ it remains to find $A_{n}$, which can easily be done using the expansion of the higher term in the Legendre polynomial into a series for $\lambda_{3}\to\lambda_{2}=0$ (see ww ):

[TABLE]

Then the condition of matching to the smooth solution gives

[TABLE]

i.e.

[TABLE]

Differentiating the expression

[TABLE]

with respect to $\lambda_{4}$ provided that $\partial x/\partial\lambda_{4}=0$ for a fixed $t$, which is equivalent to the matching condition with the group velocity at the small-amplitude edge, we obtain

[TABLE]

Substituting $\lambda_{4}$ obtained from this expression into (68) we obtain the law of motion of the small-amplitude edge in the asymptotic regime:

[TABLE]

This formula naturally reproduces the above asymptotic laws (44) for $n=2$ and (52) for $n=3$.

Expression (68) confirms the validity of formula (62) in the limit $\lambda_{2}\to 0$. Assuming that this formula also holds for a finite $\lambda_{2}$, we obtain the law of motion of the small-amplitude boundary in parametric form:

[TABLE]

It should be noted that in contrast to the theory of the KdV equations, the self-similar regime of motion of the boundaries is realized only asymptotically for long times $t\gg\rho_{0}^{n-1}$, when the velocity of motion is much higher than the velocity of sound in the background distribution. However, the limiting transition to $\rho_{0}\to 0$ is impossible in the expressions describing the wave profile since the magnitude of $m$ in elliptic functions vanishes for $\lambda_{2}=-\lambda_{1}\to 0$.

VII Conclusion

Thus, the approach developed by Gurevich and Pitaevskii makes it possible to analyze in detail the process of DSW formation during wave breaking in the Bose-Einstein condensate, the dynamics of which obeys the Gross-Pitaevskii equation. The developed theory is applicable to the initial stage of the process, in which the smooth part of the profile can be treated as a monotonic function of the coordinate. It should be noted, however, that the theory of quasi-simple waves must also describe the asymptotic stage of the evolution of a finite-duration pulse since, analogously to the simple wave theory, the initial pulse splits with time into two pulses, in each of which two of four Riemann invariants again remain constant. Therefore, the Gurevich-Pitaevskii approach supplemented with the generalized hodograph method and modern method for deriving the Whitham modulation equations remains a powerful tool for investigating dispersion shock waves, which are of considerable interest for modern nonlinear physics.

Acknowledgements.

I thank M. Isoard, S. K. Ivanov and N. Pavloff for useful discussions.