Nonlinear interactions between solitons and dispersive shocks in focusing media

TL;DR

This paper investigates how solitons interact with dispersive shocks in focusing media using the nonlinear Schrödinger equation, revealing that solitons pass through shocks with altered properties and sometimes become trapped, with analytical and numerical validation.

Contribution

It provides a detailed analytical and numerical analysis of soliton-shock interactions in focusing media, highlighting new behaviors such as property changes and trapping.

Findings

Solitons pass through shocks without radiation, with altered properties.

Solitons can become trapped inside shock structures depending on initial velocity.

Analytical predictions match numerical simulations.

Abstract

Nonlinear interactions in focusing media between traveling solitons and the dispersive shocks produced by an initial discontinuity are studied using the one-dimensional nonlinear Schrodinger equation. It is shown that, when solitons travel from a region with nonzero background towards a region with zero background, they always pass through the shock structure without generating dispersive radiation. However, their properties (such as amplitude, velocity and shape) change in the process. A similar effect arises when solitons travel from a region with zero background towards a region with nonzero background, except that, depending on its initial velocity, in this case the soliton may remain trapped inside the shock-like structure indefinitely. In all cases, the new soliton properties can be determined analytically. The results are validated by comparison with numerical simulations.

Equations

$$i{q}_t+q_{xx}+2(|q{|}^2-q_o^2)q=0.$$

$${\varphi }_x=X\varphi ,{\varphi }_t=T\varphi ,$$

$$Q(x,t)=\left(\begin{array}{cc}0& q\\ -q^{\ast }& 0\end{array}\right).$$

$$q(x,t)\to q_{\pm },x\to \pm \infty ,$$

$$V_{\mathrm{zbg}}(k)=4k_{\mathrm{re}}.$$

$$V_{\mathrm{nzbg}}(k)=2(k_{\mathrm{re}}+k_{\mathrm{im}}{\lambda }_{\mathrm{re}}/{\lambda }_{\mathrm{im}}),$$

$$q_{\mathrm{step}}(x,0)=q_{o}H(x),$$

$$4{\alpha }_{\mathrm{re}}+2(q_o^2-{\alpha }_{\mathrm{im}}^2)/{\alpha }_{\mathrm{re}}=V,$$

$$\displaystyle\big{(}\alpha_{\mathrm{re}}^{2}+(q_{o}-\alpha_{\mathrm{im}})^{2}\big{)}K(m)=(\alpha_{\mathrm{re}}^{2}-\alpha_{\mathrm{im}}^{2}+q_{o}^{2})E(m).$$

Full text

Introduction

A common way to study the response of a nonlinear system is to consider Riemann problems, i.e., the evolution of a jump discontinuity between two uniform values of the initial datum. In dispersive nonlinear media, Riemann problems can give rise to dispersive shock waves (DSWs), which are non-stationary coherent wave structures and arise in many different physical contexts, including water waves, the atmosphere, optics and Bose-Einstein condensates. As a result, considerable effort has been devoted to the study of DSW formation, propagation and interactions Whitham1976 ; Kamchatnov ; GP1974 ; elgurevich ; PLA1995v203p77 ; SJAM59p2162 ; BK2006 ; HA2007 . The topic has also attracted renewed interest in recent years BEHM2016 ; EH2016 ; prl120p144101 ; sirev59p3 ; jmp50p091406 ; pre97p032218 .

An ubiquitous tool in nonlinear physics is the nonlinear Schrödinger (NLS) equation, which is a universal model for the evolution of the envelope of weakly nonlinear dispersive wave trains JMP46p133 . The NLS equation arises in a wide variety of physical settings, including deep water waves, fiber optics, plasmas and Bose-Einstein condensates AS1981 ; Agrawal2007 ; IR2000 ; PS2003 . The NLS equation is also a completely integrable infinite-dimensional Hamiltonian system AS1981 ; NMPZ1984 ; FT1987 ; APT2004 . This means that the initial value problem can be solved by the inverse scattering transform (IST) ZS1972 ; ZS1973 .

