Modulation Instability and Associated Dissipative Rogue Waves and Breathers in Ultra-Cold Neutral Plasma With Adiabatic Trapped Electron and Warm Ions

Abstract

Nonlinear envelope dynamics of low-frequency ion-acoustic disturbances are examined in an ultra-cold, charge-neutral electron-ion plasma, in which the ions are modeled as a warm viscous fluid, and the electrons follow an adiabatically trapped degenerate equation of state. The study is devoted to clarifying how ion kinematic viscosity, together with electron degeneracy effects, control the onset and development of modulational instability and the emergence of strongly localized ion-acoustic envelopes, including breather-like configurations and rogue waves, in such ultra-cold environments. Starting from a fluid description of the plasma and adopting an adiabatic trapped model for the degenerate electrons, a dissipative nonlinear Schrödinger (DNLS) equation governing the slow evolution of the ion-acoustic carrier envelope is derived through a reductive perturbation technique (the derivative expansion method). Within this DNLS framework, the interplay between dispersion, cubic nonlinearity, and dissipation is analyzed, and the modulational instability characteristics are obtained by tracking the influence of the degenerate electron temperature T and the ion kinematic viscosity η₀ on the instability growth rate and associated stability domains of the carrier wave. To complement the analytical stability analysis, the DNLS equation is integrated numerically by means of a second-order split-step Fourier algorithm, which allows us to follow the spatiotemporal evolution of rational breather-type excitations on a finite background, including Akhmediev and Kuznetsov-Ma breathers, and the Peregrine soliton. The numerical results show that stronger ion kinematic viscosity substantially delays the development of modulationally unstable patterns and suppresses their peak amplitudes, whereas an increase in T enhances the effective nonlinearity and facilitates modulational growth. The overall picture that emerges provides a consistent description of dissipative breather and rogue wave dynamics in ultra-cold neutral plasmas and may be of relevance for interpreting localized electrostatic activity in space and astrophysical plasmas, ultra-dense laboratory configurations, and other high-energy-density plasma settings.

1. Introduction

Ion-acoustic waves (IAWs) and their nonlinear manifestations play a prominent role in the dynamics of many space and laboratory plasmas. Over the past decades, ion–acoustic solitary waves (IASWs) have been reported in a broad range of configurations and have been the subject of numerous experimental and theoretical investigations.^1–4 In these structures, the inertia is provided by the ion component, whereas the restoring mechanism is essentially associated with the electron pressure. Deviations from simple Maxwellian electrons, together with finite ion temperature and additional microscopic processes, are now known to substantially modify the existence conditions, amplitude, and polarity of ion-acoustic solitons (IASs). For instance, the combined influence of warm adiabatic ions, nonthermal electrons, and external magnetic fields on IASs in magnetized plasmas has been examined in several works, showing that nonthermal electron populations can significantly alter the shape and admissible parameter domains of IASWs.^5,6

Electron trapping and degeneracy introduce further refinements in the description of IAWs by incorporating explicit wave–particle interaction mechanisms and quantum–statistical effects. ⁷ The impact of adiabatic trapping at the microscopic level has been addressed in earlier analytical studies and is supported by experimental observations and numerical simulations, which confirm the formation of trapped particle populations around nonlinear electrostatic structures.^8–10 More recently, the role of trapping in quantum plasmas has been analyzed by Shah et al. ¹¹ following the approach proposed in Ref. 12, where IASW structures were obtained in partially and fully degenerate plasmas with small temperature corrections. These contributions highlight that the interplay among ion temperature, electron trapping, and degeneracy can strongly influence IAW propagation and the parameter regimes in which coherent nonlinear structures may form. Beyond trapped distributions, the Kaniadakis (κ) distribution has been increasingly used to model electron populations in astrophysical and laboratory plasmas. Recently, a study on dust-ion-acoustic solitary waves in an unmagnetized plasma with Kaniadakis-distributed electrons ¹³ demonstrated that the spectral index κ and the electron-to-ion temperature ratio play critical roles in determining the existence domains, amplitude, and polarity of solitary structures. While the present work focuses on adiabatically trapped degenerate electrons, the inclusion of κ-distributions represents a promising extension for future investigations.

On a larger scale, the evolution of weakly nonlinear, weakly dispersive electrostatic wave packets is controlled by slow modulation processes, which are conveniently described by envelope equations of the nonlinear Schrödinger (NLS)-type. Within this framework, modulational instability (MI) of a nearly monochromatic carrier wave provides a robust route to concentrating wave energy and to the emergence of localized envelopes, such as bright and dark solitons, breathers, and rogue waves (RWs). The MI of various electrostatic modes has been investigated in a broad range of plasma models, including multicomponent electron–ion and pair plasmas with Maxwellian or non–Maxwellian distributions, by employing reductive perturbation or multiple–scales techniques to derive the standard NLS equation.^14–16 For IAWs, MI has been analyzed, for example, in double–pair plasmas with adiabatic ions and superthermal electrons and positrons, ¹⁷ as well as in systems containing warm adiabatic ions, isothermal positrons, and two–temperature superthermal electrons.^18,19 More recently, amplitude modulation and soliton formation in dense quantum plasmas have been investigated using NLS-type frameworks, ²⁰ and field modulations of IAWs with Vasyliunas–Schamel distributed electrons have been studied in the context of MI. ²¹ These works have established that MI is a pervasive mechanism leading to energy localization in plasmas and that the stability of the carrier wave is highly sensitive to the details of the particle distributions and temperature ratios. In parallel to NLS-based envelope descriptions, much of the literature on nonlinear dispersive waves has focused on Korteweg-de Vries (KdV)-type equations, which govern weakly nonlinear long waves. For instance, the optimal homotopy asymptotic method (OHAM) has been successfully applied to the generalized KdV equation, ²² providing highly accurate approximate solutions for solitary wave profiles without requiring small parameters. Such techniques are especially valuable when exact integrability is absent or when dissipative or forcing terms break the Hamiltonian structure. While our current study employs a dissipative NLS (DNLS) framework to describe modulated envelope waves rather than KdV-type solitary waves, the OHAM approach could, in principle, be adapted to solve the stationary forms of our DNLS equation or its higher-order extensions, particularly when seeking approximate dissipative soliton solutions.

