# Asymptotic theory of not completely integrable soliton equations

A. M. Kamchatnov  Institute of Spectroscopy, Russian Academy of Sciences,
Troitsk, Moscow, 108840, Russia<EMAIL_ADDRESS>
###### Abstract

We develop the theory of transformation of intensive initial nonlinear wave
pulses to trains of solitons emerging at asymptotically large time of
evolution. Our approach is based on the theory of dispersive shock waves in
which the number of nonlinear oscillations in the shock becomes the number of
solitons at the asymptotic state. We show that this number of oscillations,
which is proportional to the classical action of particles associated with the
small-amplitude edges of shocks, is preserved by the dispersionless flow. Then
the Poincaré-Cartan integral invariant is also constant and therefore it
reduces to the quantization rule similar to the Bohr-Sommerfeld quantization
rule for linear spectral problem associated with completely integrable
equations. This rule yields a set of ‘eigenvalues’ which are related with the
asymptotic solitons’ velocities and other their characteristics. Our
analytical results agree very well with the results of numerical solutions of
the generalized nonlinear Schrödinger equation.

###### pacs:

47.35.Jk, 47.35.Fg, 02.30.Jr

> In many typical situations an intensive enough initial wave pulse evolves
> eventually to some number of solitons, and therefore it is very important to
> be able to predict velocities of these asymptotic solitons and other their
> characteristics. If the nonlinear wave equation under consideration is
> completely integrable, then this problem can be solved by finding the
> eigenvalues of the linear spectral problem associated with this equation.
> The theory becomes especially effective when the number of solitons is
> large, so that the quasiclassical asymptotical method can be applied to the
> spectral problem. However, in case of not completely integrable equations a
> similar approach was only known for a special class of initial conditions
> when they correspond to simple waves of the dispersionless approximation. In
> this paper, we generalize this theory to arbitrary smooth enough initial
> pulses described by two wave variables. We obtain the generalized Bohr-
> Sommerfeld quantization rule which defines a set of ‘eigenvalues’
> corresponding to given initial conditions and show how these eigenvalues are
> related with physical parameters of asymptotic solitons. The theory agrees
> very well with numerical simulations and sheds new light on the
> quasiclassical limit of completely integrable equations.

##  I Introduction

If a nonlinear wave system supports propagation of solitons, then an intensive
enough initial pulse evolves eventually into a certain number $N$ of solitons
and some amount of linear radiation. In view of universality of this
phenomenon, it is very important to be able to predict parameters of emerging
solitons (e.g., their velocities and amplitudes) from given initial
distributions of the wave variables. This problem admits a straightforward
solution in case of completely integrable nonlinear wave equations; see, e.g.,
Refs. nmpz-80 ; as-81 ; newell-85  and references therein. Such equations are
associated with some linear spectral problems where the wave distributions
play the role of ‘potentials’. Though these ‘potentials’ evolve with time
according to the wave equation under consideration, the spectrum of the
associated linear problem does not change, i.e. the evolution is isospectral,
and each discrete eigenvalue corresponds to a certain soliton at the
asymptotic stage of evolution at $t\to\infty$. Hence, if we find the spectrum
of the linear problem for given initial wave distributions, then the discrete
eigenvalues provide all necessary information about solitons’ parameters at
asymptotically large time $t$.

This theory considerably simplifies if the final number of solitons is large,
$N\gg 1$, and the initial distributions are smooth enough, so we can apply the
quasiclassical WKB method to the linear spectral problem. For example, in case
of Korteweg-de Vries (KdV) equation associated with the stationary Schrödinger
equation ggkm-67 , the well-known Bohr-Sommerfeld quantization rule readily
gives approximate values of the discrete eigenvalues very simply related with
solitons’ velocities and amplitudes karpman-67 ; karpman-73 . As usual in the
WKB method, the asymptotic formulas give quite accurate results even for small
quantum numbers. Similar Bohr-Sommerfeld quantization rule was derived in
Refs. JLML-99 ; kku-02  for the Zakharov-Shabat spectral problem associated
with nonlinear Schrödinger (NLS) equation zs-73 . However, when the nonlinear
wave equation under consideration is not completely integrable, such a direct
approach becomes impossible.

To find a more general method for calculation of parameters of asymptotic
solitons, we need to consider in some more detail the process of
transformation of an initially smooth pulse to a sequence of solitons. One can
distinguish in this process three typical stages. At first, the pulse’s
profile gradually steepens due to nonlinear effects, but it keeps a smooth
enough form without any oscillations, so its evolution obeys with high
accuracy the dispersionless (hydrodynamic) limit of nonlinear wave equations
under consideration. When the profile’s slope becomes very large, the
dispersion effects start to influence on the evolution of the wave, so
oscillations are generated after the so-called wave breaking moment of time.
Such regions of nonlinear oscillations are called dispersive shock waves
(DSWs) and, according to Gurevich and Pitaevskii gp-73 , they can be
approximated by modulated periodic solutions of the nonlinear wave equations.
Evolution of the modulation parameters obeys the Whitham modulation equations
whitham ; Whitham-74 . In Gurevich-Pitaevskii approximation, a DSW occupies a
finite expanding region of space. At one edge of the DSW solitons are
gradually formed and the opposite small-amplitude edge propagates with some
group velocity along a smooth part of the evolving pulse. Thus, at the second
stage of evolution the wave structure consists of one or two DSWs adjoined
with smooth parts of the pulse and the DSWs expand gradually over the whole
pulse. Numbers of oscillations in DSWs increase with time due to difference
between phase and group velocities of linear waves at the small-amplitude
edges where oscillations enter into the DSWs regions gp-87 . At the last third
stage of evolution of an initially localized pulse with a finite length its
smooth parts disappear and the number of oscillations in DSWs is stabilized,
so the solitons acquire their asymptotic values of parameters as $t\to\infty$,
when distances between solitons become much greater than their widths. Thus,
asymptotic solitons form gradually from small-amplitude oscillations which
enter into the DSW region at its small-amplitude edge, so we have to consider
propagation of this edge in some more detail.

Let the nonlinear wave dynamics under consideration be described by two
variables which for definiteness we shall call “density” $\rho$ and “flow
velocity” $u$. The periodic solutions $\rho=\rho(x,t),u=u(x,t)$ are periodic
functions of the phase $\theta$ which in non-modulated case has the form
$\theta=kx-\omega t=k(x-Vt)$, where $V$ is the phase velocity of the traveling
wave. In a slightly modulated wave, the wave number $k$ and the frequency
$\omega$ become slow functions of $x$ and $t$. Since locally the wave can
still be considered in the leading approximation as uniform, they can be
defined as

$k=\theta_{x},\qquad\omega=-\theta_{t}$ (1)

for some phase $\theta(x,t)$, and then one of the modulation equations can be
written in the form of ‘the number of waves’ conservation law whitham ;
Whitham-74

$k_{t}+\omega_{x}=0,$ (2)

where $k$ plays the role of density of waves and $\omega$ is their flux.

At the small-amplitude edge we can distinguish the short wavelength
oscillation at the length scale $\sim 2\pi/k$ and the slowly changing
“background distributions” $\rho_{b}(x,t),u_{b}(x,t)$ which change at the
length scale of the initial pulse length $\sim l$. Hence we can write
$\rho(x,t)=\rho_{b}(x,t)+\rho^{\prime}(x,t),u(x,t)=u_{b}(x,t)+u^{\prime}(x,t)$
where the smooth functions $\rho_{b}(x,t),u_{b}(x,t)$ obey the motion equation
in the dispersionless (hydrodynamic) limit obtained by neglecting term with
higher order space and/or time derivatives, and the small-amplitude
oscillations $\rho^{\prime}(x,t),u^{\prime}(x,t)$ obey the equations
linearized with respect to deviations from the background flow. Since at the
length scale $\sim 2\pi/k$ the background flow can be considered as constant,
we obtain in the main approximation from the linearized equations with
constant coefficients the harmonic wave solutions
$\rho^{\prime},u^{\prime}\propto\exp[i(kx-\omega t)]$ with the dispersion
relation $\omega=\omega(k,\rho_{b},u_{b})$. A slightly modulated linear wave
packet propagating along a slowly changing background still consists of
harmonics with the same dispersion law which satisfies the number of waves
conservation law (2), but $\rho_{b},u_{b}$ become in this case the slowly
changing solutions $\rho=\rho(x,t),u=u(x,t)$ of the dispersionless equations
(see Ref. Whitham-74 ). In Gurevich-Pitaevskii theory gp-73  of DSWs, the
small-amplitude edge is represented by such a linear wave packet and according
to the conservation law (2) oscillations enter into the DSW region where their
amplitudes grow from zero at the small-amplitude edge to soliton-like
propagation at the opposite soliton edge. Since the small-amplitude edge
propagates with the group velocity $v_{g}=\partial\omega(k,\rho,u)/\partial
k$, along its path the flux is Doppler-shifted, so the number of oscillations
inside the DSW region changes with time as (see Refs. gp-87 ; kamch-21a )