Typically, DSWs are produced either in the small dispersion, or semiclassical, limit, which arises when the dispersive effects are small compared to nonlinear ones, or in the long-time asymptotics. For the focusing NLS equation, the asymptotic behavior of solutions in the semiclassical limit with zero background (ZBG) has been studied extensively, with sech-shaped input millerkamvissis1998 ; bronskikutz1999 ; KMM2003 , generalizations thereof TVZ2004 ; lyng2012 , and box-like input JM2013 ; EKT2016 ; TE2016 . In all these cases, the focusing nonlinearity results in the formation of highly peaked oscillations in a localized region of space.

The situation is more complicated with nonzero background (NZBG), due to the presence of modulational instability (MI), namely, the fact that a constant background is unstable to long-wavelength perturbations ZO2009 . A quantitative description of the nonlinear stage of MI for generic localized perturbations of a constant background was recently obtained in BM2016 ; PRE94p060201R ; BM2017 . The corresponding behavior, which is comprised of two quiescent states separated by a central wedge with modulated periodic oscillations, was later found to arise in a broad class of NLS-type evolution equations describing a variety of focusing nonlinear media in SIREV2018 , and was recently observed experimentally in arxiv2018elrandouxsuret . Related scenarios arise from Riemann problems. Special cases of Riemann problems for the focusing NLS equation were studied in bikbaev ; BV2007 ; BKS2011 , and more general Riemann problems were recently considered in B2018 . The expanding oscillatory wedge between two uniform states can be viewed as a DSW in focusing media.

Solitons play no role in the above discussion. On the other hand, the focusing NLS equation admits a large variety of soliton solutions, both with ZBG ZS1972 and NZBG BK2014 . The combined presence of solitons and dispersive radiation in focusing media with NZBG was recently shown to produce novel phenomena such as soliton transmission, trapping and wake BLM2018 . The purpose of this work is to investigate nonlinear interactions arising when the solution contains all three of the above components, namely: NZBG, solitons and dispersive shocks. In particular, we study a practical scenario, namely, the interaction between a soliton and the oscillatory wedge formed by a discontinuity in the initial condition. We show that, when traveling from a region with nonzero background towards a region with zero background, solitons always pass through these shock structures and retain their identity, without generating dispersive radiation. Importantly, however, we also show that, even though the discrete eigenvalue in the scattering problem is of course time-independent, all of the properties of the corresponding soliton (including amplitude, velocity and shape) change once they move from a region with NZBG to one with ZBG.

A similar scenario arises when solitons travel from a region with ZBG to one with NZBG, except that now, depending its velocity, the soliton can remain trapped inside the shock-like structure indefinitely. In both cases, the new soliton properties are analytically determined by computing the long-time asymptotics of solutions, and the results are validated by comparison with extensive numerical simulations.

NLS equation, ZBG and NZBG, solitons

The cubic one-dimensional NLS equation is the partial differential equation

[TABLE]

Subscripts $x$ and $t$ denote partial differentiation, and $q(x,t)$ typically describes the complex-valued envelope of oscillations. The non-negative real parameter $q_{o}$ corresponds to a background amplitude gauge .

The IST makes crucial use of the existence of a Lax pair, namely the fact that Eq. (1) is the compatibility condition $\phi_{xt}=\phi_{tx}$ of the overdetermined linear system

[TABLE]

with $X=-ik\sigma_{3}+Q$ and $T=-i(2k^{2}+q_{o}^{2}-|q|^{2}-Q_{x})\,\sigma_{3}-2kQ$, where $\sigma_{3}=\mathop{\rm diag}\nolimits(1,-1)$ is the third Pauli matrix, and

[TABLE]

The first half of the Lax pair (2) and $q(x,t)$ are referred to as the scattering problem and the potential, respectively, while $k$ is the scattering parameter or eigenvalue of the scattering problem.

Here we study the behavior of solutions of Eq. (1) with the following step-like boundary conditions:

[TABLE]

with $q_{-}=0$, and where without loss of generality we can take $q_{+}=q_{o}$ thanks to the phase invariance of the NLS equation. We refer to Eqs. (4) as the case of one-sided nonzero background (1SNZBG). When $q_{o}=0$, they reduce to the case of ZBG ZS1972 ; AS1981 ; NMPZ1984 ; APT2004 . Conversely, when $q_{-}\neq 0$, one has a problem with a two-sided NZBG.