Once the NLS description is available, a wealth of exact and approximate solutions can be used to characterize the nonlinear stages of MI. Among the most prominent are the rational breather families of Akhmediev, Kuznetsov–Ma, and Peregrine (rogue wave (RW)), which have been widely discussed in the context of optics and hydrodynamics and transplanted to plasma physics.^23–28 The Akhmediev breather (AB) is localized in space and periodic in time, representing a transient amplification on a finite background; the Kuznetsov–Ma (KM) breather is localized in time and periodic in space, while the Peregrine soliton (PS) is doubly localized and often invoked as an archetype of RWs. Numerical and analytical studies indicate that such structures may appear ”from nowhere” and disappear without a trace, thus providing a plausible framework for understanding extreme events in dispersive nonlinear media.^23,26–28 Recent work has extended these concepts to include breather interactions in electron–ion plasmas with electron beams, ²⁹ electron acoustic Peregrine breathers in quantum plasmas with temperature anisotropy, ³⁰ and the evolution of electrostatic shocks into RW structures. ³¹ RW generation through nonlinear self-interaction of electrostatic waves in dense plasma has also been reported, ³² along with studies of breather structures in the framework of the Gardner equation in electron-positron-ion plasmas. ³³ In addition, dissipative extensions of the NLS equation have been used to explore how collisional or viscous effects alter the onset and growth of MI and the associated envelope structures in electronegative and nonthermal plasmas.^34–37 For instance, it has been shown that collisional dissipation can deform and attenuate breathers and rogue waves, while still allowing for transient energy concentration under appropriate conditions.^34,37 Chaotic excitations of RWs in pair plasmas have also been reported, ³⁸ suggesting rich dynamics beyond the regular breather regime. Furthermore, the propagation of RWs and cnoidal wave formations through low-frequency plasma oscillations has been investigated, ³⁹ highlighting the diverse manifestations of extreme wave events in plasma systems.

Despite this substantial body of work, the modulational properties of IAW packets in ultra–cold neutral plasmas with adiabatically trapped, degenerate electrons and warm viscous ions have not, to the best of our knowledge, been explored through a DNLS framework. In particular, a systematic treatment that connects the microscopic plasma parameters (degenerate electron temperature and ion kinematic viscosity) to the coefficients of a DNLS equation and uses this link to analyze MI and compute rational breather solutions appears to be missing. The present contribution aims to fill this gap.

In this work, we consider a homogeneous, ultra–cold electron–ion plasma, in which the ions form a warm, viscous, inertial fluid and the electrons are modeled as adiabatically trapped and degenerate.^40–46 The electron number density is described by a trapped degenerate distribution, whose expansion yields a set of nonlinear coefficients entering Poisson’s equation.^45,46 Using the reductive perturbation technique, we derive, from the normalized fluid and Poisson equations, a DNLS equation for the slowly varying envelope of IAWs, where the dispersion, nonlinearity, and dissipation coefficients explicitly depend on the ion temperature ratio, the degeneracy parameter, and the ion kinematic viscosity. The resulting DNLS equation is then employed to perform a detailed MI analysis and to compute dissipative rational breather and rogue wave solutions.

The paper is organized as follows. In Sec. 2, we present the normalized fluid model with adiabatically trapped degenerate electrons and warm viscous ions and derive the DNLS equation governing the ion–acoustic envelope. Section 3 is devoted to the analytical MI analysis based on the DNLS equation, where the stability criteria and growth rates are obtained as functions of the relevant plasma parameters. In Sec. 4, we perform numerical simulations of rational breather and RW solutions of the DNLS equation, using a second–order split–step Fourier method, and discuss the role of the ion kinematic viscosity and degenerate electron temperature on their spatiotemporal evolution. A concise summary of the main findings and some perspectives for future work are provided in the final section.

2. Governing Equations and Derivation of the dissipative NLS Equation

Here, we consider a one–dimensional homogeneous electron-ion plasma composed of two species, including warm inertial positive non-degenerate ions and temperature adiabatic trapped degenerate electrons. In this model, the ion kinematic viscosity is considered. Within this framework, the normalized dynamical fluid equations describing the propagation of the IAWs in this model are introduced^40–44:

$\frac{\partial }{\partial t}{n}_i+\frac{\partial }{\partial x}\left(n_i{u}_i\right)=0,$

(1)

$\frac{\partial }{\partial t}{u}_i+u_i\frac{\partial }{\partial x}{u}_i+\frac{\partial }{\partial x}\varphi +\frac{\sigma }{n_i}\frac{\partial }{\partial x}{P}_i-\frac{\eta }{n_i}\frac{{\partial }^2}{\partial x^2}{u}_i=0,$

(2)

$\frac{\partial }{\partial t}{P}_i+u_i\frac{\partial }{\partial x}{P}_i+3P_i\frac{\partial }{\partial x}{u}_i=0,$

(3)

$\frac{{\partial }^2}{\partial x^2}\varphi =n_e-n_i.$

(4)

Here, the quantities n _i and and n _e indicate the normalized number densities of both ion and electron (n _i = N _i /n_i0, n _e = N _e /n_e0), respectively, with n_i0 (n_e0) being the ion (electron) equilibrium density. The term u _i represents the re-scaled ions fluid speed (u _i = U _i /C _s ), and the normalized electrostatic potential is denoted by ϕ. The parameters σ and η correspond to the ratio of ion to electron temperature and the ion kinematic viscosity, respectively. The normalized adiabatic trapped degenerate electron number density reads^45,46:

$n_e={\left(1+\varphi \right)}^{3/2}+T^2{\left(1+\varphi \right)}^{-1/2},$

(5)