$\frac{dN}{dt}=\frac{1}{2\pi}\left|k\frac{\partial\omega}{\partial
k}-\omega\right|.$ (3)

The expression in the right-hand side coincides, up to a constant factor, with
Lagrangian of a point particle associated with the small-amplitude edge of the
DSW, if in accordance with the well-known optics-mechanical analogy (see,
e.g., Ref. lanczos ) we identify $k$ and $\omega$ with the particle’s momentum
and Hamiltonian, correspondingly. Then the number $N$ of solitons formed from
an intensive initial pulse is equal to the mechanical action produced by an
associated with the packet point-like particle,

$N=\frac{S}{2\pi}=\frac{1}{2\pi}\int(kdx-\omega dt),$ (4)

where $dx=v_{g}dt$ and integration is taken over the packet’s path. If the
initial distributions $\rho_{0}(x),u_{0}(x)$ correspond to unidirectional
simple wave propagation of the background pulse, that is there exists a
functional dependence between $\rho$ and $u$ corresponding to a constant value
of one of the Riemann invariants of the dispersionless equations LL-6 , then
the wave number $k$ can be expressed as a function of the value of $\rho$ at
the point where the small-amplitude edge is located at this moment of time,
$k=k(\rho)$. Then we have also $\omega=\omega(k,\rho)$, where $\rho=\rho(x,t)$
is a solution of the Hopf equation

$\rho_{t}+V_{0}(\rho)\rho_{x}=0$ (5)

to which the dispersionless equations can be reduced, $k=k(\rho)$ is the
solution of the equation

$\frac{dk}{d\rho}=\frac{\partial\omega/\partial\rho}{V_{0}-\partial\omega/\partial
k}$ (6)

with the initial condition

$k(0)=0,$ (7)

where it is assumed that wave breaking occurs at the point with $\rho=0$.
Equation (6) was first obtained in Ref. el-05  by means of analysis of Whitham
equations for DSWs at their small-amplitude edge and then it was derived from
the Hamilton equations for propagation of linear wave packets along simple
wave pulses in Refs. kamch-20a ; ks-21 . When the function $k=k(\rho)$ is
found, then, as was shown in Ref. kamch-19 , one can find the path $x=x(t)$ of
the small-amplitude edge of a DSW and, consequently, calculate the integral in
Eq. (4) for some particular nonlinear wave equations (see Refs. kamch-20a ;
kamch-21 ; cbk-21 ). In this case the final result can be written in the form

$N=\frac{1}{\pi}\int_{-\infty}^{\infty}k(\rho_{0}(x))dx,$ (8)

where $\rho_{0}(x)$ is the initial distribution of $\rho$. Such a formula was
suggested earlier in Refs. egkkk-07 ; egs-08  under supposition that one can
extend solution of the Whitham equations to the unmodulated yet part of the
pulse in such a way, that the total number of oscillations in DSW and
unmodulated part remains constant during evolution of the pulse.
Generalization of this rule to more general solutions
$k=k(\rho,\overline{q})$, $\overline{q}$ is an integration constant in a
solution of Eq. (6), yields also parameters of solitons at the asymptotic
stage of evolution egs-08 . This theory agrees with numerical solutions of
nonlinear wave equations and with experimental results obtained in Ref.
mfweh-20  for a viscous fluid conduit.

So far this theory was limited to an important but particular case of an
initial simple-wave pulse evolving eventually to a number of solitons and the
aim of this paper is to generalize this approach to the general case of
nonlinear wave evolutions describes by two variables $\rho(x,t),u(x,t)$ and
admitting solitonic propagation of waves. To this end, we notice that equality
of two expressions (4) and (8) can be interpreted as a special case of
preservation of a constant value of the Poincaré-Cartan integral invariants
poincare ; cartan  for a mechanical system with the phase space $(x,k)$ and
Hamiltonian $\omega(k,x,t)$. However, in this special case the tube of
trajectories preserving the invariant is defined by solutions of Eq. (5)
rather than by solutions of the Hamilton equations under consideration. In
this paper, we extend this definition of Poincaré-Cartan integral invariants
to general hydrodynamic flows for two variables $\rho,u$ not restricted by the
condition that one of the Riemann invariants is constant and arrive at
generalization of formula (8)

$N=\frac{1}{\pi}\int_{-\infty}^{\infty}k(\rho_{0}(x),u_{0}(x),q_{N})dx,$ (9)

where $\rho_{0}(x),u_{0}(x)$ are the initial distributions of these two
variables and $q_{N}$ is an integration constant in the solution of equations
which define the function $k=k(\rho,u)$ of two variables. One can say that
$q_{N}$ corresponds to the $N$-th soliton in the asymptotic distribution of
solitons. Then if we decrease $N$ by one, then we get the value $q_{N-1}$
corresponding to the $(N-1)$-th soliton, and in this way we relate each $n$-th
soliton with the corresponding value $q_{n}$ of the generalized ‘Bohr-
Sommerfeld quantization’ rule (9), which extends the above mentioned method
karpman-67 ; karpman-73 ; JLML-99 ; kku-02  of finding solitons parameters for
completely integrable wave equations to not completely integrable equations.
Relationship of the parameters $q_{n}$ with such physical parameters of
solitons as, for example, their velocities is established by means of the
Stokes remark lamb ; stokes  (see also ai-77 ; dkn-03 ; kamch-20a ) that
linear waves and exponentially small soliton tails obey the same linearized
equations, so the expressions for their phase or soliton velocities are
transformed to each other by the replacement
$k\leftrightarrow-i\widetilde{k}$, where $\widetilde{k}$ is the inverse half-
width of a soliton.

We illustrate the formulated here approach by its application to the
generalized NLS (gNLS) equation and confirm its accuracy by comparison with
numerical solutions.

##  II Poincaré-Cartan integral invariant and Bohr-Sommerfeld quantization
rule

Let the hydrodynamic (dispersionless) equations for the background flow
evolution can be written in the form of equations of the compressible fluid
dynamics,

$\rho_{t}+(\rho u)_{x}=0,\qquad u_{t}+uu_{x}+\frac{c^{2}}{\rho}\rho_{x}=0,$
(10)

where $c=c(\rho)$ is the ‘sound velocity’ related with the density $\rho$ by
means of equation of state $p=p(\rho)$, $c^{2}=dp/d\rho$, $p$ is the pressure.
The characteristic velocities of this system,

$v_{+}=u+c,\qquad v_{-}=u-c,$ (11)

correspond to the sound waves propagating upstream and downstream the
background flow with velocity $u$. Equations (10) can be transformed to the
Riemann diagonal form

$\frac{\partial r_{+}}{\partial t}+v_{+}\frac{\partial r_{+}}{\partial
x}=0,\quad\frac{\partial r_{-}}{\partial t}+v_{-}\frac{\partial
r_{-}}{\partial x}=0,$ (12)

for the variables

$r_{\pm}=\frac{u}{2}\pm\frac{1}{2}\int_{0}^{\rho}\frac{cd\rho}{\rho}$ (13)

called Riemann invariants. If we substitute here the function $\rho=\rho(c)$,
then we get the Riemann invariants in the form

$r_{\pm}=\frac{u}{2}\pm\sigma(c),\qquad\sigma(c)=\frac{1}{2}\int_{0}^{c}\frac{c}{\rho(c)}\frac{d\rho}{dc}dc,$
(14)

so that $c$ can be used as a wave variable instead of $\rho$. If Eqs. (12) are
solved, then the physical variables $u,\rho,c$ can be expressed in terms of
the Riemann invariants,

$u=r_{+}+r_{-},\quad\rho=\rho(r_{+},r_{-}),\quad c=c(r_{+},r_{-}).$ (15)

In case of unidirectional propagation of a background pulse we get a simple
wave solution with one of the Riemann invariants constant. Let such a pulse
propagate to the right through the uniform quiescent medium with constant
values of $\rho=\rho_{R}$, $c=c_{R}$, $u=u_{R}=0$. Then along this pulse the
variables $u$ and $c$ are related by the formula
$u/2-\sigma(c)=-\sigma(c_{R})$, so that $c=c(u)$ and $r_{+}=u+\sigma(c_{R})$,
$v_{+}=u+c(u)\equiv V_{0}(u)$. Consequently, the second equation (12) is
fulfilled identically and the first one reduces to the Hopf equation

$\frac{\partial u}{\partial t}+V_{0}(u)\frac{\partial u}{\partial x}=0.$ (16)

We are interested in propagation of the small-amplitude edge of a DSW and we
identify this propagation with motion of a point particle with coordinate
$x(t)$ and momentum $k(t)$ which satisfy the Hamilton equations

$\frac{dx}{dt}=\frac{\partial\omega}{\partial
k},\qquad\frac{dk}{dt}=-\frac{\partial\omega}{\partial x}.$ (17)