The IST in the case with ZBG was developed in the seminal work by Zakharov and Shabat ZS1972 . The symmetric NZBG case $|q_{-}|=|q_{+}|=q_{o}$ was developed in BK2014 , and was extended to the fully asymmetric case $|q_{\pm}|\neq 0$ and $|q_{-}|\neq|q_{+}|$ in JMP55p101505 , while the one-sided NZBG was studied in PV2015 . Recall that the IST works by associating to $q(x,t)$ time-independent scattering data via the scattering problem. Once the scattering data are obtained from the initial conditions (ICs), the solution of Eq. (1) is then reconstructed by inverting the scattering transform. The scattering data are computed in terms of the Jost eigenfunctions, which are the solutions $\phi_{\pm}(x,t,k)$ of the Lax pair (2) that reduce to plane waves as $x\to\pm\infty$, and are therefore the nonlinearization of the Fourier modes. The set of all complex values of $k$ for which the Jost eigenfunctions are defined comprises the continuous spectrum $\Sigma$ of the scattering problem. For potentials on ZBG, the continuous spectrum is simply the real $k$-axis; i.e., $\Sigma_{\mathrm{zbg}}=\mathbb{R}$ ZS1972 . For potentials on NZBG or 1SNZBG, however, the continuous spectrum acquires a subset of the imaginary axis, namely $\Sigma_{\mathrm{nzbg}}=\Sigma_{1\mathrm{snzbg}}=\mathbb{R}\cup i[-q_{o},q_{o}]$ BK2014 ; PV2015 . Moreover, for potentials on NZBG, the nonlinear analogue of the Fourier wavenumber is given by $\lambda=(q_{o}^{2}+k^{2})^{1/2}$.

The discrete spectrum of the scattering problem, when present, gives rise to soliton solutions. In particular, each discrete eigenvalue contributes one soliton to the solution. Both with ZBG and with NZBG, all the properties of the soliton are determined explicitly by the location of the discrete eigenvalue in the complex plane. These properties however differ with ZBG versus NZBG. In particular, a soliton on ZBG generated by a discrete eigenvalue $k=k_{\mathrm{re}}+ik_{\mathrm{im}}$ travels with velocity ZS1972

[TABLE]

A soliton on NZBG generated by the same discrete eigenvalue, however, travels with velocity

[TABLE]

where $\lambda=\lambda_{\mathrm{re}}+i\lambda_{\mathrm{im}}$. A contour plot of constant soliton velocity in the spectral plane for solitons on NZBG is given in Fig. 1. For solitons on ZBG, the curves of constant soliton velocity are obviously given by vertical lines, towards which the contour lines in Fig. 1 tend asymptotically as $\mathop{\rm Im}\nolimits k\to\infty$. Note also $V_{\mathrm{nzbg}}(k)>V_{\mathrm{zbg}}(k)$ for all $k$ in the first quadrant.

Interactions between solitons and DSWs: Set-up

The simplest realization of ICs consistent with the boundary conditions (4), is a “pure step” problem, namely

[TABLE]

where $H(x)$ is the Heaviside step function, defined as $H(x)=0$ for $x<0$ and $H(x)=1$ for $x>0$. The above IC results in the formation of an oscillatory “wedge” (or DSW) in the region $0<x<4\sqrt{2}q_{o}t$, to the left of which the solution is negligible and to the right of which the solution is approximately equal to the background $q_{o}$. Inside the wedge, the solution can be described as a slow modulation of the traveling wave (elliptic) solutions of the focusing NLS equation, BKS2011 (see below for details).

Here we consider situations arising from a combination of the above step ICs and a traveling soliton of the focusing NLS equation. Specifically, we consider the case of a soliton generated by a discrete eigenvalue located at $k=k_{o}$ and initially positioned at $x=X_{o}$. If $X_{o}>0$, the soliton is initially positioned to the right of the step (i.e., on a NZBG), whereas if $X_{o}<0$ it is initially positioned to the left of the step (i.e., on a ZBG). Denoting by $V_{o}$ the initial velocity of the soliton, if $\mathop{\rm Re}\nolimits(k_{o})>0$, the soliton will travel to the right (i.e., $V_{o}>0$), whereas if $\mathop{\rm Re}\nolimits(k_{o})<0$, it will travel to the left (i.e., $V_{o}<0$). Thus, if $X_{o}V_{o}>0$, the soliton will travel away from the DSW, whereas if $X_{o}V_{o}<0$ it will move towards it. The most interesting scenario is obviously the latter. We therefore consider two main cases:

   1. A left-moving soliton initially placed on NZBG (i.e., to the right of the initial discontinuity), corresponding to $V_{o}<0$ and $X_{0}>0$.

   2. A right-moving soliton initially placed on ZBG (i.e., to the left of the initial discontinuity), corresponding to $V_{o}>0$ and $X_{0}<0$.