where T denotes the normalized temperature of the degenerate electrons defined as $T=\pi T_e/2\sqrt{2}{E}_{Fe}$ , where T _e is the actual temperature of the degenerate electron gas, and $E_{Fe}={(3{\pi }^2{n}_{e0})}^{2/3}({\hslash }^2/2m_e)$ is the electron Fermi energy. ϕ is the corresponding electrostatic potential normalized by eΦ/E _Fe . Equation (5) follows from the trapped degenerate electron model derived by Shah et al. ¹¹ The first term in equation (5) is associated with the contribution of trapped electrons, whereas the second term reflects the influence of finite temperature in a partially degenerate electron population. ⁴⁶ In what follows, we consider parameter values representative of ultra-dense astrophysical environments, for which relativistic degeneracy is important, ⁴⁶ namely n_e0 = 10²⁷cm⁻³ and σ = T _i /T _e = 0.03. So, the Poisson equation can be expressed as follows according to the Taylor expansion of equation (5):

$\begin{array}{c}\frac{{\partial }^2}{\partial x^2}\varphi +n_i-(1+T^2+{\gamma }_1\varphi +{\gamma }_2{\varphi }^2+{\gamma }_3{\varphi }^3)=0,\hfill \\ {\gamma }_1=\frac{3-T^2}{2},{\gamma }_2=\frac{3+3T^2}{8},{\gamma }_3=\frac{-1-5T^2}{16}.\hfill \end{array}$

(6)

Now, to investigate the slow modulation of a weakly nonlinear ion–acoustic carrier, we employ a reductive perturbation (the derivative expansion) technique and introduce stretched independent variables ( $\zeta =ϵ\left(x-v_{g}t\right)$ and τ = ϵ²t, where v _g is the group velocity) corresponding to the envelope scales. The dependent variables are expanded around their equilibrium values in a power series of a smallness parameter, and the resulting hierarchy of equations is solved order by order. At the lowest order, one recovers the linear dispersion relation of the ion–acoustic mode. At higher orders, the compatibility conditions at the secular level yield an evolution equation for the complex envelope of the carrier wave. After a standard multiple–scale reduction, this procedure yields a dissipative nonlinear Schrödinger (DNLS) equation for the slowly varying ion–acoustic envelope, where the dispersive, nonlinear, and dissipative coefficients are expressed in terms of the underlying plasma parameters, including the ion temperature ratio, the degeneracy parameter T, and the ion kinematic viscosity η.

Accordingly, the following expansions are introduced:

$\begin{array}{ll}\hfill \left[\begin{array}{c}\hfill n_i\left(x,t\right)\hfill \\ \hfill u_i\left(x,t\right)\hfill \\ \hfill P_i\left(x,t\right)\hfill \\ \hfill \varphi \left(x,t\right)\hfill \end{array}\right]& =\left[\begin{array}{c}\hfill 1\hfill \\ \hfill 0\hfill \\ \hfill 1\hfill \\ \hfill 0\hfill \end{array}\right]+\sum _{m=1}^{\infty }{ϵ}^m\sum _{l=-m}^{\infty }\left[\begin{array}{c}\hfill {n_l}^{(m)}\left(\zeta ,\tau \right)\hfill \\ \hfill {u_l}^{(m)}\left(\zeta ,\tau \right)\hfill \\ \hfill {P_l}^{(m)}\left(\zeta ,\tau \right)\hfill \\ \hfill {{\varphi }_l}^{(m)}\left(\zeta ,\tau \right)\hfill \\ \hfill \hfill \end{array}\right]A{\left(x,t\right)}^l.\hfill \end{array}$

(7)

We assume that the above series include all overtones $A{\left(x,t\right)}^l={\mathrm{e}}^{il\left(kx-\omega t\right)}$ up to order m. These are due to nonlinear terms with the corresponding coefficients being of maximum power ϵ ^m along with the relations ${({n_l}^{(m)})}^{*}={n_{-l}}^{(m)}$ , ${({u_l}^{(m)})}^{*}={u_{-l}}^{(m)}$ , ${({P_1}^{(m)})}^{*}={P_{-l}}^{(m)}$ , and ${({{\varphi }_l}^{(m)})}^{*}={{\varphi }_{-l}}^{(m)}$ . We assumed η = ϵ² η₀. By substituting equation (7), into (1), (2), (3), and (4), the first-order harmonics are obtained in the form

$\begin{array}{c}{n_1}^{(1)}=\left(k^2+{\gamma }_1\right){{\varphi }_1}^{(1)},{u_1}^{(1)}=-\frac{{{\varphi }_1}^{(1)}k\omega }{3k^2\sigma -{\omega }^2},{P_1}^{(1)}=-3\frac{{{\varphi }_1}^{(1)}{k}^2}{3k^2\sigma -{\omega }^2}.\hfill \end{array}$

(8)

which give the dispersion relation as:

${\omega }^2=\frac{k^2}{k^2+{\gamma }_1}+3k^2\sigma .$

(9)

The dispersion relation obtained is a function of the ratio of ion temperature σ, so when (σ = 0), the same dispersion relation is found in Ref. 47. For the second order, with m = 2 and l = 0, we obtain

$\begin{array}{c}{n_0}^2=B_5{\left|{{\varphi }_1}^{(1)}\right|}^2,{u_0}^2=B_6{\left|{{\varphi }_1}^{(1)}\right|}^2,\hfill \\ {{\varphi }_0}^2=B_8{\left|{{\varphi }_1}^{(1)}\right|}^2,{P_0}^2=B_7{\left|{{\varphi }_1}^{(1)}\right|}^2,\hfill \end{array}$

(10)

with

