This paper investigates dynamical behaviors and controllability of some nonautonomous localized waves based on the Gross–Pitaevskii equation with attractive interatomic interactions. Our approach is a relation constructed between the Gross–Pitaevskii equation and the standard nonlinear Schrödinger equation through a new self-similarity transformation which is to convert the exact solutions of the latter to the former’s. Subsequently, one can obtain the nonautonomous breather solutions and higher-order rogue wave solutions of the Gross–Pitaevskii equation. It has been shown that the nonautonomous localized waves can be controlled by the parameters within the self-similarity transformation, rather than relying solely on the nonlinear intensity, spectral parameters, and external potential. The control mechanism can induce an unusual number of loosely bound higher-order rogue waves. The asymptotic analysis of unusual loosely bound rogue waves shows that their essence is energy transfer among rogue waves. Numerical simulations test the dynamical stability of obtained localized wave solutions, which indicate that modifying the parameters in the self-similarity transformation can improve the stability of unstable localized waves and prolong their lifespan. We numerically confirm that the rogue wave controlled by the self-similarity transformation can be reproduced from a chaotic initial background field, hence anticipating the feasibility of its experimental observation, and propose an experimental method for observing these phenomena in Bose–Einstein condensates. The method presented in this paper can help to induce and observe new stable localized waves in some physical systems.

The exact nonlinear localized wave solutions of nonlinear integrable equations, such as solitons, breathers, and rogue waves, can characterize nonlinear processes in diverse domains. These solutions are widely recognized as significant physical models. The rogue wave refers to the abrupt and unexpected appearance of a surface wave on a plane wave. The breather phenomenon can be considered a recurring manifestation of rogue waves. There exists evidence supporting the occurrence of rogue waves and breathers in various physical systems, including nonlinear optics, plasmas, and Bose–Einstein condensates. Hence, the effective control of localized waves can offer significant assistance in facilitating their experimental observation. This study aims to investigate the localized wave solutions of the Gross–Pitaevskii equation by using a self-similarity transformation. The self-similarity transformation parameters have been demonstrated to exert control on the localized waves and induce some unusual higher-order rogue waves. We test the stability of these localized waves by numerical experiments and confirm that they can be excited in a chaotic background field. This provides evidence for the experimental observation of these localized waves.

## I. INTRODUCTION

Localized waves are typical nonlinear waves in integrable systems such as the nonlinear Schrödinger equation (NLSE),^{1–5} including soliton, breather, and rogue wave (RW). The generations of these waves are usually thought to be the interaction of nonlinearity and dispersion.^{6} Dark and bright solitons are two types of solitons that are often found. In the context of Bose–Einstein condensates (BECs), the dark soliton represents a hole in the particle distribution, whereas the bright soliton is seen as a representation of the particle density distribution.^{7,8} Experimental evidence demonstrated that the utilization of the Feshbach resonance technique could excite the formation of stable bright or dark solitons in ultracold $ 85$Rb BECs and $ 7$Li BECs systems.^{9,10} When solitons propagate on a finite background and periodically appear peak and depression, one can obtain breathers, including the Akhmediev breather (AB)^{11,12} and the Kuznetsov–Ma breather (K–MB).^{12,13} Breathers are thought to be the origin of RWs.^{14} RWs are steep, transient, and localized, with peaks that are more than twice as high as the surrounding background.^{15,16} They were observed in a variety of physical systems, including optical fibers,^{17–20} water tanks,^{21,22} electrical lattices,^{23,24} capillary waves,^{25} and finance systems.^{26} These waves are associated with the manifestation of certain ultrashort wave phenomena. The RW is often understood as a concentrated atomic density distribution in BECs.^{27,28}