In usual situations the Hamiltonian (dispersion relation) $\omega$ is a
function of $k,x,t$, $\omega=\omega(k,x,t)$. Following the excellent textbook
gant-66 , let us define in the three-dimensional extended phase space
$(x,k,t)$ an arbitrary closed curve

$C_{0}=\left\\{x=x_{0}(\eta),\,\,k=k_{0}(\eta),\,\,t=t_{0}(\eta)\,\,|\,\,0\leq\eta\leq
1\right\\},$ (18)

where the points corresponding to $\eta=0$ and $\eta=1$ coincide with each
other. Then we can calculate along this curve the integral

$I_{0}=\oint_{C_{0}}(k\delta x-\omega\delta t),$ (19)

where $\delta x$ and $\delta t$ denote differentials of the functions
$x=x_{0}(\eta)$, $t=t_{0}(\eta)$ at the points of the contour $C_{0}$. Now, we
define the flow in the phase space by the equations

$\frac{dx}{dt}=Q(k,x,t),\qquad\frac{dk}{dt}=P(k,x,t)$ (20)

and, when solutions of these equations cross the contour $C_{0}$, we obtain a
tube of trajectories

$x=x(t,\eta),\qquad k=k(t,\eta),$ (21)

where $x(t_{0}(\eta),\eta)=x_{0}(\eta)$, $k(t_{0}(\eta),\eta)=k_{0}(\eta)$. If
we take now another function $t=t_{1}(\eta)$, then we get a new contour

$\begin{split}&C_{1}=\\{x=x_{1}(\eta)=x(t_{1}(\eta),\eta),\\\
&k=k_{1}(\eta)=k(t_{1}(\eta),\eta),\,\,t=t_{1}(\eta)\,\,|\,\,0\leq\eta\leq
1\\}\end{split}$ (22)

and a new value of the integral similar to Eq. (19),

$I_{1}=\oint_{C_{1}}(k\delta x-\omega\delta t),$ (23)

which is taken over the contour $C_{1}$ around the tube of trajectories
generated by the flow (21). Now one can ask, for which Eqs. (20) and
corresponding flows (21) the integrals (19) and (23) are equal to each other
independently of the choice of the contour $C_{1}$ (that is the choice of the
function $t=t_{1}(\mu)$). As was shown in Refs. poincare ; cartan ; gant-66 ,
this is true for

$Q=\frac{\partial\omega}{\partial k},\qquad P=-\frac{\partial\omega}{\partial
x},$ (24)

that is when the flow obeys the Hamilton equations (17). Thus, the invariance
of the Poincaré-Cartan integral

$I=\oint_{C}(k\delta x-\omega\delta t)$ (25)

means that the flow (20) is Hamiltonian.

Now we change the perspective and look at invariance of the integral (25) from
a different point of view. As was indicated in Introduction, we are interested
in invariance of this integral with respect to flows generated by the
hydrodynamic equations (10) or (16). In this case the Hamiltonian $\omega$
depends on $x$ and $t$ only via solutions $\rho=\rho(x,t)$, $u=u(x,t)$ of
these equations, i.e.,

$\omega=\omega(k,\rho,u)$ (26)

and the flow is generated by the equation

$\frac{dx}{dt}=u(x,t).$ (27)

Then invariance of the integral (25) imposes certain conditions on the
dependence of the momentum (wave number) $k$ on the background variables,

$k=k(\rho,u).$ (28)

If the condition of invariance is fulfilled, then the action (4) can be
reduced to the generalized Bohr-Sommerfeld quantization rule (9) which
provides important information about parameters of asymptotic solitons. We
shall consider the condition of invariance of the integral (25) separately for
the simple wave flow (16) and for the general solutions of Eqs. (10),

###  II.1 Simple wave background pulse

In simple wave case, the small-amplitude edge propagates along the background
pulse described by a single variable $u=u(x,t)$ which evolution obeys the Hopf
equation (16). Then we have $\omega=\omega(k,u)$ and we look for such a
function $k=k(u)$, that the Poincaré-Cartan integral (25) remains the same for
any contour $C$ around the tube generated by the flow $dx/dt=V_{0}(u)$.
Following Ref. gant-66 , we introduce a coordinate $\mu$ along the paths that
form the tube and assume that this coordinate is related with the flow by the
relations

$\frac{dx}{V_{0}(u)}=dt=\chi d\mu,$ (29)

where $\chi=\chi(t,x,k)$ is an arbitrary function: its choice fixes the choice
of the coordinate $\mu$ along the flow trajectories. Thus, for any fixed value
of $\mu=\mathrm{const}$ we get a point on every trajectory, so these points
define the closed contour around the tube. As a result, the integral (25),
generally speaking, becomes a function of $\mu$, but we want to find such a
function $k=k(u)$, that this integral does not depend on $\mu$ for any choice
of $\chi$. Differentiation of Eq. (25) with respect to $\mu$ gives

$\begin{split}dI=&\oint\Bigg{\\{}\left[\frac{dk}{du}\delta
x-\left(\frac{\partial\omega}{\partial
k}\frac{dk}{du}+\frac{\partial\omega}{\partial u}\right)\delta t\right]\\\
&\times\left(\frac{\partial u}{\partial t}dt+\frac{\partial u}{\partial
x}dx\right)-\delta x\cdot dx+\delta\omega\cdot dt\Bigg{\\}}=0,\end{split}$

where we have integrated the terms $kd\delta x-\omega d\delta t$ by parts over
the closed contour. Now we substitute $dx=V_{0}\chi d\mu$, $dt=\chi d\mu$ and
obtain after simple transformations with account of Eq. (16) the expression

$\Bigg{\\{}\oint\left[\left(\frac{\partial\omega}{\partial
k}-V_{0}\right)\frac{dk}{du}+\frac{\partial\omega}{\partial
u}\right]\left(\frac{\partial u}{\partial t}\delta t+\frac{\partial
u}{\partial x}\delta x\right)\chi\Bigg{\\}}\cdot d\mu=0.$

This integral must vanish for any choice of $\chi$, so we get the equation

$\frac{dk}{du}=\frac{\partial\omega/\partial
u}{V_{0}(u)-\partial\omega/\partial k}.$ (30)

This equation was obtained earlier from analysis of Whitham equations at the
small-amplitude edge of DSWs el-05  and from Hamilton equations for
propagation of wave packets along a simple wave background pulse kamch-20a ;
ks-21 . Presented here derivation gives important new information that the
background evolution preserves the value of the Poincaré-Cartan integral (25).

Let us assume now that we have found such a solution

$k=k(u,\overline{q}),$ (31)

($\overline{q}$ being an integration constant) of Eq. (30) that the closed
contour $C$ corresponds to the path of the small-amplitude edge gone first in
one direction with negative $k$ and then in the opposite direction with
positive $k$ according to the Hamilton equations (17). Then we can define
“number of oscillations”

$N(\overline{q})=\frac{1}{4\pi}\oint_{C(\overline{q})}(k\delta x-\omega\delta
t)$ (32)

in the DSW corresponding to this value $\overline{q}$ (additional factor $1/2$
compared with Eq. (4) is introduced for taking into account that in the case
of a closed contour the edge goes its path twice). This number
$N(\overline{q})$ remains the same, when we transform the contour to the
initial state at $t=0$ ($\delta t=0$):

$N(\overline{q})=\frac{1}{4\pi}\oint_{C_{0}(\overline{q})}k\delta
x=\frac{1}{2\pi}\int_{x_{1}(\overline{q})}^{x_{2}(\overline{q})}k(u_{0}(x),\overline{q})dx,$
(33)

where $x_{1,2}(\overline{q})$ are two ‘turning point’ at which

$k(u_{0}(x_{1,2}(\overline{q})),\overline{q})=0$ (34)

and $u_{0}(x)$ is the initial distribution of the background pulse. The
asymptotic expression (8) for the number of solitons assumes that these
turning point go to infinities, but less formal consideration makes it clear
at once that there exists the maximal integer value of $N$ for which the
expression (33) with a given distribution $u_{0}(x)$ still has the
corresponding ‘eigenvalue’ $\overline{q}_{N}$ but for $N+1$ such an eigenvalue
does not exist anymore. In this case the integration interval
$(x_{1}(\overline{q}_{N}),x_{2}(\overline{q}_{N}))$ is so wide for intensive
initial pulses with $u_{0}(x)\to 0$ as $|x|\to\infty$, that it can be replaced
with good enough accuracy by the interval $(-\infty,\infty)$. We can say that
$\overline{q}_{N}$ corresponds to the $N$-th soliton in the asymptotic state.
If we decrease $N$ by one, then we get $N(\overline{q}_{N-1})=N-1$, that is
$\overline{q}_{N-1}$ corresponds to $(N-1)$-th soliton, and so on. As a
result, we arrange correspondence between the $n$-th soliton and the parameter
$\overline{q}_{n}$ which is to be determined from the Bohr-Sommerfeld
quantization rule