It should be noted that, numerically, the ICs corresponding to the above situations are realized in a different way depending on whether $X_{o}$ is positive or negative. If $X_{o}<0$, one can simply add a one-soliton solution of the focusing NLS equation with ZBG to the pure step IC. If $X_{o}>0$ instead, one should multiply the pure step IC by the one-soliton solution of the focusing NLS equation with NZBG. Once the appropriate ICs have been realized, the time integration of Eq. (1) was performed using an eighth-order Fourier split-step method. As in BM2016 ; PRE94p060201R ; BM2017 ; SIREV2018 , the computed time evolution is only accurate up to the time at which round-off error grows to $O(1)$ accuracy .

Interaction between solitons and DSWs: Results

Figure 2 presents the results of numerical simulations corresponding to case 1 above (left-moving soliton initially placed to the right of the step), while Fig. 3 corresponds to case 2 (right-moving soliton initially placed to the left of the step).6 In each case, the left column shows density plots of the solution amplitude $|q(x,t)|$. while the right column shows the difference between the solution to the left and that obtained from the pure step IC (7), which provides a direct visual illustration of the nonlinear interaction effects. Each row corresponds to a different choice of discrete eigenvalue.

In case 1 (soliton initially on a NZBG and traveling leftward toward the DSW) we observe that, for all choices of discrete eigenvalue, the soliton is transmitted through the DSW and emerges as a soliton on ZBG. Importantly, however, Fig. 2 clearly shows that all of the properties of the soliton (that is, its amplitude, width, velocity, and its breather-like versus traveling-wave nature) are changed after it has traveled through the oscillatory structure.

The change of the soliton features may be surprising, since the properties of the soliton are completely determined by the location of the discrete eigenvalue, and both the continuous and discrete spectrum of the scattering problem are independent of time. As we show below, however, these changes are not a numerical artifact, and are indeed reflective of the true nonlinear dynamics of the system.

As shown in Fig. 3, a similar outcome arises in case 2 (soliton initially on a ZBG and traveling rightward toward the DSW). In addition, however, here there are cases (e.g., second and third row) in which the soliton does not escape the oscillatory wedge, and remains trapped there forever, as is easily seen by comparing the soliton trajectory with that of the wedge boundary escape .

Long-time asymptotics

We emphasize that none of the velocities in the trajectories shown in the right column of Figs. 2 and 3 were determined numerically, and all of them are determined analytically instead. Indeed, we next show that the numerical results of Figs. 2 and 3 can be fully characterized by studying the long-time asymptotics of solutions of the focusing NLS equation with nonzero background.

A full calculation of the long-time asymptotics is beyond the scope of this work, so here we limit ourselves to presenting the essential details. Recall that, in the inverse problem of IST, the solution of the NLS equation is obtained from that of a suitable matrix Riemann-Hilbert problem (RHP) defined in terms of the reflection coefficient and, if present, the discrete spectrum. A key part of the RHP is a controlling phase function appearing in the jump conditions. In the case of ZBG, this phase function is $\theta_{\mathrm{zbg}}(x,t,k)=k(x-2kt)$ APT2004 , whereas with NZBG, $\theta_{\mathrm{nzbg}}(x,t,k)=\lambda\,(x-2kt)$ instead BK2014 , with $\lambda(k)$ as above. Indeed, it is precisely by looking along directions $x=Vt$ and setting $\mathop{\rm Im}\nolimits[\theta(x,t,k_{o})]=0$ that one finds the soliton velocities corresponding to a discrete eigenvalue at $k=k_{o}$ with ZBG and NZBG.