Integrable systems, such as the NLSE,^{1–5} the Ginzburg–Landau equation,^{29–31} and the Gross–Pitaevskii equation (GPE),^{28,32–35} are fundamental partial differential equations that incorporate both dispersion and nonlinearity. These equations enable the generation of localized wave excitations. The well-known GPE, a class of the NLSE deformations, contains an external potential and variable nonlinear coefficient.^{7,36,37} In BECs, nonlinearities of the GPE are due to interatomic interactions, which are defined by the $s$-wave scattering length.^{38–40} Some comprehensive representation of various phenomena, including electric fields in optical fibers,^{41} Langmuir waves in plasmas,^{42–45} deformation waves in seas,^{15} and BECs,^{46–48} has benefited from the key role of NLSE/GPE in mathematical physics in recent decades. Since the discovery of BECs in 1995,^{49} the GPE has been widely employed as a theoretical framework for investigating the experiment of atomic matter waves in BECs. This includes various phenomena such as bright solitons,^{9,10} dark solitons,^{9,50} gap solitons,^{51} gray solitons,^{52} RWs,^{27} vortices,^{53,54} quantized vortices,^{55} vortex lattices,^{56} Bloch oscillations,^{57,58} Landau–Zener tunneling,^{59} Mott insulator,^{60} and compression of condensation.^{61} Theoretical investigations on the GPE demonstrated that BECs could help in the excitation of localized waves.^{27,48} Therefore, there is a growing interest in the exact and numerical localized wave solutions to the GPE.^{28,32–35}

The scattering length, a significant element for all BEC quantities, directly affects the structure of localized waves. Recent experiments have shown evidence that the Feshbach resonance can alter the scattering length,^{62,63} while also offering an efficient means of manipulating the waveform through the Feshbach resonance technique.^{64} In relevant mathematical physical models, including the NLSE/GPE, the nonlinear coefficients represent the scattering length.^{62–64} Several findings indicate that the spectral parameter plays a crucial role in determining the propagation characteristics of solitons and breathers, and the finite background height and frequency dominate the waveforms of the RW.^{65,66} Hence, the amplitude, width, and waveform of the localized waves are predictable when the spectral parameter and background are fixed. The higher-order RW usually includes the tightly bound and the loosely bound. The tightly bound higher-order RW is symmetrical before and after the peak, but the loosely bound higher-order RW is frequently made up of a set number of lower-order RWs, each with the same structure.^{66–68} The number of first-order RWs in the loosely bound higher-order RW is contingent on its order and typically constant.^{66} Due to the variable nonlinear coefficient, the localized waves in the GPE will show some nonautonomous structures, which may lose stability and undergo compression during propagation.^{34,35} Nevertheless, the conventional analytical method cannot easily change the previous findings. The control of the localized wave is inseparable from these conclusions.

To improve the limitations, we propose a new method to control them for a GPE. Our idea is to propose a new transformation, which makes the localized waves no longer depend on the spectral parameter, background, and nonlinear intensity (or scattering length) as the NLSE, and it also makes the model studied accompanied by some new localized wave phenomena. Highly localized matter waves may exhibit particularly robust dynamics in condensate’s density and energy in the context of BECs. Therefore, the adequate management of localized waves is worth studying, and their management methods are also a great challenge. It would be beneficial to investigate condensate dynamics with more extensive control of localized waves in the GPE.

Based on a GPE in this paper, we propose a new control of the localized waves and give some novel localized wave dynamical structures via the self-similarity transformation and generalized Darboux transformation (DT). In Sec. II, the model, analytical method, and exact solutions are given. The exact breather solutions and RW solutions of the GPE are discussed in Secs. III and IV, respectively. The numerical simulations of exact solutions and applications in experimental observations are presented in Sec. V. Section VI gives a summary of this paper.

## II. PRELIMINARIES

### A. Model

^{7,69,70}

^{28,71}Thus, the 3D GPE (1) can reduce to the following dimensionless quasi-1D GPE with time-dependent:

^{71}

^{,}

^{71}

^{,}

^{72–75}Equation (2) can describe the dynamics of BECs with the attractive interaction in an expulsive potential, and it is Lax integrable with a Lax pair as

^{76}For Eq. (2), the dynamics of one-breather and first-order RW had been reported.

^{28,71}However, the higher-order localized wave phenomena have yet to be reported. Therefore, in Secs. III–V, we will focus on distinct localized wave structures of Eq. (2) and investigate some unreported localized wave dynamics.

### B. Analytical methods

This part will present the computation method for the exact solutions to Eq. (2). Despite the fact that the Lax pair (4) exists in Eq. (2) and can be solved using the DT, the variable coefficients make it difficult to find the fundamental solution. We so put forward the following theorem:

In this instance, transformation (8) is the self-similarity transformation, which can be used to convert exact solutions of Eq. (7) to ones of Eq. (2). In the transformation, we introduce several free parameters ( $ c 0$, $ c 1$, $ c 2$, $ c 3$, and $ c 4$). Even though their values may be to zero, we decided to preserve them since they provide new control over localized waves. The self-similarity transformation is a widely used auxiliary method for solving low- and high-dimensional one-component and multi-component GPEs,^{77–86} such as 1D GPEs,^{77–80} two-dimensional GPEs,^{81–83} and 3D GPEs.^{84–86} However, none of these existing studies have regulated the new phenomena of breather and RWs for GPE (2). It is worth mentioning that the self-similarity transformation can be applied to any one- and high-dimension equation, and it can theoretically induce new nonlinear localized wave phenomena. The find in this paper indicates that GPE (2) can give rise to new nonlinear wave phenomena under the influence of self-similarity transformation. Sections III and IV will explore how these parameters in the self-similarity transformation might affect localized waves.

^{4}

Theorem 2 is a generalized DT in the matrix form.^{66} This is a new generalized DT in the NLSE and is different from generalized DT in the iterative form.^{5} The localized wave solutions of Eq. (7), including breather, and higher-order RW, can be obtained. Then, via the transformation (8), we can give the exact solutions of Eq. (2).

### C. Rogue wave and breather solutions

**Case 1: Rogue wave solutions**

For $N=3$, the generalized $(1,2)$-fold DT and transformation (8) are used to achieve the third-order RW solution of Eq. (2). Due to the overlong expressions, we omit them here.

**Case 2. Breather solutions**

^{87,88}In this work, due to the physical structure of the solution, we believe that it is also a breather.

Various localized wave solutions $q(\xi ,\tau )$ for Eq. (7) can be derived from the generalized $(m,N\u2212m)$-fold DT. Then, the exact localized wave solutions of Eq. (2) may be obtained by combining the self-similarity transformation (8).

Here is emphasized that $ |\psi |$, $ | \psi 1 |$, $ | \psi 2 |$, and $ | \psi 3 |$ in all after figures stand for $ |\psi (x,t) |$, $ | \psi 1(x,t) |$, $ | \psi 2(x,t) |$, and $ | \psi 3(x,t) |$, respectively.

## III. DYNAMICS OF BREATHER SOLUTIONS

### A. One-breather solutions

**Akhmediev breather.**When $ \zeta 1= 1 2 i$, $ C 1=\u2212 C 2=1$, the one-breather solution (34) reads

Increase values of $ c 0$, $\lambda $, or $ | c 4 |$, the AB moves along the negative $t$-axis, as shown in Line $3$ ( $ c 0=\u221220$ to $ c 0=\u22125$), Line $4$ ( $\lambda =0.01$ to $\lambda =0.03$), and Line $5$ ( $ c 4=\u22120.07$ to $ c 4=\u22120.25$), respectively.

Increase the value of $ c 1$, the phase of the AB moves along $x$-axis, and the structure is invariant as shown in Line $2$ ( $ c 1=0$ to $ c 1=2$).

Increasing the value of $ c 0$ or reducing the value of $\lambda $ results in an increase in the period of the AB, as shown in Line $3$ and Line $4$, respectively.

Increasing the value of $ c 0$ or reducing values of both $ | c 4 |$ and $\lambda $ results in a decrease in the amplitude of the AB, as shown in Line $3$, Line $5$, and Line $4$, respectively.

When the value of $ c 3$ is changed from $0$ to $ 1 8$, the resulting AB is shown in Fig. 1(a2), exhibiting a distinct trace compared to Fig. 1(a1). There is no need to modify the nonlinear strength, dispersion, and background frequency to control the trajectory of the breather. Based on Eqs. (36) and (37), it can be seen that the parameter $ c 2$ does not have an impact on the AB. Furthermore, the AB fails to hold when $ c 0\u22650$. Consequently, these cases are not considered for further discussion.

**Kuznetsov–Ma breather.**When $ \zeta 1= 4 3 i$, $ C 1=\u2212 C 2=1$, the one-breather solution (34) reads

Increase the value of $ c 0$, the K–MB moves along the negative $t$-axis, as seen in Fig. 2(b2) ( $ c 0=\u221220$ to $ c 0=0$).

Increase the value of $ c 1$, the K–MB propagates from the negative $x$-axis to the center of the field, as seen in Fig. 2(b3) ( $ c 1=0$ to $ c 1=2$).

Increase the value of $ c 3$, the K–MB propagates to the positive $x$-axis, as seen in Fig. 2(b4) ( $ c 3=0$ to $ c 3= 1 30$).