$\begin{array}{c}{B}_5=B_8{\gamma }_1+2{\gamma }_2,\hfill \\ B_6=-\frac{-9k^6{\sigma }^2{v}_g+6k^5\omega \sigma {v_g}^2+\left(-3{\omega }^2\sigma +9{\sigma }^2{B}_8+9{\sigma }^2{\gamma }_1\right)v_g{k}^4+C_4}{\left(9k^4{\sigma }^2-6k^2{\omega }^2\sigma +{\omega }^4\right)\left(-{v_g}^2+3\sigma \right)},\hfill \\ B_7=\frac{27k^6{\sigma }^2-18k^5\omega \sigma v_g+\left(9{\omega }^2\sigma -27{\sigma }^2{B}_8-27{\sigma }^2{\gamma }_1\right)k^4+C_5}{\left(-{v_g}^2+3\sigma \right){\left(3k^2\sigma -{\omega }^2\right)}^2},\hfill \\ B_8=\frac{\left(-2{\omega }^3{v_g}^2+6{\omega }^3\sigma -18\omega {\sigma }^2{\gamma }_1\right)k^3+\left(-2{\omega }^3{\gamma }_1{v_g}^2+6{\omega }^3\sigma {\gamma }_1\right)k-6{\omega }^4\sigma {\gamma }_2{v}_g+C_6+C_7}{v_g{\left(3k^2\sigma -{\omega }^2\right)}^2\left(-{\gamma }_1{v_g}^2+3\sigma {\gamma }_1+1\right)},\hfill \end{array}$

(11)

where

$\begin{array}{c}{C}_4=6k^3\omega \sigma {\gamma }_1{v_g}^2+\left(\left(-6{\omega }^2\sigma B_8-9{\omega }^2\sigma {\gamma }_1+{\omega }^2\right)v_g+{\omega }^2\sigma \right)k^2+{\omega }^4{B}_8{v}_g,\hfill \\ C_5=-18k^3\omega \sigma {\gamma }_1{v}_g+\left(18{\omega }^2\sigma B_8+27{\omega }^2\sigma {\gamma }_1-{\omega }^2{v}_g-3{\omega }^2\right)k^2-3{\omega }^4{B}_8,\hfill \\ C_6=\left(-12{\omega }^2\sigma {\gamma }_2{v_g}^3+\left(36{\omega }^2{\sigma }^2{\gamma }_2+9{\omega }^2\sigma {\gamma }_1-{\omega }^2\right)v_g-{\omega }^2\sigma \right)k^2+2{\omega }^4{\gamma }_2{v_g}^3,\hfill \\ C_7=9k^6{\sigma }^2{v}_g-18k^5\omega {\sigma }^2+\left(18{\sigma }^2{\gamma }_2{v_g}^3+\left(-54{\sigma }^3{\gamma }_2+3{\omega }^2\sigma -9{\sigma }^2{\gamma }_1\right)v_g\right)k^4.\hfill \end{array}$

(12)

At the same order, i.e, m = 2 for l = 1, solutions exist under the compatibility condition

$\begin{array}{c}{v}_g=-2\frac{k\omega \left(9k^4{\sigma }^2-6k^2{\omega }^2\sigma +{\omega }^4-{\omega }^2\right)}{9k^6{\sigma }^2+\left(-6{\omega }^2\sigma +9{\sigma }^2{\gamma }_1+3\sigma \right)k^4+R},\hfill \\ R={\omega }^2\left({\omega }^2-6\sigma {\gamma }_1+1\right)k^2+{\omega }^4{\gamma }_1.\hfill \end{array}$

(13)

The relations (13) represent the group velocity of IAWs, which strongly depends on the ratio of ion temperature σ. For m = 2 and l = 2, we find the second harmonic quantities ${{\varphi }_2}^{(2)}$ , ${n_2}^{(2)}$ , ${u_2}^{(2)}$ , ${v_2}^{(2)}$ and ${w_2}^{(2)}$ , in term of ${{\varphi }_1}^{(1)}$ in the form

$\begin{array}{c}{n_2}^{(2)}=B_1{\left|{{\varphi }_1}^{(1)}\right|}^2,{u_2}^{(2)}=B_2{\left|{{\varphi }_1}^{(1)}\right|}^2,\hfill \\ {P_2}^{(2)}=B_3{\left|{{\varphi }_1}^{(1)}\right|}^2,{{\varphi }_2}^{(2)}=B_4{\left|{{\varphi }_1}^{(1)}\right|}^2.\hfill \end{array}$

(14)

with

$\begin{array}{c}{B}_1=\left(4k^2+{\gamma }_1\right)B_4+{\gamma }_2,\hfill \\ B_2=-\frac{\left(\left(18k^4{\sigma }^2-12k^2{\omega }^2\sigma +2{\omega }^4\right)B_4+9k^6{\sigma }^2+C_3\right)k\omega }{2\left(9k^4{\sigma }^2-6k^2{\omega }^2\sigma +{\omega }^4\right)\left(3k^2\sigma -{\omega }^2\right)},\hfill \\ B_3=-\frac{k^2\left(27k^6{\sigma }^2-9\sigma \left({\omega }^2-6\sigma \left(B_4+1/2{\gamma }_1\right)\right)k^4+C_2\right)}{2{\left(3k^2\sigma -{\omega }^2\right)}^3},\hfill \\ B_4=-\frac{27k^8{\sigma }^2+\left(54{\sigma }^3{\gamma }_2-15{\omega }^2\sigma +27{\sigma }^2{\gamma }_1\right)k^6+k^5\omega \sigma +C_1}{2{\left(3k^2\sigma -{\omega }^2\right)}^2\left(\left(-4k^2-{\gamma }_1\right){\omega }^2+k^2\left(12k^2\sigma +3\sigma {\gamma }_1+1\right)\right)},\hfill \end{array}$

(15)

where

$\begin{array}{c}{C}_1=\left(-54{\omega }^2{\sigma }^2{\gamma }_2+2{\omega }^4-15{\omega }^2\sigma {\gamma }_1+{\omega }^2\right)k^4+\left(18{\omega }^4\sigma {\gamma }_2+2{\omega }^4{\gamma }_1\right)k^2-2{\omega }^6{\gamma }_2,\hfill \\ C_2=-\left(-3+36\left(B_4+1/4{\gamma }_1\right)\sigma \right){\omega }^2{k}^2+k{\omega }^3+6{\omega }^4{B}_4,\hfill \\ C_3=\left(-3{\omega }^2\sigma +9{\sigma }^2{\gamma }_1\right)k^4+k^3\omega \sigma +\left(-3{\omega }^2\sigma {\gamma }_1+{\omega }^2\right)k^2.\hfill \end{array}$