$\int_{x_{1}(\overline{q}_{n})}^{x_{2}(\overline{q}_{n})}k(u_{0}(x),\overline{q}_{n})dx=2\pi
n,\qquad n=1,2,\ldots,N,$ (35)

where $k=k(u,\overline{q})$ is the solution (31) of Eq. (30). At last, we
notice that Eq. (35) can be interpreted as a quasiclassical quantization rule
for eigenvalues of a linear wave equation for waves propagating along the
“potential” $u_{0}(x)$. As is well known from quantum mechanics (see, e.g.,
LL-3 ), a more accurate description of wave solutions in vicinities of the
turning points leads to the replacement $n\mapsto n-1/2$, and this is very
general phenomenon of the short wavelength asymptotic behavior at $k\to\infty$
(see, e.g., Appendix 11 in Ref. arnold-89 ). We will make such a replacement
in practical applications of this theory without derivation just as a
heuristic approximation confirmed by comparison with numerical results.

We will relate the parameters $\overline{q}_{n}$ obtained from the Bohr-
Sommerfeld rule (35) with the physical parameters of the asymptotic solitons
in the next Section and now we turn to generalization of this quantization
rule to not simple-wave background flows.

###  II.2 General flow of background pulse

Now we suppose that the dispersionless flow is described by the solution
$\rho=\rho(x,t)$, $u=u(x,t)$ of Eqs. (10) with some initial distributions
$\rho=\rho_{0}(x),u=u_{0}(x)$. The dispersion relation in the Hamilton
equations (17) for the packet’s propagation depends on the values $\rho$ and
$u$ at the point $x$ of the packet’s instant location at the moment $t$, so
the function

$\omega=\omega(k,\rho,u)$ (36)

is known from the linearized equations. We look again for such a function
$k=k(\rho,u)$ that the Poincaré-Cartan integral (25) is preserved by the flow
in the extended phase space $(x,k,t)$ provided the contour $C$ is taken around
a tube of streamlines defined by the equation $dx/dt=u(x,t)$. As in the simple
wave case of the background flow, we define a coordinate $\mu$ along
streamlines by the equations

$\frac{dx}{u}=dt=\chi d\mu,$ (37)

where $\chi=\chi(x,t,k)$ is an arbitrary function, and we demand the integral
does not depend neither on $\chi$ (that is the choice of the contour’s form),
no on $\mu$ (that is on position of this contour on the tube). Differentiation
of Eq. (25) with respect to $\mu$ gives with account of Eqs. (10) and (17) the
expression

$\begin{split}dI=\oint&\Bigg{\\{}\left[\frac{\partial
k}{\partial\rho}\left(\frac{\partial\omega}{\partial
k}-u\right)-\frac{c^{2}}{\rho}\frac{\partial k}{\partial
u}+\frac{\partial\omega}{\partial\rho}\right]\rho_{x}\\\
&+\left[\frac{\partial k}{\partial u}\left(\frac{\partial\omega}{\partial
k}-u\right)-{\rho}\frac{\partial
k}{\partial\rho}+\frac{\partial\omega}{\partial u}\right]u_{x}\Bigg{\\}}\\\
&\times(\delta x-u\delta t)\chi\cdot d\mu=0.\end{split}$

Since $\chi$ is an arbitrary function, the expression in curly brackets must
vanish. Moreover, the expressions in square brackets are only functions of
$\rho$ and $u$, whereas the local values of $\rho_{x}$ and $u_{x}$ depend on
the choice of the initial distributions, so they can be considered as
arbitrary functions, too. Hence, the expressions in square brackets vanish
separately and we arrive at the equations

$\begin{split}&\frac{\partial
k}{\partial\rho}=\frac{(v_{g}-u)\frac{\partial\omega}{\partial\rho}+\frac{c^{2}}{\rho}\frac{\partial\omega}{\partial
u}}{c^{2}-(v_{g}-u)^{2}},\\\ &\frac{\partial k}{\partial
u}=\frac{(v_{g}-u)\frac{\partial\omega}{\partial
u}+\rho\frac{\partial\omega}{\partial\rho}}{c^{2}-(v_{g}-u)^{2}},\end{split}$
(38)

which should determine the function $k=k(\rho,u)$ for given dispersion
relation (36) and equation of state $c=c(\rho)$. Eqs. (38) have recently been
derived in Ref. sk-23  from the Hamilton equations (17) under the same
supposition that $k$ is only a function of $\rho$ and $u$.

For existence of the function $k=k(\rho,u)$, its derivatives defined by Eqs.
(38) must satisfy the compatibility condition

$\frac{\partial}{\partial\rho}\left(\frac{\partial k}{\partial
u}\right)=\frac{\partial}{\partial u}\left(\frac{\partial
k}{\partial\rho}\right).$ (39)

If this condition is not fulfilled identically, we can confine ourselves to
the limit of large values of $k$ since just this limit corresponds to the high
quantum numbers in the Bohr-Sommerfeld quantization rule. To this end, we
expand the right-hand sides of Eqs. (38) with respect to a small parameter
$\sim c/k$ and obtain

$\begin{split}\frac{\partial
k}{\partial\rho}=R_{1}(k,\rho,u),\qquad\frac{\partial k}{\partial
u}=R_{2}(k,\rho,u),\end{split}$ (40)

where in $R_{1,2}$ only the terms are held that satisfy the condition (39).
Then these equations allow one to restore the function

$k=k(\rho,u,q),$ (41)

where $q$ is an integration constant determined by the initial value
$k=k_{0}\gg c$. It is convenient to define $q$ in such a way that the
condition $k\sim k_{0}\gg c$ corresponds to $q\gg c$.

The solution (41) (exact of approximate) allows one to define the Poincaré-
Cartan integral

$N(q)=\frac{1}{4\pi}\oint_{C}(k\delta x-\omega\delta t)$ (42)

for contours $C$ around a tube of streamlines of the dispersionless flow. If
we transform such a contour to the initial state at $t=0$, we obtain

$N(q)=\frac{1}{4\pi}\oint_{C_{0}(q)}k\delta
x=\frac{1}{2\pi}\int_{x_{1}(q)}^{x_{2}(q)}k[\rho_{0}(x),u_{0}(x),q]dx,$ (43)

where $x_{1,2}(q)$ are the turning points defined by the equation