Reducing the value of $ c 4$ or $\lambda $ results in an increase in the period of K–MB, as seen in Fig. 2(b5) ( $ c 4= 1 20$ to $ c 4= 4 25$) and Fig. 2(b6) ( $\lambda = 1 100$ to $\lambda = 1 150$).

The K–MB remains unaffected by changes in parameter $ c 2$.

### B. Two-breather solutions

We can also give two-breather solution via generalized $(1,N\u22121)$-fold DT. The two-breather solution (29) can be calculated by Eq. (29) with $ q 0=c e i [ a \xi + ( c 2 \u2212 1 2 a 2 ) \tau ]$, $ \zeta 1=2 i$, $ \varphi 1 ( 0 )=\u2212 i c n \u2212 e A$, $ \gamma 1 ( 0 )=( e m\u2212 e \u2212 m) e \u2212 A$, $ \varphi 1 ( 1 )= i c 2 c + 1(b+ i) n + e A$, $ \gamma 1 ( 1 )= 1 2 c + 1( e m+ e \u2212 m)b e \u2212 A$, here $ n \xb1=(c+ 2 c + 1+1) e \u2212 m\xb1(c\u2212 2 c + 1+1) e m$, $m= 2 c + 1[( ic+ i\u2212a)\tau +\xi ]$, $A= i 4( a 2\tau \u22122 c 2\tau \u22122a\xi )$, $b=( iac+ ia+ c 2+4c+2)\tau \u2212 i(c+1)\xi $, which consists of rational polynomials and exponential functions. Figure 1(a4) shows the structure of the two-breather solution (29) for $c=1$, $a=0$, $ C 2=\u2212 C 1=\u22121$, $ c 0=\u22122$, $ c 1= c 2= c 3=0$, $ c 4=\u2212 1 5$, $ g 0=2$, $\lambda = 1 25$. It forms two parallel periodic waves after the interaction, which can be regarded as the non-interacting breather bound state.

## IV. DYNAMICS OF ROGUE WAVE SOLUTIONS

### A. First-order rogue wave solutions

Figures 3(a1), 3(a2), and 3(a3) show the structures of the first-order RW solution (41). The first-order RW often exhibits transitory behavior and has a limited duration. In Fig. 3(a1), a localized wave structure with an eye-shaped pattern is seen. This structure features a peak and two depressions. The parameters used in this observation are as follows: $c=1$, $ c 1= c 2= c 3=0$, $ c 4= 1 2$, $ g 0=2$, $a=0$, $ c 0=\u221214$, $\lambda = 1 50$. The second row in Table I displays the peak and depression positions, amplitude, width, and lifespan. When the aspect ratio $\lambda $ is decreased to $\lambda = 1 110$, the structure of the first-order RW is displayed in Fig. 3(a2). The RW exhibits an amplitude $A=0.90$, a width $W=4.86$, and a lifespan $t\u224830$. This indicates that the first-order RW moves in the negative t-axis direction, while experiencing an increase in width and lifespan, and a decrease in amplitude. These data of the RW structures are presented in the second and third rows of Table I. When the background frequency is changed from $a=0$ to $a=\u2212 3 2$, the trajectory of the RW transforms from a straight line to a curve, as shown by comparing Figs. 3(a1) and 3(a3). The physical characteristics of first-order RW are shown in the fourth row of Table I. When the value of $\u2212\lambda c 0$ decreases, the width of RW rises. Consequently, there is a decrease in amplitude and an extension of the lifespan. When $ c 0\u22650$, solution (41) is no longer an RW, and this case will not be further discussed in this context.

Rogue waves . | Peaks points (x, t)
. | Depression points (x, t)
. | Amplitudes . | Widths . | Lifespan . |
---|---|---|---|---|---|

Fig. 3(a1) | $ ( 5 14 7 , 50 ln \u2061 2 7 5 )$ | $ ( 5 7 14 ( 1 \xb1 3 ) , 50 ln \u2061 2 7 5 )$ | $ 3 10 5 7 4$ | $ 5 7 21$ | t ≈ 9 |

≈(0.94, 2.83) | ≈(2.58, 2.83),( − 0.69, 2.83) | ≈1.09 | ≈3.27 | ||