(16)

Finally, by retaining the third-harmonic contributions corresponding to m = 3 and l = 1, the algebraic system obtained from the preceding orders can be combined and simplified, which leads to the following DNLS equation for the modulated IAWs:

$\begin{array}{c}i\frac{\partial }{\partial \tau }\varphi +P\frac{{\partial }^2}{\partial {\zeta }^2}\varphi +Q{\left|\varphi \right|}^2\varphi +i\gamma \varphi =0,\hfill \end{array}$

(17)

where the dispersion coefficient P and the nonlinear coefficient Q are given by:

$\begin{array}{c}P=\left[{\omega }^7+2k{\omega }^6{v}_g+\left(-9k^2\sigma -1\right){\omega }^5+\left(-18k^3\sigma v_g+2k{v}_g\right){\omega }^4+\right.\hfill \\ \left(27{\sigma }^2{k}^4-k^2{v_g}^2-9k^2\sigma \right){\omega }^3+\left(54k^5{\sigma }^2{v}_g+18k^3\sigma v_g\right){\omega }^2+\hfill \\ \left.\left(-27k^6{\sigma }^3-9k^4\sigma {v_g}^2\right)\omega -54k^7{\sigma }^3{v}_g\right]\times \frac{1}{\mathrm{\Upsilon }}.\hfill \\ Q=\left[\left(-2B_4{\gamma }_2-2B_8{\gamma }_2-3{\gamma }_3\right){\omega }^7+\left(\left(-k^3-k{\gamma }_1\right)B_2+\left(k^3+k{\gamma }_1\right)B_6\right){\omega }^6+\right.\hfill \\ \left(18k^2\sigma {\gamma }_2{B}_8+\left(2k^4+18k^2\sigma {\gamma }_2+2k^2{\gamma }_1\right)B_4+k^2{B}_5+27k^2\sigma {\gamma }_3+k^2{B}_1\right){\omega }^5+\hfill \\ \left(\left(3k^5\sigma +3k^3\sigma {\gamma }_1+k^3\right)B_2+\left(-9k^5\sigma -9k^3\sigma {\gamma }_1+k^3\right)B_6\right){\omega }^4+\hfill \\ \left(\left(-12k^6\sigma -54k^4{\sigma }^2{\gamma }_2-12k^4\sigma {\gamma }_1\right)B_4+C_8\right){\omega }^3+\hfill \\ \left(\left(9k^7{\sigma }^2+9k^5{\sigma }^2{\gamma }_1-3k^5\sigma \right)B_2+\left(27k^7{\sigma }^2+27k^5{\sigma }^2{\gamma }_1-3k^5\sigma \right)B_6\right){\omega }^2+\hfill \\ \left(\left(18k^8{\sigma }^2+54k^6{\sigma }^3{\gamma }_2+18k^6{\sigma }^2{\gamma }_1\right)B_4+C_9\right)\omega +\hfill \\ \left.\left(-27k^9{\sigma }^3-27k^7{\sigma }^3{\gamma }_1\right)B_2+\left(-27k^9{\sigma }^3-27k^7{\sigma }^3{\gamma }_1\right)B_6\right]\times \frac{1}{\mathrm{\Upsilon }}.\hfill \\ \gamma =\left[3k^6{\omega }^2\sigma -k^4{\omega }^4\right]\times \frac{{\eta }_0}{\mathrm{\Upsilon }},\hfill \end{array}$

with

$\begin{array}{c}{C}_8=k^6+\left(\left(-54B_8{\gamma }_2-81{\gamma }_3\right){\sigma }^2+\left(-3B_1+B_2-9B_5+B_6\right)\sigma +{\gamma }_1\right)k^4,\hfill \\ C_9=54k^6{\sigma }^3{\gamma }_2{B}_8+81k^6{\sigma }^3{\gamma }_3-3k^6{\sigma }^2{B}_2+18k^6{\sigma }^2{B}_5-3k^6{\sigma }^2{B}_6,\hfill \end{array}$

(18)

and $\mathrm{\Upsilon }=\left(3k^2\sigma -{\omega }^2\right)9k^6{\sigma }^2-6k^4{\omega }^2\sigma +9k^4{\sigma }^2{\gamma }_1+k^2{\omega }^4-6k^2{\omega }^2\sigma {\gamma }_1+3k^4\sigma +{\omega }^4{\gamma }_1+k^2{\omega }^2$ .

It is shown that the dispersion and nonlinear coefficients P and Q, respectively, are real functions and independent of the ion kinematic viscosity η₀. The term +iγϕ in equation (17) represents linear damping (energy loss), where the positive coefficient (γ > 0) is proportional to the ion kinematic viscosity η₀. If the coefficient γ vanishes, the DNLS equation is reduced to the standard NLS equation found in Ref. 47.

3. Stability analysis of IAWs

Once the DNLS equation governing the modulated envelope IAWs has been obtained, the modulational properties of a nearly monochromatic carrier can be examined by means of a standard stability analysis. In particular, one considers a continuous–wave (CW) background solution of the DNLS equation and superposes a small modulation with a given perturbation wave number. Linearizing with respect to the modulation amplitude yields a dispersion relation for the perturbations, from which growth rates and stability domains can be determined as functions of the system parameters, including T and the effective viscous coefficient η₀. For this purpose the perturbed amplitude ϕ can be written as ⁴⁸ :

$\varphi =\left({\varphi }_0+\delta \varphi \right){\mathrm{e}}^{-i{\int }_0^{\tau }g\left(t\right)\mathrm{d}t-\gamma \tau },$

(19)

where g denotes the shift of nonlinear frequency, ϕ₀ (real constant) represents the pump carrier wave amplitude, and δϕ indicates the small perturbation of this amplitude. The substitution of equation (19) into equation (17) and collecting both zeroth and first order terms, we ultimately get:

$g\left(\tau \right)=-Q{\left|{\varphi }_0\right|}^2{\mathrm{e}}^{-2\gamma \tau },$

(20)

$i\frac{\partial }{\partial \tau }\delta \varphi +P\frac{{\partial }^2}{\partial {\xi }^2}\delta \varphi +(\delta \varphi +\delta {\varphi }^{*})Q{\left|{\varphi }_0\right|}^2{\mathrm{e}}^{-2\gamma \tau }=0.$

(21)

We make use of the transformation δϕ = c + id, with

$\begin{array}{c}c=c_0{\mathrm{e}}^{iK\zeta -i\mathrm{\Omega }\tau }+cc,d=d_0{\mathrm{e}}^{iK\zeta -i\mathrm{\Omega }\tau }+cc.\hfill \end{array}$

(22)

Upon inserting these expressions into equation (17) and performing a linearization about the underlying plane–wave state, one arrives at the nonlinear dispersion relation previously derived by Kengne et al. ³⁷ :

${\mathrm{\Omega }}^2=K^4{P}^2\left(1-\frac{{K_c}^2(\tau )}{K^2}\right).$

(23)

Here, K is the modulation wavenumber and K _c is the time-dependent critical wavenumber delimiting the MI band:

${K_c}^2(\tau )=2\frac{Q{\left|{\varphi }_0\right|}^2{\mathrm{e}}^{-2\gamma \tau }}{P}.$

(24)

The local-instability growth rate of the wave packet whose amplitude follows equation (17) basically depends on the sign of the ratio P/Q.

Equation (23) shows that the carrier waves become stable when P/Q < 0 and $K^2>K_c^2(\tau )$ , for which Ω is always real and the modulated wave does not grow. Whereas, for P/Q > 0 and $K^2<K_c^2(\tau )$ , the carrier waves become modulationally unstable and Ω becomes imaginary and the instability is a purely growing mode, which is materialized by the local-instability growth rate ³⁷ :

$\mathrm{I}\mathrm{m}(\mathrm{\Omega })=|P|K^2\sqrt{\frac{K_c^2(\tau )}{K^2}-1}.$

(25)

We can observe from relation (25) that, there is a limited time to propagate the modulated waves at τ = τ_max, where the evolution of the wave breaks and the wave frequency Ω = 0, whereas at τ > τ_max, Ω² > 0, then the wave becomes modulationally stable and the modulation instability period is given as^37,48:

$\tau \ge {\tau }_{\mathrm{m}\mathrm{a}\mathrm{x}}=\frac{1}{2\gamma }\mathrm{l}\mathrm{n}\left(2\frac{Q{\left|{\varphi }_0\right|}^2}{{PK}^2}\right).$

(26)

The total growth rate of the modulation during the unstable period reads^48,49:

${\mathrm{\Gamma }}_g={\int }_0^{{\tau }_{\max }}Im\mathrm{\Omega }\left({\tau }^{\prime }\right)\mathrm{d}{\tau }^{\prime },$

(27)

${\mathrm{\Gamma }}_g=\frac{{|P|K}^2}{\gamma }\left(\sqrt{2\frac{Q{\left|{\varphi }_0\right|}^2}{{PK}^2}-1}-\arctan \sqrt{2\frac{Q{\left|{\varphi }_0\right|}^2}{{PK}^2}-1}\right),$

(28)

where the total growth rate is seen to reach its maximum value for $K=Q{\left|{\varphi }_0\right|}^2/P$ given by ${\mathrm{\Gamma }}_g(max)=\left(1-\pi /4\right)Q{\left|{\varphi }_0\right|}^2/\gamma$ .

The evolution of a fundamental wave whose amplitude follows the DNLS equation depends on both P and Q, which are also dependent on the relevant plasma parameters. When P and Q have the same sign, the evolution of the IAWs amplitude is modulationally unstable in the presence of external perturbations. Otherwise, the IAWs are modulationally stable under external perturbations.

The dispersion relation (23) makes it clear that the ion temperature ratio σ and the degeneracy contribution encoded in γ₁ control the basic linear properties of the ion–acoustic mode. In particular, Figure 1(a) shows that the wave frequency ω grows with the carrier wavenumber k at small k and then approaches an almost saturated value as k increases further, whereas the group velocity v _g exhibits a maximum and then decreases with k, as seen in Figure 1(b). An increase in the degeneracy temperature T strengthens the electron pressure term in the dispersion relation, raising the frequency and modifying the group velocity curve; this behavior can be interpreted as an enhancement of the restoring force associated with the degenerate electron component.

Figure 1.

Variation of the wave frequency ω and the group velocity v _g against the wave number k in panels (a) and (b) respectively. The influence of the degeneracy temperature T is studied in both cases by setting σ = 0.01.

The influence of T and of the ion kinematic viscosity η₀ on modulational stability is summarized in the stability maps of Figure 2. In the (T, k)–plane [Figure 2(a)], one observes that higher values of T reduce the domain where the carrier is modulationally stable and enlarge the region where MI takes place, with the unstable band extending towards larger wavenumbers. This is consistent with the way T enters the expansion coefficients γ₁, γ₂, γ₃ and, through them, the nonlinear coefficient Q, thereby reinforcing the effective nonlinearity and favoring the onset of MI. By contrast, the (η₀, k)–diagram in Figure 2(b) shows that the low–k part of the spectrum remains stable for all considered values of η₀, whereas an unstable band appears only beyond a certain k and its extent is not substantially modified when η₀ is increased. Since the dissipative coefficient γ in the DNLS equation is proportional to η₀, larger viscosity mainly enhances damping and reduces the MI growth rate, without significantly shifting the range of carrier wavenumbers where instability is possible. This result agrees with the findings reported by El-Tantawy et al. ⁵⁰

Figure 2.

Diagram of stability/instability versus the degeneracy temperature T in panel (a), and the ion kinematic viscosity η₀ in panel (b). The other parameters are: σ = 0.01, η₀ = 0.01, and T = 0.5.