$k[\rho_{0}(x_{1,2}(q),u_{0}(x_{1,2}(q)),q]=0.$ (44)

If we only know the asymptotic solution $k\gg c$ of Eq. (40), then Eq. (44)
give only approximate solution for the turning points. Nevertheless, Eq. (43)
gives accurate enough value of the integral, since along the most part of the
integration interval we have $k\gg c$ and vicinities of the turning points
give negligibly small contribution to the integral. As in the simple-wave
background pulse case, we arrange correspondence between the soliton’s number
$n$ and the value $q_{n}$ determined by the Bohr-Sommerfeld rule

$\int_{x_{1}(q_{n})}^{x_{2}(q_{n})}k[\rho_{0}(x),u_{0}(x),q_{n})dx=2\pi
n,\quad n=1,2,\ldots,N.$ (45)

In practical applications of this formula we will make a replacement $n\mapsto
n+1/2$, as it was explained below Eq. (35).

Now we have to relate the constants $\overline{q}_{n}$ (simple-wave case) or
$q_{n}$ (general case) with parameters of asymptotic solitons emerging
eventually from a given initial pulse.

##  III Parameters of asymptotic solitons

Now, when each soliton is labeled by a specific value of the parameter
$\overline{q}_{n}$ or $q_{n}$, we need to relate it with soliton’s physical
parameters, as, for example, its velocity $V_{n}$. This can be done with the
use of Stokes remark lamb ; stokes  (see also ai-77 ; dkn-03 ) that the small
amplitude tails $\propto\exp[\pm\widetilde{k}(x-Vt)]$ of a soliton obey the
same linearized equations as a small-amplitude harmonic wave
$\propto\exp[i(kx-\omega t)]$. Consequently, the soliton’s velocity $V_{n}$
can be expressed in the form

$V_{n}=\widetilde{\omega}(\widetilde{k}_{n})/\widetilde{k}_{n},\quad\text{where}\quad\widetilde{\omega}(\widetilde{k})=i\omega(-i\widetilde{k}),$
(46)

and $\widetilde{k}$ denotes the soliton’s inverse half-width (see Refs. el-05
; kamch-20a ; kamch-21 ). This means analytical continuation of the function
$\omega=F(k^{2})$ from the interval of positive values of its argument $k^{2}$
to the interval of its negative values of the argument
$k^{2}=-\widetilde{k}^{2}<0$. In practice it means that we consider the same
function $F$ at different intervals of its argument $k^{2}$.

We generalize this Stokes observation to other functions of $k^{2}$ that can
also be continued to the region $k^{2}<0$ where they express some assertions
about the inverse half-widths $\widetilde{k}$ of solitons. For example,
suppose we have found the exact solution (31) of Eq. (30) and it can be
expressed in the form

$k^{2}=\overline{K}(u,\overline{q}).$ (47)

This formula can be continued to the solitonic region to give

$\widetilde{k}^{2}=-\overline{K}(u,\overline{q}).$ (48)

Suppose we have found from the Bohr-Sommerfeld rule the value
$\overline{q}_{n}$ for the $n$-th soliton and at asymptotically large time
this soliton propagates along zero background $u=0$. Then its inverse half-
width is given by the formula

$\widetilde{k}_{n}^{2}=-\overline{K}(0,\overline{q}_{n})$ (49)

and, consequently, its velocity is readily obtained from Eq. (46). This
approach reproduces in a different form the theory developed for simple wave
pulses in Refs. egs-08 ; mfweh-20 .

We can extend the above approach to the case of exact solutions of Eqs. (38),
so we imply that the condition (39) is fulfilled. Again we express the exact
solution in the form

$k^{2}=K(\rho,u,q)$ (50)

and analytically continue it to the soliton region to obtain

$\widetilde{k}^{2}=-K(\rho,u,q).$ (51)

If we found for the $n$-th soliton the eigenvalue $q_{n}$ from the Bohr-
Sommerfeld rule (45) and if the asymptotic soliton propagates along the
uniform background with $\rho=\rho_{R},u=0$, then its inverse half-width is
given by the formula

$\widetilde{k}_{n}^{2}=-K(\rho_{R},0,q_{n})$ (52)

and its velocity is to be found from the Stokes rule (46).

However such an analytical continuation becomes impossible in case of
asymptotic solutions (41) which are only correct for $k\gg c\sim u$.
Nevertheless, we can overcome this difficulty in the following way. Suppose we
have found the asymptotic solution in the form (50). Naturally, this general
solution remains correct in case of simple-wave pulses when there is the
functional relationship between $\rho$ and $u$, say,
$\rho=\rho_{\mathrm{sw}}(u)$, so we get the asymptotic formula

$k^{2}=K(\rho_{\mathrm{sw}}(u),u,q),\qquad k\gg u.$ (53)

Now, the asymptotic soliton trains propagate along simple-wave solutions, so
their inverse half-widths are asymptotic specifications of the formula (48)
obtained by analytic continuation of Eq. (47). If we take its series expansion
with respect to small parameter $u$, then the resulting formula

$k^{2}=\overline{K}_{\text{asymp}}(u,\overline{q})$ (54)

must coincide with Eq. (53). This matching condition allows one to find the
relationship between the parameters $q$ and $\overline{q}$:

$\overline{q}=\overline{q}(q).$ (55)

Consequently, the inverse half-widths of asymptotic solitons are given by Eq.
(49),

$\widetilde{k}_{n}^{2}=-\overline{K}(0,\overline{q}(q_{n})),$ (56)

where the parameters $q_{n}$ are to be obtained from the Bohr-Sommerfeld rule
(45) and we assumed that the asymptotic solitons propagate along the
background with $u=0$. At last, the soliton’s velocities are to be found
according to the Stokes rule (46).

Let us illustrate this theory by examples.

##  IV Solitons evolved from a large deep in the generalized NLS equation
theory

Here we shall consider evolution of a pulse according to the defocusing
generalized NLS equation (generalized Gross-Pitaevskii equation for a Bose-
Einstein condensate of atoms with repulsive interaction)

$i\psi_{t}+\frac{1}{2}\psi_{xx}-f(|\psi|^{2})\psi=0,$ (57)

where the nonlinearity function is positive: $f(\rho)>0$, $f(0)=0$. It is well
known that by means of the substitution

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

this equation can be transformed to the hydrodynamic-like system

$\begin{split}&\rho_{t}+(\rho u)_{x}=0,\\\
&u_{t}+uu_{x}+\frac{c^{2}}{\rho}\rho_{x}+\bigg{(}\frac{\rho_{x}^{2}}{8\rho^{2}}-\frac{\rho_{xx}}{4\rho}\bigg{)}_{x}=0,\end{split}$
(59)

where

$c^{2}=\rho f^{\prime}(\rho).$ (60)

Linearization of this system with respect to small deviations from a uniform
flow yields the Bogoliubov dispersion relation

$\omega=k\left(u\pm\sqrt{c^{2}+\frac{k^{2}}{4}}\right),$ (61)

that is $c$ has the meaning of the sound velocity in the dispersionless limit
$k\to 0$. Equations (59) reduce in the same limit to the standard Euler
equations (10).

Now we can check whether the derivatives (38) commute. In this case it is
convenient to replace $\rho$ by the variable $c$ where the function
$\rho=\rho(c)$ is obtained by inversion of the function $c=c(\rho)$ (see Eq.
(60)) so Eqs. (38) are cast to the form

$\begin{split}&\frac{\partial k}{\partial
c}=-\frac{c[(2+c\rho^{\prime}/\rho)k^{2}+4(1+c\rho^{\prime}/\rho)c^{2}]}{k(k^{2}+3c^{2})},\\\
&\frac{\partial k}{\partial
u}=-\frac{\sqrt{k^{2}+4c^{2}}[k^{2}+2(1+\rho/(c\rho^{\prime}))c^{2}]}{k(k^{2}+3c^{2})}.\end{split}$
(62)

The straightforward calculation yields (see also sk-23 )

$\begin{split}&\frac{\partial}{\partial c}\left(\frac{\partial k}{\partial
u}\right)-\frac{\partial}{\partial u}\left(\frac{\partial k}{\partial
c}\right)=\frac{\sqrt{k^{2}+4c^{2}}}{k\rho\rho^{\prime}(k^{2}+3c^{2})^{2}}\\\
&\times\left[(k^{2}+6c^{2})\rho^{\prime
2}(\rho^{\prime}-2c)+2(k^{2}+3c^{2})\rho^{2}(c\rho^{\prime\prime}-\rho^{\prime})\right],\end{split}$
(63)

and we see that the compatibility condition is only fulfilled for the case of
the Kerr-like nonlinearity with $f(\rho)=\rho$ and $\rho=c^{2}$. Otherwise, it
is fulfilled only in the limit $k\to\infty$, so it is natural to consider
these two situations separately.

###  IV.1 Kerr-like nonlinearity $f(\rho)=\rho$

NLS equation with Kerr-like nonlinearity is completely integrable zs-73 , so
the Bohr-Sommerfeld quantization rule for solitons parameters can be derived
from the Zakharov-Shabat linear spectral problem JLML-99 ; kku-02 . It is
instructive to see, how these known results follow from our general approach
not based on the complete integrability property. We assume here that the
initial pulse is represented by a deep in the density or sound velocity
distribution $c=c_{0}(x)$ and there is also some distribution of the initial
flow velocity $u=u_{0}(x)$. Far enough from the initial pulse we have
$c_{0}(x)\to c_{R}$, $u_{0}(x)\to 0$ as $|x|\to\infty$. Naturally, these two
distributions are equivalent to some initial distributions $r_{\pm}^{(0)}(x)$
of the dispersionless Riemann invariants (see (14))

$r_{\pm}=\frac{u}{2}\pm c.$ (64)

In this case with $\rho=c^{2}$ Eqs. (62) can be written in the form

$\begin{split}\frac{\partial k}{\partial c}=-4\frac{c}{k},\qquad\frac{\partial
k}{\partial u}=-\frac{\sqrt{k^{2}+4c^{2}}}{k},\end{split}$ (65)

and they have the exact solution

$k^{2}=(q-u)^{2}-4c^{2},$ (66)

where $q$ is an integration constant. If we write it in the form

$k^{2}=4\left(\frac{q}{2}-\frac{u}{2}-c\right)\left(\frac{q}{2}-\frac{u}{2}+c\right)=4(\lambda-
r_{+})(\lambda-r_{-}),$ (67)

where we have defined $\lambda=q/2$, then the Bohr-Sommerfeld rule (45) reads

$\begin{split}\int_{x_{1}(\lambda_{n})}^{x_{2}(\lambda_{n})}\sqrt{(\lambda_{n}-r_{+}^{(0)}(x))(\lambda_{n}-r_{-}^{(0)}(x))}\,dx=n\pi,\\\
n=1,2,\dots,N.\end{split}$ (68)

It coincides in the main approximation in the limit $n\to\infty$ with the
asymptotic formula for the linear spectral problem eigenvalues (see JLML-99 ;
kku-02 ). It is worth, however, notice that in better approximation developed
in Ref. JLML-99  the integer $n$ should be replaced by $n+1/2$ with
$n=1,2,\ldots,N$ in agreement with the replacement $n\mapsto n+1/2$ mentioned
below Eq. (35).

Since Eq. (66) represents the exact solution of Eqs. (65), it is correct for
any values of $k$ and can be continued to the soliton region, so we obtain

$\widetilde{k}^{2}=-4(\lambda-r_{+})(\lambda-r_{-}).$ (69)

In the asymptotic region, when solitons propagate along background with
$r_{\pm}=\pm c_{R}$, we get for the half-width of the $n$-th soliton the
expression $\widetilde{k}_{n}^{2}=4(c_{R}^{2}-\lambda_{n}^{2})$. Consequently,
according to the Stokes rule (46), its velocity is equal to

$V_{n}=\frac{\widetilde{\omega}(\widetilde{k}_{n})}{\widetilde{k}_{n}}=\sqrt{c_{R}^{2}-\frac{\widetilde{k}^{2}}{4}}=\lambda_{n}$
(70)

in agreement with the known results (see, e.g., Ref. kku-02 ).

The possibility to find the exact solution of Eqs. (62) in this case is
related, apparently, with complete integrability of the NLS equation (57) with
Kerr nonlinearity $f(\rho)=\rho$. We will clarify this point in Section V.

###  IV.2 Non-Kerr nonlinearity

In the limit of large $k\gg c$ the right-hand side of Eq. (63) tends to zero
as $\propto k^{-2}$, that is the condition (39) is fulfilled in this limit and
we can find the asymptotic expression for the function $k=k(c,u)$. Eqs. (62)
in the limit of large $k$ can be written as

$\frac{\partial k^{2}}{\partial
c^{2}}=-2\left(1+\frac{c^{2}}{\rho}\frac{d\rho}{dc^{2}}\right),\qquad\frac{\partial
k}{\partial u}=-1.$ (71)

The second equation gives $k=q-u+F(c^{2})$, where in the main approximation
$k\approx q\sim k_{0}\gg|u|,c$. We suppose that $F(c^{2})\sim c^{2}$, so
$k^{2}\approx(q-u)^{2}+2qF(c^{2})$, where we have neglected small terms
$uF\sim c^{3}$ and $F^{2}\sim c^{4}$. Then the first equation (71) gives

$F(c^{2})=-\frac{1}{q}\left(c^{2}+\int_{0}^{\rho(c)}\frac{c^{2}}{\rho}d\rho\right)\sim\frac{c^{2}}{q}$

in agreement with our supposition about the order of magnitude of $F$. Thus,
we obtain the asymptotic solution

$k^{2}=(q-u)^{2}-2\left(c^{2}+\int_{0}^{\rho(c)}\frac{c^{2}}{\rho}d\rho\right).$
(72)

For a particular case of the nonlinearity function

$f(\rho)=\frac{1}{\gamma-1}\rho^{\gamma-1},\quad\text{when}\quad
c^{2}=\rho^{\gamma-1},\quad\rho=c^{2/(\gamma-1)},$ (73)

we get

$k^{2}=(q-u)^{2}-\frac{2\gamma}{\gamma-1}c^{2}.$ (74)

If $\gamma=2$, that is $f(\rho)=\rho$, this asymptotic formula reproduces the
exact solution (66). The function $k=k(c,u)$ defined by Eq. (74) can be used
in the Bohr-Sommerfeld rule (45), so we get

$\begin{split}&\int_{x_{1}(q_{n})}^{x_{2}(q_{n})}\sqrt{(q_{n}-u_{0}(x))^{2}-\frac{2\gamma}{\gamma-1}c_{0}^{2}(x)}\,\,dx\\\
&=2\pi\left(n+\frac{1}{2}\right),\qquad n=1,2,\ldots,N,\end{split}$ (75)

where we have made the replacement $n\mapsto n+1/2$. As a result, we obtain a
set of parameters $q_{n}$ for solitons emerging from the pulse with given
initial distributions $c=c_{0}(x),u=u_{0}(x)$. To relate $q_{n}$ with
velocities of these solitons, we should turn to the simple-wave solution of
the corresponding equation (30).

The problem of evolution of a step-like discontinuity of the simple-wave type
was studied in much detail in Ref. hoefer-14  by the method of Ref. el-05 .
The solution of Eq. (30) found there can be written in the form

$k^{2}=4c^{2}(\alpha^{2}(c,\overline{q})-1),$ (76)

where the function $\alpha=\alpha(c,\overline{q})$ is defined in implicit form
in our notation by the formula

$\frac{c}{\overline{q}}=\left(\frac{2}{1+\alpha}\right)^{\frac{\gamma-1}{3\gamma-5}}\left(\frac{\gamma+1}{3-\gamma+2(\gamma-1)\alpha}\right)^{\frac{2(\gamma-2)}{3\gamma-5}}.$
(77)

In the limit $c\ll\overline{q}$ we have $\alpha\gg 1$, so we obtain the series
expansion

$\frac{c}{\overline{q}}=\beta(\gamma)\left(\frac{1}{\alpha}-\frac{1}{\gamma-1}\cdot\frac{1}{\alpha^{2}}+\frac{\gamma^{2}-3\gamma+6}{4(\gamma-1)^{2}}\cdot\frac{1}{\alpha^{3}}+\ldots\right),$
(78)

where

$\beta(\gamma)=2^{\frac{3-\gamma}{3\gamma-5}}\left(\frac{\gamma+1}{\gamma-1}\right)^{\frac{2(\gamma-2)}{3\gamma-5}}.$
(79)

Inversion of this series and substitution of the result into Eq. (76) yield
with necessary accuracy

$\begin{split}k^{2}\approx&(2\overline{q}\beta)^{2}-\frac{4}{\gamma-1}(2\overline{q}\beta)c\\\
&+\left(\frac{4}{(\gamma-1)^{2}}-\frac{2}{\gamma-1}-2\right)c^{2},\qquad
c\ll\overline{q}.\end{split}$ (80)

This asymptotic expression should be matched with Eq. (74) for the simple wave
pulses.

The dispersionless Riemann invariants (14) in our case (73) are equal to

$r_{\pm}=\frac{u}{2}\pm\frac{c}{\gamma-1}.$ (81)

We consider the asymptotic train of right-propagating solitons along the
background with

$r_{-}=\frac{u}{2}-\frac{c}{\gamma-1}=-\frac{c_{R}}{\gamma-1}=\mathrm{const},$

where $c\to c_{R}$, $u\to 0$ as $|x|\to\infty$. Hence, we have to substitute

$u=\frac{2}{\gamma-1}(c-c_{R})$ (82)

into Eq. (74) to obtain

$\begin{split}k^{2}=&\left(q+\frac{2c_{R}}{\gamma-1}\right)^{2}-\frac{4}{\gamma-1}\left(q+\frac{2c_{R}}{\gamma-1}\right)c\\\
&+\left(\frac{4}{(\gamma-1)^{2}}-\frac{2}{\gamma-1}-2\right)c^{2},\qquad c\ll
q.\end{split}$ (83)

This expression coincides with Eq. (80) for

$\overline{q}=\frac{1}{2\beta}\left(q+\frac{2c_{R}}{\gamma-1}\right).$ (84)

In the solitonic region Eq. (76) gives

$\widetilde{k}^{2}=4c^{2}(1-\alpha^{2}(c,\overline{q})),$ (85)

so the asymptotic velocities of solitons propagating along a uniform
background with $c=c_{R}$, $u=0$ are equal to

$\begin{split}V_{n}&=\frac{\widetilde{\omega}(\widetilde{k}_{n})}{\widetilde{k}_{n}}=\sqrt{c_{R}^{2}-\frac{\widetilde{k}_{n}^{2}}{4}}\\\
&=c_{R}\alpha\left(c_{R},\frac{1}{2\beta}\left(q_{n}+\frac{2c_{R}}{\gamma-1}\right)\right),\end{split}$
(86)

where $q_{n}$ are to be obtained from the Bohr-Sommerfeld quantization rule
(75).

Soliton’s velocity $V$ is related with minimal value $\rho_{m}$ of the density
at its center by the formula (see Ref. ks-09 )

$V^{2}=\frac{2\rho_{m}}{(\rho_{R}-\rho_{m})^{2}}\int_{\rho_{m}}^{\rho_{R}}\left[f(\rho_{R})-f(\rho)\right]d\rho,$
(87)

so when the soliton’s velocity is known, we can easily find the amplitude of
this soliton.

Let us illustrate the theory by a concrete example $\gamma=3$, when
$\beta=\sqrt{2}$ and the function $\alpha=\alpha(c)$ can be found in the
explicit form (see Ref. hoefer-14 ),

$\alpha(c,\overline{q})=\frac{1}{2}\left(\sqrt{1+8\left(\frac{\overline{q}}{c}\right)^{2}}-1\right).$
(88)

As a result, we obtain a simple formula

$V_{n}=\frac{1}{2}\left(\sqrt{c_{R}^{2}+(c_{R}+q_{n})^{2}}-c_{R}\right).$ (89)

Figure 1: Distribution of $c(x)$ in the right-propagating soliton train
evolved from the pulse with the initial distributions (90) at $t=1000$.

We compared these analytical predictions with results of the exact numerical
solution of the gNLS equation (57) with $f(\rho)=\rho^{2}/2$ ($\gamma=3$) for
the initial distributions

$c_{0}(x)=1-\frac{0.9}{\cosh^{2}(x/20)},\qquad u_{0}(x)=0.$ (90)

In this case the initial pulse evolves into two symmetrical right- and left-
propagating dark soliton trains and a typical distribution of $c(x)$ at
$t=1000$ is shown in Fig. 1 for the right-propagating solitons. If we assume
that all the solitons started their motion at $x=0$ at the moment $t=0$ with
asymptotic velocities $V_{n}$, then their coordinates at the moment $t$ are
equal to $x_{n}(t)=t(\sqrt{1+(1+q_{n})^{2}}-1)/2$. This formula allows one to
find $q_{n}$ from numerical values of the coordinates,

$q_{n}^{\text{num}}=\sqrt{\left(2x_{n}(t)/t+1\right)^{2}-1}-1.$ (91)

These values can be compared with the values $q_{n}^{\text{B.-S.}}$ found from
the Bohr-Sommerfeld rule (75) in a particular case of the initial
distributions (90), and the results are shown in Fig. 2, where crosses
correspond to $q_{n}^{\text{B.-S.}}$ and dots to $q_{n}^{\text{num}}$. As one
can see, the agreement is quite good even for small values of $n$.

Figure 2: Parameters $q_{n}^{\text{B.-S.}}$ obtained from the Bohr-Sommerfeld
rule (crosses) and $q_{n}^{\text{num}}$ obtained from the numerical solution
(red dots).

In case of $f(\rho)=\rho^{2}/2$ Eq. (87) gives the expression for the minimal
density at the $n$-th soliton

$\rho_{m}^{(n)}=\left[\rho_{R}^{2}+\frac{3}{4}\left(\sqrt{\rho_{R}^{2}+(q_{n}+\rho_{R})^{2}}-\rho_{R}\right)^{2}\right]^{1/2}-\rho_{R}.$
(92)

This formula also agrees very well with the minimal values of density at
centers of solitons shown in Fig. 1.

##  V Connection with the theory of completely integrable equations

As we saw in the preceding Section, there is sharp difference between the
cases with $\gamma=2$ and $\gamma\neq 2$: if $\gamma=2$ the derivatives (62)
commute and we get the exact solution (66) of these equations, whereas if
$\gamma\neq 2$ we can only obtain the asymptotic solution (74) correct in the
limit $q\gg|u|,c$. It seems very plausible that such a difference between two
situations is related with the complete integrability of the NLS equation (57)
with $f(\rho)=\rho$, ($\gamma=2$). Here we shall consider this relationship
for the class of completely integrable equations which belong to the Ablowitz-
Kaup-Newell-Segur (AKNS) scheme akns-74 .

Usually the AKNS scheme is formulated in $2\times 2$-matrix form, but for
discussion of its quasiclassical limit it is convenient to use its scalar form
(see Ref. ak-02 ). So we assume that an integrable equation under
consideration can be written as a compatibility condition of two linear
equations

$\phi_{xx}=\mathcal{A}\phi,\qquad\phi_{t}=-\frac{1}{2}\mathcal{B}_{x}\phi+\mathcal{B}\phi_{x}$
(93)

for the function $\phi$, where $\mathcal{A}$ and $\mathcal{B}$ depend on the
wave variables (let them be $\rho$ and $u$) and the spectral parameter
$\lambda$. Then the compatibility condition $(\phi_{xx})_{t}=(\phi_{t})_{xx}$
leads to the equation

$\mathcal{A}_{t}-2\mathcal{B}_{x}\mathcal{A}-\mathcal{B}\mathcal{A}_{x}+\frac{1}{2}\mathcal{B}_{xxx}=0,$
(94)

which must be fulfilled for any values of $\lambda$, so we get the nonlinear
equations under consideration. This provides a number of useful relations for
solutions without actual finding them. In particular, the first Eq. (93) has
two basis solutions $\phi_{+},\phi_{-}$, so we define the function

$g=\phi_{+}\phi_{-},$ (95)

which satisfies the equations

$g_{xxx}-2\mathcal{A}_{x}{g}-4\mathcal{A}{g}_{x}=0,\qquad{g}_{t}=\mathcal{B}{g}_{x}-\mathcal{B}_{x}{g}.$
(96)

The first one is readily integrated to give

$\frac{1}{2}{g}{g}_{xx}-\frac{1}{4}{g}_{x}^{2}-\mathcal{A}{g}^{2}={P},$ (97)

where $P$ is an integration constant. Then the second Eq. (96) yields

$\left(\frac{\sqrt{P}}{{g}}\right)_{t}=\left(\frac{\sqrt{P}}{{g}}\mathcal{B}\right)_{x},$
(98)

where we have introduced the constant factor $\sqrt{P}$ under the
differentiation signs to fix the asymptotic behavior of $g$ in the limit
$\lambda\to\infty$. Eq. (98) can serve as a generating function of
conservation laws of our nonlinear wave equations and averaging of this
generating function provides a convenient method of derivation of the Whitham
modulation equations kamch-94 ; kamch-04 . Solutions of the linear equations
(93) can also be expressed in terms of the function $g$ (see, e.g., kku-02 ;
ak-01 )

$\phi_{\pm}=\sqrt{g}\exp\left(\pm i\int^{x}\frac{\sqrt{P}}{g}\,dx\right).$
(99)

Now we can turn to the discussion of the quasiclassical limit of these
equations.

Let us study transformation of a large-scale pulse to soliton trains. This
transformation can be represented as formation and evolution of dispersive
shock waves (see, e.g., kamch-21a ). According to Gurevich and Pitaevskii
gp-73 , such an evolution can be represented as a slow change of modulation
parameters in the periodic solution of equations under consideration what
leads to the slow variation of $P$ as well as parameters in $g$. Nevertheless,
following to the Whitham method whitham ; Whitham-74  in Krichever’s
formulation krichever-88  (see also Ref. dn-93 ), we assume that expressions
(98) and (99) remain correct during evolution of $\rho$ and $u$ from the
initial smooth distributions to the asymptotic soliton trains. When $\rho$ and
$u$ are still smooth functions, the function $g$ is also smooth, so the terms
with derivatives in Eq. (97) can be neglected and we obtain

$\frac{\sqrt{P}}{g}\approx\sqrt{-\overline{\mathcal{A}}}\equiv\overline{k}(\rho,u,\lambda),$
(100)

where $\overline{\mathcal{A}}=\overline{\mathcal{A}}(\rho,u)$ is to be
obtained from $\mathcal{A}$ in the same approximation with omitted derivatives
of $\rho$ and $u$. Then we arrive at a quasiclassical limit of Eq. (99),

$\phi_{\pm}\approx\sqrt{g}\exp\left(\pm
i\int^{x}\overline{k}(\rho,u,\lambda)\,dx\right).$ (101)

The condition that $\phi_{\pm}$ are single-valued functions of $x$ gives the
Bohr-Sommerfeld quantization rule (see, e.g., karpman-73 ; JLML-99 ; kku-02 )

$\begin{split}\int_{x_{1}(\lambda_{n})}^{x_{2}(\lambda_{n})}\overline{k}(\rho_{0}(x),u_{0}(x),\lambda_{n})\,dx=\pi\left(n+\frac{1}{2}\right),\\\
n=1,2,\dots,\end{split}$ (102)

where ${x_{1}(\lambda_{n})},{x_{2}(\lambda_{n})}$ are the turning points.
Comparison with Eq. (45) gives the expression for the wave number $k$ in our
approach,

$k(\rho,u)=2\overline{k}(\rho,u)=2\sqrt{-\overline{\mathcal{A}}(\rho,u)},$
(103)

as well as some relationship between the integration constant $q$ in our
theory and the spectral parameter $\lambda$ in the AKNS scheme. Consequently,
in case of nonlinear wave equations, which are completely integrable in the
AKNS scheme, we obtain the expression for the function $k=k(\rho,u)$ without
solving Eqs. (38).

Substitution of $\sqrt{P}/g=k/2$ into Eq. (98) gives the equation

$k_{t}-(k\overline{\mathcal{B}})_{x}=0,$ (104)

where $\overline{\mathcal{B}}(\rho,u)$ is obtained from $\mathcal{B}$ in the
same quasiclassical approximation. This equation must coincide with the number
of waves conservation law (2), so we get

$\frac{\omega}{k}=-\overline{\mathcal{B}}.$ (105)

As we see, the quasiclassical limits of the functions $\mathcal{A}$ and
$\mathcal{B}$ are related with the dispersion relation
$\omega=\omega(k,\rho,u)$ of linear waves and the function $k=k(\rho,u,q)$
which can be obtained by solving Eqs. (38) if these derivatives commute. From
this point of view, the condition of commutativity of derivatives in Eqs. (38)
can be considered as an “integrability test” of nonlinear equations in
framework of the AKNS scheme. If this test is fulfilled, then equations

$\begin{split}\overline{\mathcal{A}}(\rho,u,q)&=-\frac{1}{4}k^{2}(\rho,u,q),\\\
\overline{\mathcal{B}}(\rho,u,q)&=-\frac{\omega(k(\rho,u,q),\rho,u)}{k(\rho,u,q)}\end{split}$
(106)

give us the quasiclassical limit of the functions $\mathcal{A}$, $\mathcal{B}$
in the Lax pair (93) and $q$ plays the role of the spectral parameter.

Let us apply the above consideration to the NLS equation (57) with
$f(\rho)=\rho$. In this case we have (see Ref. alber-93 )

$\begin{split}\mathcal{A}&=-\left(\lambda+\frac{i\psi_{x}}{2\psi}\right)^{2}+\psi^{*}\psi-\left(\frac{\psi_{x}}{2\psi}\right)_{x},\\\
\mathcal{B}&=-\lambda+\frac{i}{2}\frac{\psi_{x}}{\psi}.\end{split}$ (107)

The substitution of Eq. (58) and neglecting derivatives of $\rho$ and $u$
yield

$\overline{\mathcal{A}}=-\left(\lambda-\frac{u}{2}\right)^{2}+\rho,\qquad\overline{\mathcal{B}}=-\lambda-\frac{u}{2}.$
(108)

Then Eq. (103) gives $k^{2}=4\left[(\lambda-u/2)^{2}-\rho\right]$ in agreement
with Eq. (66) for $q=2\lambda$. Exclusion of $\lambda$ from this equation and
from Eq. (105), that is from $\omega=k(\lambda+u/2)$, reproduces the
dispersion relation (61).

##  VI Conclusion

The asymptotic method developed in this paper provides a simple way to finding
parameters of solitons evolved from given initial distributions of wave
variables for a wide class of nonintegrable nonlinear equations. Besides that,
it presents a simple method of derivation of the asymptotic Bohr-Sommerfeld
generalized rule for a linear spectral problem in framework of the AKNS
scheme. Condition of commutativity of derivatives defined in Eqs. (38) or
their generalizations (see Ref. sk-23 ) can be used as an ‘integrability
test’: if they commute, then it is quite plausible that the nonlinear wave
equations under consideration are completely integrable and Eqs. (106) give
useful information about quasiclassical limit of functions defining the Lax
pair in the scalar representation of the AKNS scheme. At last, the results of
the paper can be used in discussions of problems about propagation of high-
frequency wave packets along non-uniform and evolving with time background, as
it was demonstrated in Ref. sk-23 .

###### Acknowledgements.

I thank L. F. Calazans de Brito and D. V. Shaykin for useful discussions. This
research is funded by the research project FFUU-2021-0003 of the Institute of
Spectroscopy of the Russian Academy of Sciences and by the Foundation for the
Advancement of Theoretical Physics and Mathematics “BASIS”.

## References

  * (1) V. E. Zakharov, S. V. Manakov, S. P. Novikov, and L. P. Pitaevskii, The Theory of Solitons: The Inverse Scattering Method, (Nauka, Moscow, 1980) (translation: Consultants Bureau, 1984).
  * (2) M. J. Ablowitz, H. Segure, Solitons and the Inverse Scattering Transform, (SIAM, Philadelphia, 1981).
  * (3) A. C. Newell, Solitons in Mathematics and Physics, (SIAM, Philadelphia, 1985).
  * (4) S. C. Gardner, J. M. Greene, M. D. Kruskal, R. M. Miura, Phys. Rev. Lett., 19, 1095 (1967).
  * (5) V. I. Karpman, Phys. Lett. A, 25, 708 (1967).
  * (6) V. I. Karpman, Non-Linear Waves in Dispersive Media, (Nauka, Moscow, 1973) (English translation: Pergamon Press, Oxford, 1975).
  * (7) S. Jin, C. D. Levermore, D. W. McLaughlin, Comm. Pure Appl. Math., 52, 613 (1999).
  * (8) A. M. Kamchatnov, R. A. Kraenkel, B. A. Umarov, Phys. Rev. E 66, 036609 (2002).
  * (9) V. E. Zakharov, A. B. Shabat, Zh. Eksp. Teor. Fiz., 64, 1627 (1973) [Sov. Phys. JETP, 37, 823 (1973)].
  * (10) A. V. Gurevich and L. P. Pitaevskii, Zh. Eksp. Teor. Fiz., 65, 590 (1973) [Sov. Phys.-JETP, 38, 291 (1974)].
  * (11) G. B. Whitham, Proc. Roy. Soc. London, A  283, 238 (1965).
  * (12) G. B. Whitham, Linear and Nonlinear Waves, (Wiley Interscience, New York, 1974).
  * (13) A. V. Gurevich and L. P. Pitaevskii, Zh. Eksp. Teor. Fiz., 93, 871 (1987) [Sov. Phys.-JETP, 93, 871 (1987)].
  * (14) A. M. Kamchatnov, Usp. Fiz. Nauk., 191, 52-87 (2021) [Phys.–Uspekhi, 64, 48-82 (2021)].
  * (15) C. Lanczos, The Variational Principles of Mechanics, (University of Toronto Press, Toronto, 1962).
  * (16) L. D. Landau and E. M. Lifshitz, Fluid Mechanics, (Pergamon, Oxford, 1987).
  * (17) G. A. El, Chaos, 15, 037103 (2005).
  * (18) A. M. Kamchatnov, Chaos 30, 123148 (2020).
  * (19) A. M. Kamchatnov, D. V. Shaykin, Phys. Fluids, 33, 052120 (2021).
  * (20) A. M. Kamchatnov, Phys. Rev. E 99, 012203 (2019).
  * (21) A. M. Kamchatnov, Zh. Eksp. Teor. Fiz. 159, 76 (2021) [JETP, 132, 63 (2021)].
  * (22) L. F. Calazans de Brito, A. M. Kamchatnov, Phys. Rev. E 104, 054203 (2021).
  * (23) G. A. El, A. Gammal, E. G. Khamis, R. A. Kraenkel, A. M. Kamchatnov, Phys. Rev. A 76, 053813 (2007).
  * (24) G. A. El, R. H. J. Grimshaw, N. F. Smyth, Physica D, 237, 2423 (2008).
  * (25) M. D. Maiden, N. A. Franco, E. G. Webb, G. A. El, and M. A. Hoefer, J. Fluid Mech. 883, A10 (2020).
  * (26) H. Poincaré, Les méthodes nouvelles de la Mécanique céleste, t. III, (Paris, Gauthier-Villar, 1899).
  * (27) É. Cartan, Leçons sur les invariants intégraux, (Hermann, Paris, 1922).
  * (28) H. Lamb, Hydrodynamics, (Cambridge University Press, Cambridge, 1994).
  * (29) G. G. Stokes, Mathematical and Physical Papers, Vol. V, p. 163 (Cambridge University Press, Cambridge, 1905).
  * (30) O. Akimoto and K. Ikeda, J. Phys. A: Math. Gen. 10, 425 (1977).
  * (31) S. A. Darmanyan, A. M. Kamchatnov, M. Neviére, Zh. Eksp. Teor. Fiz., 123, 997 (2003) [JETP, 96, 876 (2003)].
  * (32) F. R. Gantmacher, Lectures in Analytical Mechanics, (Moscow, Nauka, 1966) (English translation: Mir, Moscow, 1975).
  * (33) L. D. Landau and E. M. Lifshitz, Quantum Mechanics, Pergamon, Oxford (1981).
  * (34) V. I. Arnold, Mathematical Methods of Classical Mechanics, (N. Y., Springer,1989).
  * (35) D. V. Shaykin, A. M. Kamchatnov, Phys. Fluids (accepted, 2023), preprint arXiv:2303.16592 (2023).
  * (36) M. A. Hoefer, J. Nonlinear Sci., 24, 525 (2014).
  * (37) A. M. Kamchatnov and M. Salerno, J. Phys. B: At. Mol. Opt. Phys. 42, 185303 (2009).
  * (38) M. J. Ablowitz, D. J. Kaup, A. C. Newell, H. Segur, Stud. Appl. Math. 53, 249 (1974).
  * (39) A. M. Kamchatnov and R. A. Kraenkel, J. Phys. A: Math. Gen., 35, L13 (2002).
  * (40) A. M. Kamchatnov, Phys. Lett. A, 186, 387 (1994).
  * (41) A. M. Kamchatnov, Physica D, 188, 247 (2004).
  * (42) A. M. Kamchatnov, J. Phys. A: Math. Gen., 34, L441 (2001).
  * (43) I. M. Krichever, Funk. Analiz. Prilozh. 22, 37 (1988) [Funct. Anal. Appl., 22, 200 (1988)].
  * (44) B. A. Dubrovin, S. P. Novikov, Sov. Sci. Rev. C. Math. Phys. 9, 1 (1993).
  * (45) S. J. Alber, Complex deformations of integrable Hamiltonians over generalized Jacobi varieties, in Nonlinear Processes in Physics, eds A. S. Fokas, D. J. Kaup, A. C. Newell, and V. E. Zakharov, p.6 (Berlin, Springer, 1993).