Fig. 3(a2) | $ ( 385 14 , 110 ln \u2061 2 385 55 )$ | $ ( 385 14 ( 1 \xb1 3 ) , 110 ln \u2061 2 385 55 )$ | $ 3 110 55 3 4 7 4$ | $ 1 7 3 385$ | t ≈ 30 |

≈(1.40, − 37.13) | ≈(3.83, − 37.13),( − 1.03, − 37.13) | ≈0.90 | ≈4.86 | ||

Fig. 3(a3) | $ ( 5 2 2 , 50 ln \u2061 2 5 )$ | $ ( 5 2 2 ( 1 \xb1 3 ) , 50 ln \u2061 2 5 )$ | $ 3 5 20 2 3 4$ | $5 6$ | t ≈ 200 |

≈(3.54, − 63.16) | ≈(9.66, − 63.16),( − 2.59, − 63.16) | ≈0.56 | ≈12.25 |

Rogue waves . | Peaks points (x, t)
. | Depression points (x, t)
. | Amplitudes . | Widths . | Lifespan . |
---|---|---|---|---|---|

Fig. 3(a1) | $ ( 5 14 7 , 50 ln \u2061 2 7 5 )$ | $ ( 5 7 14 ( 1 \xb1 3 ) , 50 ln \u2061 2 7 5 )$ | $ 3 10 5 7 4$ | $ 5 7 21$ | t ≈ 9 |

≈(0.94, 2.83) | ≈(2.58, 2.83),( − 0.69, 2.83) | ≈1.09 | ≈3.27 | ||

Fig. 3(a2) | $ ( 385 14 , 110 ln \u2061 2 385 55 )$ | $ ( 385 14 ( 1 \xb1 3 ) , 110 ln \u2061 2 385 55 )$ | $ 3 110 55 3 4 7 4$ | $ 1 7 3 385$ | t ≈ 30 |

≈(1.40, − 37.13) | ≈(3.83, − 37.13),( − 1.03, − 37.13) | ≈0.90 | ≈4.86 | ||

Fig. 3(a3) | $ ( 5 2 2 , 50 ln \u2061 2 5 )$ | $ ( 5 2 2 ( 1 \xb1 3 ) , 50 ln \u2061 2 5 )$ | $ 3 5 20 2 3 4$ | $5 6$ | t ≈ 200 |

≈(3.54, − 63.16) | ≈(9.66, − 63.16),( − 2.59, − 63.16) | ≈0.56 | ≈12.25 |

Figure 3(a4) is a phase diagram showing the position of first-order RW under different parameters $ c i$ and $\lambda $,

Increase values of $ c 1$ and $ c 3$, the RW moves in space along the negative and positive $x$-axis, respectively.

Increase values of $ c 0$ and $ c 4$, the RW moves along both the positive $x$-axis and negative $t$-axis.

Increase the value of $\lambda $, the RW moves along the positive $(x,t)$-axis.

The RW remains unaffected by changes in parameter $ c 2$.

### B. Second-order rogue wave solutions

### C. Third-order rogue wave solutions

When $N=3$, the third-order RW solution can be calculated by the generalized $(1,2)$-fold DT and transformation (8) and is omitted due to its complexity. Figure 5 plots the tightly and loosely bound third-order RW. Figure 5(a1) exhibits the tightly bound third-order RW with parameters $a= c 1= c 2= c 3=0$, $c=1$, $ g 0=2$, $\lambda = 1 50$, $ c 0=\u221214$, $ c 4= 2 3$, $ b 1= d 1= b 2= d 2=0$. When we change $ c 0$ and $ c 4$, a long-lifespan tightly bound third-order RW is plotted in Fig. 5(a2). These structures are similar to the second-order RW, but the number of low-energy peaks around the peaks is different. When we choose the nonzero $ b 1$, $ b 2$, $ d 1$, or $ d 2$ (e.g., $ b 1=100$, $ d 1= b 2= d 2=0$) the tightly bound third-order RW will divide into six loosely bound low-order RWs [see Fig. 5(a3)]. Change the parameters $ c 0$ and $ c 4$ on the basic of the parameters in Fig. 5(a3), the tightly bound third-order RW can divide into five first-order RWs [e.g., $ c 0=\u221210$ and $ c 4= 5 3$, Fig. 5(a4)], four first-order RWs [e.g., $ c 0=\u22123$ and $ c 4= 8 3$, Fig. 5(a5)], three first-order RWs [e.g., $ c 0=2$ and $ c 4= 14 3$, Fig. 5(a6)]. As time increases, the single RWs in the loosely bound RW have a compression, increasing in amplitude and decreasing in width and lifespan. Due to the external potential, it is no longer guaranteed that each single RW has the same structure concerning each other.