The effect of these parameters on the total growth rate Γ _g is displayed in Figure 3. For fixed k in the regime P/Q > 0, panel 3(a) reveals that the peak value of Γ _g increases with T, and the corresponding MI band in K becomes broader. This confirms that degeneracy promotes nonlinear focusing by amplifying the cubic nonlinearity and enlarging the range of modulationally unstable perturbations. In contrast, panel 3(b) shows that increasing η₀ reduces the maximum of Γ _g and narrows the effective MI band in K, meaning that viscous dissipation shortens the time window over which MI can act and limits the achievable amplification. These features are in full agreement with the qualitative expectations based on equations (25)–(28) and with earlier findings on dissipative NLS–type models.^19,37 However, inside the modulational unstable region, i.e. for P/Q > 0, it is possible for a contingent perturbation of the amplitude to grow and therefore lead to the creation of the solitonic structures. ⁵¹

Figure 3.

Growth rate of MI, versus the wavenumber K. The influence of the degeneracy temperature T in panel (a), and the ion kinematic viscosity η₀ in panel (b). The other parameters are: k = 2, ϕ₀ = 0.2, σ = 0.01, η₀ = 0.02, and T = 0.5.

4. Numerical simulation of rational breathers: solution of dissipative NLS equation

In general, linear stability analysis is limited because it can only detect the onset of instability and therefore doesn’t reveal what kinds of dynamical patterns one might observe in the system as instability persists over time.

The NLS equation without dissipation is known to be integrable and to exhibit a broad class of localized envelope structures, both stationary and pulsating. Among these are rational breathers for P/Q > 0 (MI region) with a continuous wave background. ²⁶ A generic solution of the derived NLS equation (17), (with η₀ = 0) with a continuous wave background can be compactly expressed as follows ⁵² :

$\varphi (\xi ,\tau )={\varphi }_0\left[1+\frac{2(1-2a)cosh(2P\beta \tau )+i\beta sinh(2P\beta \tau )}{\sqrt{2a}cos(\phi \xi )-cosh(2P\beta \tau )}\right]e^{\mathit{2iP\tau }}.$

(29)

A perturbation wave number is represented by the independent variable parameter φ, which is a family of solutions to this equation. Here, ϕ₀ and a represent the continuous wave field’s starting amplitude and breather parameter, respectively. Depending on these parameters (ϕ₀ and a), the independent variable φ and the growing and decaying rate β are determined by ⁵³ :

$\begin{array}{c}\beta =\sqrt{8a(1-2a)},{\phi }^2={\varpi }^2(1-2a),{\varpi }^2=2{\varphi }_0^2(Q/P).\hfill \end{array}$

(30)

For a = 0.25, or $\phi =\varpi /\sqrt{2}$ , the gain reaches its maximum at β = 1. Keep in mind that the solution indicates an AB when 0 < a < 0.5.

For a → 0.5, the gain becomes algebraic rather than exponential, the modulation frequency approaches infinity, and the PS solution describes the wave dynamics.^54,55

For a > 0.5, the KM breather characterizes the modulated wave dynamics. ⁵⁶

For numerical integration of the DNLS equation (17), we cannot directly use the analytical expression (29) as an initial condition when η₀ ≠ 0 (dissipative case), because (29) is an exact solution only for γ = 0. Instead, we use a perturbed continuous wave as the initial condition, which triggers the MI-induced growth of the desired breather pattern. The general form of the initial condition is:

$\varphi (\xi ,\tau =0)={\varphi }_0\left[1+\vartheta cos(\phi \xi )\right],$

(31)

where $\vartheta =\psi \left(1+i(2P\beta /{\phi }^2)\right)$ , is the modulation parameter and ψ being a real parameter.

We are going to use the well known second order split-step Fourier scheme ⁵⁷ with periodic boundary conditions to accomplish the numerical simulation with ξ ∈ [−L, L], L = 15, the number of spatial grid points N = 1000, the temporal step size Δτ = 0.01, and the spatial step size Δξ = 2 L/N.

The spatiotemporal evolution of AB for different values of η₀ is seen in Figure 4. This solution is periodic in time and localized in space. It is found that the behavior of the kinematic viscosity parameter (η₀) leads to a reduction of the nonzero continuous wave backgrounds (see panels 4(a)-(d)). This indicates that the distinctive features of the AB are influenced by the kinematic viscosity.

Figure 4.

Computation of AB. In panels (a) and (b) η₀ = 0, whereas in panels (c) and (d) η₀ = 0.09. The other parameters are: k = 0.4, ϕ₀ = 0.45, a = 0.4, ψ = 0.001, σ = 0.03, and T = 0.4.

An illustration of the evolution of dissipative KM breathers due to the influence of kinematic viscosity is shown in Figure 5. The results demonstrate how the dissipative KM breather exhibits a decrease in amplitude as it propagates. These waves are extremely sensitive to changes in the η₀ value, as evidenced by the fact that the maximum amplitude of the KM breather for η₀ = 0 is larger than that of η₀ = 0.1. As a result, the system naturally dissipates sufficient energy when η₀ grows, which causes the amplitudes of the dissipative KM breather to decrease.

Figure 5.

Computation of KM breather. In panels (a) and (b) η₀ = 0, whereas in panels (c) and (d) η₀ = 0.1. The other parameters are: k = 0.5, ϕ₀ = 0.3, a = 0.7, ψ = 0.5, σ = 0.03, and T = 0.4.

The PS has an endless period in both space and time since it is a doubly confined structure. Mouhammdoul et al. conducted a thorough investigation of the dynamics of RWs. They have explained how the nonthermal parameter and the magnetic field affect the RWs’ amplitude, but they have used an analytical solution of the DNLS equation. The dissipative PS evolution for different values of η₀ is shown in Figure 6. Figure 6(a) and (b) are obtained for η₀ = 0, and exhibit RW with higher amplitude. It is observed that the PS has concentrated the energy of the IAWs in a small area, as shown in these panels. In contrast to the traditional PS, the kinematic viscosity parameter η₀ causes the backgrounds of the dissipative PS to diminish.