Importantly, however, in the long-time asymptotics of solutions with symmetric NZBG in the presence of a small disturbance initially placed near $x=0$, the governing phase function gets modified in the wedge $x\in(-4\sqrt{2}q_{o}t,4\sqrt{2}q_{o}t)$ BM2016 ; BM2017 . For those values of $x$, $\theta(x,t,k)$ is replaced by a new phase function $h(x,t,k)$ defined in terms of certain Abelian integrals BM2017 . It was also shown in BLM2018 that it is by setting $\mathop{\rm Im}\nolimits[h(x,t,k_{o})]=0$ that one determines the velocity of a soliton inside the wedge. For the one-sided NZBG studied here, the same considerations apply in the half wedge $x\in(0,4\sqrt{2}q_{o}t)$.

Classification of interactions

For a left-moving soliton starting on a NZBG, the plots in the right column of Fig. 2 clearly show that, for the choices of $k_{o}$ considered, the soliton velocity is given by $V_{\mathrm{nzbg}}(k_{o})$ before entering the oscillation structure, and by $V_{\mathrm{zbg}}(k_{o})$ after exiting it.

The situation is more complex for a right-moving soliton starting on a ZBG. In this case, one observes different outcomes depending on the precise location of the discrete eigenvalue. Consider a discrete eigenvalue $k_{o}$ in the first quadrant of the complex plane, and recall the contour plot of soliton velocity with NZBG in Fig. 1. One can distinguish four domains: in $\smash{D_{1}^{+}}$ (purple region) and $\smash{D_{4}^{+}}$ (blue region), one has $V_{\mathrm{nzbg}}(k_{o})>4\sqrt{2}q_{o}$, whereas in in $\smash{D_{2}^{+}}$ (yellow region) and $\smash{D_{3}^{+}}$ (gray region), one has $V_{\mathrm{nzbg}}(k_{o})<4\sqrt{2}q_{o}$. Figure 3 shows the results obtained from a discrete eigenvalue located in each of these domains. (The domains $\smash{D_{1}^{-}},\dots,\smash{D_{4}^{-}}$ used in Fig. 2 are the symmetric counterparts to $\smash{D_{1}^{+}},\dots,\smash{D_{4}^{+}}$ relative to the imaginary axis.) The difference between $\smash{D_{1}^{+}}$ and $\smash{D_{4}^{+}}$ is that the latter collects eigenvalues close to the branch cut, and results in broader solitons compared to the former BLM2018 . The difference between $\smash{D_{2}^{+}}$ and $\smash{D_{3}^{+}}$ originates from the long-time asymptotics.

Recall that the elliptic solutions of the focusing NLS equation are determined (up to translations and phase invariance) by four complex conjugate constants, which are the branch points of the associated Riemann surface in the IST Kamchatnov . For the boundary conditions (4) with $q_{-}=0$, two of these branch points are fixed at $\pm iq_{o}$, whereas the other two, $\alpha=\alpha_{\mathrm{re}}+i\alpha_{\mathrm{im}}$ and its conjugate, are given by the system of modulation equations elgurevich

[TABLE]

Here, $m={4\alpha_{\mathrm{im}}q_{o}}/{|\alpha-q_{o}|^{2}}$, while $K(m)$ and $E(m)$ are the complete elliptic integrals of the first and second kind, respectively NIST . The trajectory described by $\alpha$ in the complex plane as $V$ varies between 0 and $4\sqrt{2}q_{o}$ is shown by the blue dashed curve in Fig. 1 It is this curve that defines the boundary between $\smash{D_{2}^{+}}$ and $\smash{D_{3}^{+}}$. As shown in elgurevich ; BM2016 ; PRE94p060201R , the same slow modulation of the traveling wave solutions of the focusing NLS equation also describes the nonlinear stage of MI induced by localized perturbations of a constant background.

Importantly, each of the domains $D_{1}^{+},\dots,D_{4}^{+}$ results in a different outcome for the nonlinear interaction. The four cases shown in Fig. 3 correspond to a choice for the discrete eigenvalue $k_{o}$ in each of the four domains (the precise location is identified by the red circles in Fig. 1).

In all four cases, the soliton initially travels towards the oscillatory structure with velocity $V_{\mathrm{zbg}}(k_{o})$. The simplest case is that of $k_{o}\in\smash{D_{1}^{+}}$ (top row of Fig. 3). Here the soliton is transmitted through the DSW, and emerges as a soliton on a nonzero background. However, its velocity is different after the interaction, and is given by $V_{\mathrm{nzbg}}(k_{o})>4\sqrt{2}q_{o}$.