The following is a concise overview of the content covered in this particular section: (i) The parameters $ c i$ ( $i=0,1,3,4$) in the self-similarity transformation (8) can regulate the breather and RW. The background, nonlinear strength (scattering length), and external potential need no modification. (ii) Each first-order RW in the loosely bound higher-order RW has an increase in amplitude and a compression in width with the increasing of time. (iii) The quantity of these loosely bound $N$th-order RWs is not constant within the field, including at least $N$ and at most $ N 2 + N 2$ first-order RWs, different from the existent conclusion.^{66}

## V. NUMERICAL SIMULATIONS AND APPLICATIONS

### A. Dynamical stability

We build a numerical simulation to test and validate the steady propagation of the previously obtained breather and RW solutions of the GPE. In this study, we will use a second-order finite difference method to analyze the wave propagation characteristics of these localized waves. It is noteworthy to emphasize that our numerical experimental indicate the validity of the central difference for the variable coefficient system. In the next discussion, we will examine the following cases:

* Case 1. One-AB.* We choose a waveform $ |\psi (x,170) |$ in a spatial domain as the initial condition, where the parameters of $ |\psi (x,t) |$ are the same as in Fig. 1(a1). Figure 6(a1) exhibits the numerical evolution result of $ |\psi (x,170) |$. When a low-amplitude white noise ( $\u22483\xd7 10 \u2212 2$) is added to the initial conditions, the corresponding numerical result is shown in Fig. 6(b1), which exhibits a great oscillation, and indicates the AB is unstable.

* Case 2. Two-breather.* Choose waveforms $ |\psi (x,\u221212) |$ and $ |\psi (x,12) |$ as the initial conditions for the spatial domain, with the parameters of $ |\psi (x,t) |$ being consistent with those shown in Fig. 1(a4). Subsequently, Figs. 6(a2) and 6(a3) display the numerical evolution results for $ |\psi (x,\u221212) |$ and $ |\psi (x,12) |$, respectively. A low-amplitude white noise ( $\u22486\xd7 10 \u2212 2$) is introduced to the initial condition $ |\psi (x,\u221212) |$, and the numerical results are shown in Fig. 6(b2). Another low-amplitude white noise ( $\u22483\xd7 10 \u2212 2$) is introduced to the initial condition $ |\psi (x,12) |$. The numerical results are shown in Fig. 6(b3). The findings suggest that the long period region exhibits stability for a minimum duration of $t=35$ [ranging from $t=\u221212$ to $t=23$, as shown in Fig. 6(b2)]. Conversely, the short period region can only maintain stability for a duration of less than $t=14$ [ranging from $t=12$ to $t=26$, as shown in Fig. 6(b3)], even when subjected to a lower perturbation. The dynamics of the one-K–MB exhibit identical characteristics. Therefore, increasing the period of one- and two-breather can ensure the stability of the breather. Reducing $\lambda $ and $ | c 4 |$ can enhance the stability of the breather.

* Case 3. Loosely bound second-order RW.* Consider the waveforms $ |\psi (x,\u221215) |$ and $ |\psi (x,\u221250) |$ in a spatial domain as the initial condition, where $ |\psi (x,t) |$ is the loosely bound second-order RW, parameters are, respectively, the same as in Figs. 3(a7) and 3(a8), then, Figs. 7(a1) and 7(a2) exhibit the numerical evolution results of initial states $ |\psi (x,\u221215) |$ and $ |\psi (x,\u221250) |$, respectively. When a low-amplitude white noise ( $\u22481\xd7 10 \u2212 2$) is added to these initial conditions, the shorter-lived RW can maintain the stability of less than $t=19$ [ranging from $t=\u221221$ to $t=\u22122$, see Fig. 7(b1)], but the longer-lived RW can maintain the stability at least $t=40$ [ranging from $t=\u221250$ to $t=\u221210$, see Fig. 7(b2)].