Figure 6.

Computation of Peregrine soliton. In panels (a) and (b) η₀ = 0, whereas in panels (c) and (d) η₀ = 0.8. The other parameters are: k = 0.05, ϕ₀ = 0.15, a = 0.505, ψ = 0.9, σ = 0.03, and T = 0.4.

All these phenomena are due to the manifestation of MI in the system. Generally, this happens in the same domain of parameters where localized structures such as bright solitons, breathers, and rogue waves are found, suggesting that MI is a progenitor of soliton formation, as extensively discussed across a broad range of physical settings.

5. Conclusion

In this work, we have examined the nonlinear evolution of low-frequency modulated ion-acoustic waves in an ultra-cold neutral electron-ion plasma, where the ions are treated as warm and viscous, and the electrons are modeled as adiabatically trapped and degenerate. Starting from the basic fluid equations and Poisson’s equation and adopting a reductive perturbation technique (the derivative expansion method), a dissipative nonlinear Schrödinger (DNLS) equation governing the slow modulation of the ion–acoustic envelope has been derived. The coefficients of this equation explicitly incorporate the effects of the degenerate electron temperature and the ion kinematic viscosity, thus providing a clear connection between the microscopic plasma parameters and the large-scale envelope dynamics.

Based on the DNLS equation, we have carried out a detailed modulational instability analysis. The signs and relative magnitudes of the dispersion and nonlinearity coefficients, encapsulated in the ratio P/Q, determine whether the carrier wave is modulationally stable or unstable: for P/Q > 0 the carrier is prone to modulational instability, whereas for P/Q < 0 the modulated wave remains stable. The results indicate that an increase in the degenerate electron temperature enhances the system’s effective nonlinearity, broadens the instability domain, and increases the modulation growth rate. On the other hand, the ion kinematic viscosity lowers the growth rate and narrows the unstable region so that viscous effects act as a damping agent and hinder the development of strongly localized envelope structures.

To complement the analytical treatment, we have integrated the DNLS equation numerically using a second-order split-step Fourier scheme. For appropriate choices of the plasma parameters and initial conditions, the simulations reproduce the emergence of Akhmediev breathers, Kuznetsov–Ma breathers, Peregrine solitons, and higher–order rogue–wave patterns on a finite background. In the absence of dissipation, these structures reach large peak amplitudes and display the characteristic focusing–defocusing cycles associated with modulational instability. Considering the ion kinematic viscosity in the plasma model results in clear attenuation of the background, a progressive reduction in the maximum amplitudes, and a delay in the onset of rogue waves and breathers. This confirms that viscous effects play a decisive role in limiting the intensity and lifetime of extreme envelope structures in such plasmas.

For quick reference, the key parametric thresholds governing modulation instability and dissipative breather/rogue wave formation are summarized below:

• MI condition: P/Q > 0 (dispersion and nonlinearity coefficients have the same sign).

• Critical wavenumber for MI: $K_c^2(\tau )=2Q{\left|{\varphi }_0\right|}^2{\mathrm{e}}^{-2\gamma \tau }/P$ . MI occurs when $K^2<K_c^2(\tau )$ ; stable propagation when $K^2>K_c^2(\tau )$ .

• Maximum instability growth rate: ${\mathrm{\Gamma }}_g(max)=\left(1-\pi /4\right)Q{\left|{\varphi }_0\right|}^2/\gamma$ achieved at modulation wavenumber $K=Q{\left|{\varphi }_0\right|}^2/P$ .

• Dissipation-induced stability time: ${\tau }_{\mathrm{m}\mathrm{a}\mathrm{x}}=1/2\gamma \mathrm{l}\mathrm{n}\left(2Q{\left|{\varphi }_0\right|}^2/{PK}^2\right)$ . After τ > τ_max, MI is suppressed and the wave returns to a stable state.

Taken as a whole, the analytical and numerical results clarify the respective roles of adiabatic trapping, electron degeneracy, and ion kinematic viscosity in controlling modulational stability and the subsequent evolution of ion-acoustic envelopes in ultra-cold neutral plasmas. The DNLS description obtained in this study provides a concrete framework for analyzing localized electrostatic structures in ultracold laboratory plasmas and in dense astrophysical settings, and can be naturally extended to include additional ingredients such as external magnetic fields, multicomponent plasma mixtures, or fractional transport models. Moreover, our analysis is restricted to one spatial dimension. While this captures the essential physics of modulational instability along the direction of wave propagation, it can be extended to transverse instabilities and two- or three-dimensional envelope dynamics. Finally, our derivation relies on the reductive perturbation method, which assumes weak nonlinearity; consequently, the resulting DNLS equation describes only small-amplitude waves and cannot account for strongly nonlinear structures.

The present study suggests several natural extensions that may be pursued in future work:

• It would be of interest to construct and analyze a fractional counterpart of the dissipative NLS equation derived in this study by introducing time- and/or space-fractional derivatives to model memory effects and nonlocal transport in ultra-cold plasmas.

• The present model could be generalized to include an external magnetic field and to allow for oblique propagation of low-frequency ion-acoustic wave packets. In that case, one may derive anisotropic (possibly vector) DNLS or fractional-DNLS equations and investigate how the combined action of magnetic-field–induced anisotropy, adiabatic trapping, and viscosity modifies the modulational instability thresholds and the morphology of breather and rogue-wave structures.

• Another natural extension is to consider multi-species configurations, for example, by adding positrons or heavy ion species while retaining the adiabatically trapped degenerate electron component and viscous effects. The resulting envelope equations may exhibit competing nonlinearities or multiple characteristic scales, leading to new classes of coupled breathers, mixed-polarity rogue waves, or multi-component localized structures.

• From a methodological perspective, it would be useful to develop and benchmark advanced numerical schemes specifically adapted to dissipative and fractional NLS-type models, extending the second-order split-step Fourier method used here.

• Finally, a meaningful step forward would be to confront the theoretical predictions with either laboratory experiments on ultra-cold neutral plasmas or large-scale kinetic simulations.