When $k_{o}\in\smash{D_{2}^{+}}$ (second row of Fig. 3), we have $V_{\mathrm{zbg}}(k_{o})<V_{\mathrm{nzbg}}(k_{o})<4\smash{\sqrt{2}}q_{o}$. Here the soliton is not transmitted through the DSW, and remains trapped inside the wedge. The velocity $V_{s}(k_{o})$ of the trapped soliton is obtained by solving the equation $h_{\mathrm{im}}(k_{o},V)=0$, with $h(k,V)=h(Vt,t,k)/t$, which has a unique root $V=V_{s}(k_{o})$ for $k_{o}\in\smash{D_{2}^{+}}$ BLM2018 .

When $k_{o}\in\smash{D_{3}^{+}}$ (third row of Fig. 3), we also have $V_{\mathrm{zbg}}(k_{o})<V_{\mathrm{nzbg}}(k_{o})<4\sqrt{2}q_{o}$, and the soliton is again trapped inside the wedge. Here however the equation $h_{\mathrm{im}}(k_{o},V)=0$ has two solutions, for $V=V_{\ast}(k_{o})$ and $V=V_{s}(k_{o})$ with $V_{\ast}(k_{o})<V_{s}(k_{o})$ BLM2018 . The second of these roots corresponds to the trapped soliton, the first to the soliton-generated wake BLM2018 ; limitations .

Finally, when $k_{o}\in\smash{D_{4}^{+}}$ (fourth row of Fig. 3), we have $V_{\mathrm{zbg}}(k_{o})<4\sqrt{2}q_{o}$ but $V_{\mathrm{nzbg}}(k_{o})>4\sqrt{2}q_{o}$. Here the soliton is transmitted through the wedge and eventually re-emerges with speed $V_{\mathrm{nzbg}}(k_{o})$ limitations . However, the equation $h_{\mathrm{im}}(k_{o},V)=0$ has a solution for $V=V_{\ast}(k_{o})$ which gives rise to a soliton-generated wake BLM2018 .

Discussion

We emphasize that there is no data fitting in Figs. 2–3. Thus, the figures demonstrate excellent agreement between the numerically computed soliton velocities and those obtained from the long-time asymptotics. In terms of the inverse problem in the IST, the properties of the soliton depend on whether the controlling phase function is the one with ZBG or that with NZBG. Soliton trapping by a dispersive shock was also recently discussed in Refs. prl120p144101 ; BLM2018 , while soliton trapping by a rarefaction wave was considered in Ref. jmp50p091406 .

The dynamics remain virtually unchanged if the sharp discontinuity in the IC is replaced by a smooth function which interpolates between the asymptotic values of the potential, demonstrating the robustness of the results. Note also that the dynamical behavior produced by a pure step IC is markedly different in the focusing and defocusing case. Namely, instead of the oscillatory wedge arising here, in the defocusing case the pure step IC (7) gives rise to Gibbs-like phenomena difrancomclaughlin ; biondinitrogdon .

The numerical results do not allow us to draw any conclusions regarding the possible presence of a soliton-generated wake in interactions between left-moving solitons and the DSW. We were unable to observe a wake, but we cannot exclude its existence a priori. A definitive answer can be obtained through a rigorous calculation of the long-time asymptotics, which however is beyond the scope of this work.

The results of this work open up a number of interesting avenues for further research. From a theoretical point of view, an obvious open problem is the rigorous validation of the results of this work by explicit computation of the long-time asymptotic behavior of solutions, for example using the Deift-Zhou nonlinear steepest descent method for oscillatory Riemann-Hilbert problems DZ1993 . From a practical point of view, another obvious question is whether the results of this work can be generalized to other NLS-type models, as was previously done in SIREV2018 for the results obtained in BM2016 . Finally, from an even more practical point of view, an obvious challenge will be the experimental observation of the results of this work, perhaps in nonlinear optical fiber experimts, similarly to those recently conducted in arxiv2018elrandouxsuret .

Acknowledgement

We are grateful to D. Mantzavinos and S. Li for many insightful discussions. This work was supported in part by the National Science Foundation under grant DMS-1614623.