* Case 4. Tightly bound third-order RW.* We choose the waveforms $ |\psi (x,\u221220) |$ and $ |\psi (x,\u221260) |$ in a spatial domain as the initial condition, where $ |\psi (x,t) |$ is the tightly bound third-order RW, parameters are, respectively, the same as in Figs. 5(a1) and 5(a2), then, Figs. 7(a3) and 7(a4) exhibit the numerical evolution results of $ |\psi (x,\u221220) |$ and $ |\psi (x,\u221260) |$, respectively. When a low-amplitude white noise ( $\u22485\xd7 10 \u2212 2$) is added to the initial conditions, the shorter-lived third-order RW can maintain the stability of less than $t=15$ [ranging from $t=\u221220$ to $t=\u22125$, see Fig. 7(b3)], but the longer-lived RW can maintain stability of at least $t=30$ [ranging from $t=\u221260$ to $t=\u221230$, see Fig. 7(b4)].

* Case 5. Loosely bound third-order RW.* Choose the initial conditions $ |\psi (x,\u221210) |$ and $ |\psi (x,\u2212145) |$, where $ |\psi (x,t) |$ is the loosely bound third-order RW, parameters are, respectively, the same as in Figs. 5(a3) and 5(a6), then, Figs. 7(a5) and 7(a6) exhibit the numerical evolution results of unperturbed initial conditions $ |\psi (x,\u221210) |$ and $ |\psi (x,\u2212145) |$, and Figs. 7(b5) and 7(b6) show the numerical evolution results of perturbed initial conditions by a low-amplitude white noise ( $\u22481\xd7 10 \u2212 2$). The shorter-lived RW can maintain the stability of less than $t=15$ [ranging from $t=\u221210$ to $t=5$, see Fig. 7(b5)], but the longer-lived RW can maintain the stability of at least $t=25$ [ranging from $t=\u2212145$ to $t=\u2212120$, see Fig. 7(b6)].

Numerical simulations show that producing multiple stable condensates in the field is possible, and they can be highly localized and characterized by extremely high energies. These novel RW phenomena provide a new perspective on the manipulation of the number, energy, and phase of the condensate. These results show that the dynamical stability of the tightly bound RW is better than that of the loosely bound RW, while the loosely bound RW is usually more unstable, and the results are the same as previous studies.^{66,89,90} However, when the lifespan of the RW prolongs, either tightly or loosely bound, the dynamical stability can be improved, allowing the loosely bound RW to have stable dynamics. Hence, adjusting parameters $ c i$ to prolong the lifespan of the RW can enhance the stability of RW, which can realize the unstable loosely bound higher-order RW exists steadily.

### B. Applications

To be closer to real experimental environments, we numerically confirm that loosely bound higher-order RWs with unusual quantities can be excited in a chaotic background field via the split-step Fourier method^{91,92} (see the Appendix for details). For this purpose, we put a low-amplitude harmonic waves noise on the plane-wave $ \psi 0(x,t)= q 0 e i \mu + \nu $ and $ q 0=c e i [ a \xi + ( c 2 \u2212 1 2 a 2 ) \tau ]$ [ $\xi $, $\tau $, $\mu $, $\nu $ have been listed in Eq. (9)] as the initial condition, which can be expressed as $ \psi i n i(x,t)= \psi 0(x,t)[1+ a 0cos\u2061( w 0x)]$, where $ a 0$ is intensity and $ w 0$ is frequency. Then consider the following two cases:

Under parameters as in Fig. 3(a7), i.e., $c=1$, $a=0$, $ c 0=\u221214$, $ c 1= c 2= c 3=0$, $ c 4= 2 3$, $ g 0=2$, $\lambda = 1 50$, this chaotic background ( $ a 0= 10 \u2212 2$ and $ w 0=2$) then develops into a “sea” of abundant nonlinear waves, as plotted in Figs. 8(a1) and 8(b1), where the typical loosely bound higher-order RW structure (three first-order RWs) can be singled out at a certain area (red circles). For further comparison, Fig. 8(c1) plots the profiles of the numerical solution chosen in the red circle at $t=99.5$ and $t=104.5$ in the chaotic field, which coincides extremely with the exact solution.

Under parameters as in Fig. 3(a8), i.e., only changing the $ c 0=\u22122$ and $ c 4= 1 2$, this chaotic background develops into a “sea” of abundant nonlinear waves, as shown in Figs. 8(a2) and 8(b2). It is worth noting that the energy of the K–MB at $x=0$ in the previous case disappears, therefore an unusual loosely bound higher-order RW structure (two first-order RWs) can be singled out at a certain time ( $t=104$, red circle). Figure 8(c2) plots the profiles of the numerical solution chosen in the red circle at $t=104$ in the chaotic field, which coincides extremely with the exact solution.

^{9,10,46}The parameters in $ \psi 0(x,t)$ can manage the dynamics of localized waves in an established experimental setup.

## VI. CONCLUSIONS

This paper investigated an integrable GPE. A self-similarity transformation was constructed to establish a connection between the GPE and the NLSE. The localized wave solutions were obtained through the combination of the DT and the self-similarity transformation. The self-similarity transformation described in this research introduced new control and dynamical properties for localized waves:

We plotted nonautonomous breathers, including the K–MB, the AB, and two-K–MB. The period of peak of the K–MB gradually decreased during propagation. Numerical simulations showed that the K–MB of the GPE had strong dynamical stability, and the longer period had stronger stability than the smaller period. Decreasing $\lambda $ increased its stability. The AB had weak dynamical stability. The breathers could be modulated by the parameters $ c i$, $i=0,1,3,4$.

The RW showed that the fundamental first-order RW could be controlled without changing the nonlinear strength (or scattering length). The proposed self-similarity transformation allowed for the compression and amplification of tightly and loosely bound higher-order RW, which supported the idea that there could be compression in RW. We could increase the stability of RW by prolonging its lifespan controlled by the parameters $ c i$, $i=0,1,3,4$ in the self-similarity transformation.

We found that $ c i$ would change the number of loose $N$th-order RWs. The asymptotic analysis of RW suggested that this reduction in quantity was an energy transfer. We could also stabilize the existence of loosely bound higher-order RWs by reducing the number of loosely bound RWs. We demonstrated numerically that this unusual higher-order RW could be excited from a chaotic background field. It supported the feasibility of experimental observation.

It is predicted that the application of the self-similarity transformation will yield more novel nonlinear wave dynamics in other nonlinear models, including but not limited to the nonlocal coupled GPE, high-dimensional GPE, etc. Our future investigation will prioritize the examination of experimental observations pertaining to nonlinear waves based on some nonlinear models. In addition, we will continue in the exploration of further nonlinear phenomena inside these nonlinear integrable models.

## ACKNOWLEDGMENTS

This work was supported by the Beijing Major and Special Project (Grant No. Z231100006623006), the National Key Research and Development Program of China (Grant No. 2022YFA1604202), the Beijing Natural Science Foundation (Grant No. JQ21019), the National Natural Science Foundation of China (Grant Nos. 11975001, 12075034, and 12261131495), the Hebei Key Laboratory of Physics and Energy Technology (Grant No. HBKLPET2023 03), and the BUPT Excellent Ph.D. Students Foundation (Grant No. CX2022238).

## AUTHOR DECLARATIONS

### Conflict of Interest

The authors have no conflicts to disclose.

### Author Contributions

**Haotian Wang:** Investigation (equal); Methodology (equal); Writing – original draft (equal). **Hujiang Yang:** Investigation (equal); Writing – review & editing (equal). **Ye Tian:** Investigation (equal); Writing – review & editing (equal). **Wenjun Liu:** Funding acquisition (equal); Writing – review & editing (equal).

## DATA AVAILABILITY

The data that support the findings of this study are available within the article.

### APPENDIX: THE DETAILED PROCESS OF SPLIT-STEP FOURIER METHOD

Following the above steps, we can use the split-step Fourier method to explore the RW excitation from the initial perturbed plane wave, whose excitation mechanism is the existence of modulation instability.^{91,92} Therefore, when the perturbed plane wave with the same parameters of the exact solution is as the initial condition, this exact solution can be observed in the numerical evolution results.

## REFERENCES

*et al.*,

*Emergent Nonlinear Phenomena in Bose-Einstein Condensates*

*et al.*, “

*et al.*, “

*Rogue Waves in the Ocean*

*et al.*, “

*et al.*, “

*et al.*, “

*Bose-Einstein Condensation and Superfluidity*

*Fundamentals and New Frontiers of Bose-Einstein Condensation*

*Solitons in Optical Communications*

*Nonlinear Physics of Plasmas*

*et al.*, “

*et al.*, “

*Nonlinear Waves in Integrable and Nonintegrable Systems*