Investigating the Dynamics of a Unidirectional Wave Model: Soliton Solutions, Bifurcation, and Chaos Analysis

Abstract

The current work investigates a recently introduced unidirectional wave model, applicable in science and engineering to understand complex systems and phenomena. This investigation has two primary aims. First, it employs a novel modified Sardar sub-equation method, not yet explored in the literature, to derive new solutions for the governing model. Second, it analyzes the complex dynamical structure of the governing model using bifurcation, chaos, and sensitivity analyses. To provide a more accurate depiction of the underlying dynamics, they use quantum mechanics to explain the intricate behavior of the system. To illustrate the physical behavior of the obtained solutions, 2D and 3D plots, along with a phase plane analysis, are presented using appropriate parameter values. These results validate the effectiveness of the employed method, providing thorough and consistent solutions with significant computational efficiency. The investigated soliton solutions will be valuable in understanding complex physical structures in various scientific fields, including ferromagnetic dynamics, nonlinear optics, soliton wave theory, and fiber optics. This approach proves highly effective in handling the complexities inherent in engineering and mathematical problems, especially those involving fractional-order systems.

1. Introduction

Nonlinear partial differential equations (NLPDEs) are widely used to model various modern phenomena. Researchers in this field are interested in gaining deeper insights into these phenomena, driving the demand for accurate solutions to these models. The interactions between terms that result in behaviors like shocks, solitons, turbulence, and chaos make NLPDEs much more complex. In fields like physics, biology, finance, and engineering, where complex dynamics are frequently missed by linear assumptions, these equations are essential in modeling real-world phenomena [1,2,3,4]. Consequently, a variety of analytical and numerical techniques have been developed, including the inverse scattering transform and others [5]. However, no single method universally provides exact solutions for all NLPDEs; a technique that is effective for one model may be unsuitable for another [6]. These approaches, introduced in the literature by various researchers, include the extended modified tanh method [7], the Sardar sub-equation method [8], the Adomian decomposition method [9], the variational iteration method [10], the unified method [11], the Darboux transform method [12], the Hirota bilinear method [13,14], the ${\displaystyle \frac{G^{\prime }}{G^2}}$-expansion method [15,16], the F-expansion method [17], phase portrait analysis [18], the $\Psi \left(\xi \right)$ method [19], the bilinear network technique [20], the biliear residual network technique [21], and many more [22,23,24]. In addition, the growing application of fractional calculus has opened up new pathways for the investigation of the dynamics of nonlinear systems, allowing researchers to model more complex physical phenomena with greater accuracy.

Recent research on NLPDEs has led to the discovery of several novel types of solitary wave solutions. For instance, a study of soliton molecules, Y-type waves, and complex multiple soliton solutions for the extended (3+1)-D model is presented in [25]. Perturbations of the Gerdjikov–Ivanov model and the derivation of soliton solutions through the Backlund transformation are detailed in [26]. Research on Painlevé integrability and the nonlinear characteristics of a (3+1)-dimensional Boussinesq-type equation in fluid dynamics, including lump and multiple soliton/shock solutions, can be found in [27]. The optimal system, invariant solutions, and soliton dynamics for the Wazwaz–Benjamin–Bona–Mahony equation are explored in [28]. Furthermore, [29] investigates exact soliton solutions and the role of time-dependent coefficients in the Boussinesq equation, highlighting its theoretical implications and applications in mathematical physics. In [30], the authors use wave models to study different behaviors of solutions. In [31], the authors extract different dynamical structures of wave models in quantum mechanics using analytical techniques. Some works on dynamical features and soliton solutions using different approaches include [32,33,34,35]. Motivated by these valuable contributions, this paper analyzes the complex dynamics of the following recently introduced model [36], whose dynamics have not yet been explored in the literature:

$$\mathcal{U}t+{\displaystyle \frac{3}{2}}\alpha \mathcal{U}\mathcal{U}x+{\displaystyle \frac{1}{6}}\beta \mathcal{U}xxx+{\displaystyle \frac{15}{32}}{\alpha }^2{\mathcal{U}}^2\mathcal{U}x+\mathrm{a}\mathcal{U}x+\mathrm{b}\mathcal{U}y+\mathrm{c}{\mathcal{U}}_z=0,$$

(1)

where x, y, and z represent spatial variables; t is the temporal variable; and $\alpha$, $\beta$, $\mathrm{a}$, $\mathrm{b}$, and $\mathrm{c}$ are constants. Because it can represent intricate wave interactions that are crucial in applications such as plasma physics and soliton theory, the new unidirectional wave model introduced in this research has great significance in nonlinear dynamics. In contrast to the current models, this equation places a strong emphasis on unidirectional propagation, which makes the analysis easier while preserving crucial nonlinear features. Since it captures realistic circumstances in media that enable solitary wave structures, this approach is especially useful in systems where wave dispersion or shock production is common.

By concentrating on wave behavior in a single direction, unidirectional wave models simplify complex wave dynamics. Wave breaking is studied using a variety of models, such as the Korteweg–de Vries equation and non-local wave models. These models help to forecast the evolution and kinematics of waves in irregular wave trains. To comprehend wave behavior and its impacts in coastal regions, the characteristics of unidirectional wave spectra are also investigated. To learn more about wave propagation, computational simulations are used to study solitary wave solutions in unidirectional wave models.

This article is organized in the following manner. Section 2 covers a brief overview of the MSSE method. Section 3 examines the governing model by converting it into a dynamical system and studies the phase portrait behaviors by using bifurcation, chaos, and sensitivity analysis. Section 4 illustrates how to apply the MSSE method to Equation (1) and obtain soliton solutions. Section 5 explores the physical explanation of the developed solutions. Lastly, Section 6 synthesizes the conclusions.

2. Description of Modified Sardar Sub-Equation Method

By including additional terms and clusters in the ansatz for the solution, the MSSE method broadens the scope of the original Sardar sub-equation method in the nonlinear phenomena. This method has been successfully used in physics and mathematics NLPDEs on multiple occasions. The generic form of NLPDEs is

$$\begin{array}{c}\hfill \mathcal{T}(\mathcal{U},{\mathcal{U}}_x,{\mathcal{U}}_t,{\mathcal{U}}_{xx},{\mathcal{U}}_{xz},{\mathcal{U}}_{xt},...)=0.\end{array}$$

(2)

Step 1. 

Utilize the transformation of the traveling wave solution

$$\begin{array}{c}\hfill \mathcal{U}=\mathcal{V}\left(\omega \right)\omega ={{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t.\end{array}$$

The NLODEs are obtained from Equation (2) as

$$\mathcal{R}(\mathcal{V},{\mathcal{V}}^{\prime },{\mathcal{V}}^{\prime \prime },...)=0.$$

(3)

Step 2. 

The given form describes the general solution of Equation (3), as per the method.

$$\mathcal{V}\left(\omega \right)={\mathcal{P}}_0+\sum _{N=1}^n{\mathcal{P}}_N{\mathcal{K}}^N\left(\omega \right),{\mathcal{P}}_N\ne 0,$$

(4)

where $\mathcal{V}=\mathcal{V}\left(\omega \right)$ ensures that

$${\mathcal{K}}^{\prime }\left(\omega \right)=\sqrt{q_2\mathcal{K}{\left(\omega \right)}^4+q_1\mathcal{K}{\left(\omega \right)}^2+q_0},$$

(5)

where the integers are $q_0$, $q_1$, and $q_2$, which are constants to be determined. ${\mathcal{K}}_0$ and ${\mathcal{K}}_1$ are calculated, and ${\mathcal{K}}_n$ is invertible when it is zero. Using the principle of balance, we determine the value of n. The clusters corresponding to Equation (5) are listed below.

Cluster 1. 

When $q_0=0,q_1>0\mathrm{and}{q}_2<0$, we acquire

$${\mathcal{K}}_1\left(\omega \right)=\sqrt{-{\displaystyle \frac{q_1}{q_2}}}\mathrm{sech}\left(\sqrt{q_1}(\omega +w)\right),$$

$${\mathcal{K}}_2\left(\omega \right)=\sqrt{-{\displaystyle \frac{q_1}{q_2}}}\mathrm{csch}\left(\sqrt{q_1}(\omega +w)\right).$$

Cluster 2. 

For constants $r_{1}and{r}_2$, when $q_0=0,q_1>0$ and $q_2=+4r_1{r}_2$, we acquire

$${\mathcal{K}}_3\left(\omega \right)={\displaystyle \frac{4r_1\sqrt{q_1}}{(4r_1^2-q_2)sinh\left(\sqrt{q_1}(\omega +w)\right)+(4r_1^2-q_2)cosh\left(\sqrt{q_1}(\omega +w)\right)}}.$$

Cluster 3. 

For constants $e_{1}and{e}_2$, when $q_0={\displaystyle \frac{q_1^2}{4q_2}},q_1<0and{q}_2>0$, we acquire

$${\mathcal{K}}_4\left(\omega \right)=\sqrt{-{\displaystyle \frac{q_1}{2q_2}}}tanh\left(\sqrt{-{\displaystyle \frac{q_1}{2}}}(\omega +w)\right),$$

$${\mathcal{K}}_5\left(\omega \right)=\sqrt{-{\displaystyle \frac{q_1}{2q_2}}}coth\left(\sqrt{-{\displaystyle \frac{q_1}{2}}}(\omega +w)\right),$$

$${\mathcal{K}}_6\left(\omega \right)=\sqrt{-{\displaystyle \frac{q_1}{2q_2}}}\left(tanh\left(\sqrt{-{\displaystyle \frac{q_1}{2}}}(\omega +w)\right)+i\mathrm{sech}\left(\sqrt{-2q_1}(\omega +w)\right)\right),$$

$${\mathcal{K}}_7\left(\omega \right)=\sqrt{-{\displaystyle \frac{q_1}{8q_2}}}\left(tanh\left(\sqrt{-{\displaystyle \frac{q_1}{8}}}(\omega +w)\right)+coth\left(\sqrt{-{\displaystyle \frac{q_1}{8}}}(\omega +w)\right)\right),$$

$${\mathcal{K}}_8\left(\omega \right)={\displaystyle \frac{\sqrt{-\frac{q_1}{2q_2}}(\sqrt{e_1^2+e_2^2}-e_{1}cosh(\sqrt{-2q_1}(\omega +w)\left)\right)}{e_{1}sinh\left(\sqrt{-2q_1}(\omega +w)\right)+e_2}},$$

$${\mathcal{K}}_9\left(\omega \right)={\displaystyle \frac{\sqrt{-\frac{q_1}{2q_2}}cosh\left(\sqrt{-2q_1}(\omega +w)\right)}{sinh\left(\sqrt{-2q_1}(\omega +w)\right)+i}}.$$

Cluster 4. 

When $q_0=0,q_1<0and{q}_2\ne 0$ are constants, we acquire

$${\mathcal{K}}_{10}\left(\omega \right)=\sqrt{-{\displaystyle \frac{q_1}{q_2}}}sec\left(\sqrt{-q_1}(\omega +w)\right),q_1<0,$$

$${\mathcal{K}}_{11}\left(\omega \right)=\sqrt{-{\displaystyle \frac{q_1}{q_2}}}csc\left(\sqrt{-q_1}(\omega +w)\right),q_1<0.$$

Cluster 5. 

When $q_0={\displaystyle \frac{q_1^2}{4q_2}},q_1>0$ and $q_2>0$ and $e_1^2-e_2^2>0$ are constants, we acquire

$${\mathcal{K}}_{12}\left(\omega \right)=\sqrt{{\displaystyle \frac{q_1}{2q_2}}}tan\left(\sqrt{{\displaystyle \frac{q_1}{2}}}(\omega +w)\right),$$

$${\mathcal{K}}_{13}\left(\omega \right)=-\sqrt{{\displaystyle \frac{q_1}{2q_2}}}cot\left(\sqrt{{\displaystyle \frac{q_1}{2}}}(\omega +w)\right),$$

$${\mathcal{K}}_{14}\left(\omega \right)=-\sqrt{{\displaystyle \frac{q_1}{2q_2}}}(tan(\sqrt{2q_1}(\omega +w))-sec(\sqrt{2q_1}(\omega +w)\left)\right),$$

$${\mathcal{K}}_{15}\left(\omega \right)=\sqrt{{\displaystyle \frac{q_1}{8q_2}}}(tan(\sqrt{{\displaystyle \frac{q_1}{8}}}(\omega +w))-cot(\sqrt{{\displaystyle \frac{q_1}{8}}}(\omega +w)\left)\right),$$

$${\mathcal{K}}_{16}\left(\omega \right)={\displaystyle \frac{\sqrt{\frac{q_1}{2q_2}}\left(\sqrt{e_1^2-e_2^2}-{\mathfrak{A}}_{1}cos\left(\sqrt{2q_1}(\omega +w)\right)\right)}{e_2+{\mathfrak{A}}_{1}sin\left(\sqrt{2q_1}(\omega +w)\right)}},$$

$${\mathcal{K}}_{17}\left(\omega \right)={\displaystyle \frac{\sqrt{\frac{q_1}{2q_2}}cos\left(\sqrt{2q_1}(\omega +w)\right)}{sin\left(\sqrt{2q_1}(\omega +w)\right)-1}}.$$

Cluster 6. 

When $q_0=0,q_1>0,$ we acquire

$${\mathcal{K}}_{18}\left(\omega \right)={\displaystyle \frac{4q_1{\mathrm{e}}^{\sqrt{q_1}(\omega +w)}}{{\mathrm{e}}^{2\sqrt{q_1}(\omega +w)}-4q_1{q}_2}},$$

$${\mathcal{K}}_{19}\left(\omega \right)={\displaystyle \frac{4q_1{\mathrm{e}}^{\sqrt{q_1}(\omega +w)}}{1-4q_1{q}_2{\mathrm{e}}^{2\sqrt{q_1}(\omega +w)}}}.$$

Cluster 7. 

When $q_0=0,q_1=0and{q}_2>0$, we acquire

$${\mathcal{K}}_{20}\left(\omega \right)={\displaystyle \frac{1}{\sqrt{q_2}(\omega +w)}},$$

$${\mathcal{K}}_{21}\left(\omega \right)={\displaystyle \frac{i}{\sqrt{q_2}(\omega +w)}}.$$

Step 3. 

After applying Equation (4) to Equation (3) and the second-order derivatives of Equation (4) using Equation (5), a polynomial with a power of $\mathcal{K}\left(\omega \right)$ is acquired.

Step 4. 

Collect all the $\mathcal{K}\left(\omega \right)$ coefficients with the same power and set them to zero; then, we acquire the algebraic system of the equation for ${\mathcal{P}}_0,{\mathcal{P}}_n$, where $n=1,2,3,...$.

Step 5. 

At the final step, utilize Wolfram Mathematica or any other computational software to solve the algebraic systems of the equations and determine the values of the parameters. The solutions to Equation (1) are acquired by plugging these parameter values into Equation (3).

3. Dynamical System Governed by Proposed Equation

This section includes the conversion of the proposed Equation (1) into a system of ODEs. First, we use traveling wave transformation to convert the considered PDE into the ODE. Then, we convert the proposed ODE into a system of ODEs using the Galilean transformation. Therefore, consider the following:

$$\mathcal{U}\left(x,y,z,t\right)=\mathcal{V}\left(\omega \right),$$

(6)

where $\omega ={{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t$. Inserting Equation (6) into Equation (1), we derive the following nonlinear ODE.

$$\left(\mathrm{a}{{\rm Y}}_1+\mathrm{b}{{\rm Y}}_2+\mathrm{c}{{\rm Y}}_3-{{\rm Y}}_4+{\displaystyle \frac{3}{2}}{{\rm Y}}_1\alpha \mathcal{V}\left(\omega \right)+{\displaystyle \frac{15}{32}}{{\rm Y}}_1{\alpha }^2\mathcal{V}{\left(\omega \right)}^2\right){\mathcal{V}}^{\prime }\left(\omega \right)+{\displaystyle \frac{1}{6}}{{\rm Y}}_1^3\beta {\mathcal{V}}^{\left(3\right)}\left(\omega \right)=0,$$

(7)

On performing the integration of Equation (7) with respect to the $\omega$ only once and supposing the constant of integration to be zero, we obtain the result presented below:

$$\left(\mathrm{a}{{\rm Y}}_1+\mathrm{b}{{\rm Y}}_2+\mathrm{c}{{\rm Y}}_3-{{\rm Y}}_4\right)\mathcal{V}\left(\omega \right)+{\displaystyle \frac{3}{2}}{{\rm Y}}_1\alpha \mathcal{V}{\left(\omega \right)}^2+{\displaystyle \frac{5}{32}}{{\rm Y}}_1{\alpha }^2\mathcal{V}{\left(\omega \right)}^3+{\displaystyle \frac{1}{6}}{{\rm Y}}_1^3\beta {\mathcal{V}}^{\prime \prime }\left(\omega \right)=0,$$

(8)

Here, making use of the Galilean transformation, Equation (8) gives rise to the following system of ODEs:

$$\begin{array}{c}\hfill \{\begin{array}{c}{\displaystyle \frac{d{\mathcal{S}}_1\left(\omega \right)}{d\omega }}=S_2\left(\omega \right),\hfill \\ {\displaystyle \frac{d{\mathcal{S}}_2\left(\omega \right)}{d\omega }}={\mathcal{G}}_1{\mathcal{S}}_1\left(\omega \right)+{\mathcal{G}}_2{\mathcal{S}}_1{\left(\omega \right)}^2+{\mathcal{G}}_3{\mathcal{S}}_1{\left(\omega \right)}^3,\hfill \end{array}\end{array}$$

(9)

where

$$\begin{array}{ccc}\hfill {\mathcal{G}}_1& =& -{\displaystyle \frac{6\left(\mathrm{a}{{\rm Y}}_1+\mathrm{b}{{\rm Y}}_2+\mathrm{c}{{\rm Y}}_3-{{\rm Y}}_4\right)}{{{\rm Y}}_1^3\beta }},\hfill \\ \hfill {\mathcal{G}}_2& =& -{\displaystyle \frac{9\alpha }{{{\rm Y}}_1^2\beta }},\hfill \\ \hfill {\mathcal{G}}_3& =& -{\displaystyle \frac{15{\alpha }^2}{16{{\rm Y}}_1^2\beta }}.\hfill \end{array}$$

3.1. Analysis and Graphical Visualization of Bifurcation, Chaos, and Other Behaviors of Equation (9)

This portion of the present work is devoted to the comprehensive study and analysis of the bifurcations, chaotic dynamics, sensitivity, and other analyses of the proposed system of ODEs.

3.1.1. Analysis of Bifurcations

Here, we aim to present the analysis of bifurcations. Further, we graphically demonstrate the types of fixed points that the system has. First, let us consider the Hamiltonian function of Equation (9) as presented below:

$$\begin{array}{c}\hfill H({\mathcal{S}}_1,{\mathcal{S}}_2)={\displaystyle \frac{1}{2}}{\mathcal{S}}_2^2+{\displaystyle \frac{{\mathcal{G}}_1}{2}}{\mathcal{S}}_1^2+{\displaystyle \frac{{\mathcal{G}}_2}{3}}{\mathcal{S}}_1^3+{\displaystyle \frac{{\mathcal{G}}_3}{4}}{\mathcal{S}}_1^4=\mathsf{h},\end{array}$$

where $\mathsf{h}$ denotes the Hamiltonian constant. In order to obtain the equilibrium points, let us consider Equation (9), as follows:

$$\begin{array}{c}\hfill \left\{\begin{array}{c}{\mathcal{S}}_2=0\hfill \\ {\mathcal{G}}_1{\mathcal{S}}_1+{\mathcal{G}}_2{\mathcal{S}}_1^2+{\mathcal{G}}_3{\mathcal{S}}_1^3=0.\hfill \end{array}\right.\end{array}$$

On solving the above system for ${\mathcal{S}}_1$ and ${\mathcal{S}}_2$, we obtain the following stationary points:

$${\mathcal{P}}_1=\left(0,0\right),{\mathcal{P}}_2=\left({\displaystyle \frac{-{\mathcal{G}}_2-\sqrt{{\mathcal{G}}_2^2-4{\mathcal{G}}_1{\mathcal{G}}_3}}{2{\mathcal{G}}_3}},0\right),{\mathcal{P}}_3=\left({\displaystyle \frac{-{\mathcal{G}}_2+\sqrt{{\mathcal{G}}_2^2-4{\mathcal{G}}_1{\mathcal{G}}_3}}{2{\mathcal{G}}_3}},0\right).$$

The determinant of the Jacobian of (9) is

$$\begin{array}{ccc}\hfill \mathcal{D}\left({\mathcal{S}}_1,{\mathcal{S}}_2\right)& =& \left|\begin{array}{cc}0& 1\\ {\mathcal{G}}_1+2{\mathcal{G}}_2{\mathcal{S}}_1+3{\mathcal{G}}_3{\mathcal{S}}_1^2& 0\end{array}\right|=-({\mathcal{G}}_1+2{\mathcal{G}}_2{\mathcal{S}}_1+3{\mathcal{G}}_3{\mathcal{S}}_1^2).\hfill \end{array}$$

We know that

• $({\mathcal{S}}_1,{\mathcal{S}}_2)$ will be a saddle point, when $\mathcal{D}({\mathcal{S}}_1,{\mathcal{S}}_2)<0$;
• $({\mathcal{S}}_1,{\mathcal{S}}_2)$ will be a center point, when $\mathcal{D}({\mathcal{S}}_1,{\mathcal{S}}_2)>0$;
• $({\mathcal{S}}_1,{\mathcal{S}}_2)$ will be a cuspid point, when $\mathcal{D}({\mathcal{S}}_1,{\mathcal{S}}_2)=0$.

The possible cases following the different conditions are presented below.

Case I: 

${\mathcal{G}}_1>0,{\mathcal{G}}_2>0$ and ${\mathcal{G}}_3>0$

By choosing particular values for the parameters in the form $\mathrm{a}=1,{{\rm Y}}_1=1,\mathrm{b}=1,{{\rm Y}}_2=1,c=1,{{\rm Y}}_3=1,{{\rm Y}}_4=1,\beta =-2,\alpha =1$, we observe that there is a single real stationary point, which is $(0,0)$, as visualized in Figure 1a. Clearly, it can be seen that $\left(0,0\right)$ is a saddle point.

Case II: 

${\mathcal{G}}_1<0,{\mathcal{G}}_2>0$ and ${\mathcal{G}}_3>0$

Using the values of the parameters as $\mathrm{a}=1,{{\rm Y}}_1=-2,\mathrm{b}=1,{{\rm Y}}_2=1,\mathrm{c}=1,{{\rm Y}}_3=1,{{\rm Y}}_4=1,\beta =-1,\alpha =1$, we find that there exist three stationary points $\left(0,0\right)$, $\left(-5.39333,0\right)$ and $\left(0.593326,0\right)$, out of which $\left(0,0\right)$ behaves as a cuspid point, as demonstrated in Figure 1b. Moreover, $\left(-5.39333,0\right)$ and $\left(0.593326,0\right)$ are saddle points.

Case III: 

${\mathcal{G}}_1<0,{\mathcal{G}}_2<0$ and ${\mathcal{G}}_3>0$

By selecting the parameters $\mathrm{a}=1,{{\rm Y}}_1=-2,\mathrm{b}=1,{{\rm Y}}_2=1,\mathrm{c}=1,{{\rm Y}}_3=1,{{\rm Y}}_4=1,\beta =-1,\alpha =-1$, we find that there exist three non-complex stationary points, $\left(0,0\right)$, $\left(-0.593326,0\right)$ and $(0,5.39333)$, as shown in Figure 1c. Obviously, one can see that the stationary point $(0,0)$ acts as a cuspid point, with $\left(-0.593326,0\right)$ and $\left(0,5.39333\right)$ corresponding to saddle points.

Case IV: 

${\mathcal{G}}_1<0,{\mathcal{G}}_2<0$ and ${\mathcal{G}}_3<0$

By selecting the parameters $\mathrm{a}=1,{{\rm Y}}_1=2,\mathrm{b}=1,{{\rm Y}}_2=1,\mathrm{c}=1,{{\rm Y}}_3=1,{{\rm Y}}_4=1,\beta =1,\alpha =1$, we recognize that the non-complex real stationary point is $\left(0,0\right)$, as depicted in Figure 1d. Evidently, $\left(0,0\right)$ is a center point.

Case IV: 

${\mathcal{G}}_1>0,{\mathcal{G}}_2<0$ and ${\mathcal{G}}_3<0$

By selecting the parameters $\mathrm{a}=1,{{\rm Y}}_1=2,\mathrm{b}=1,{{\rm Y}}_2=1,\mathrm{c}=1,{{\rm Y}}_3=1,{{\rm Y}}_4=1,\beta =1,\alpha =1$, we obtain three real stationary points, which are $(0,0)$, $(0.593326,0)$ and $(-5.39333,0)$, as presented in Figure 1e. Obviously, $(0,0)$ is a saddle, $(0.593326,0)$ is a cuspid and $(-5.39333,0)$ is a center point.

3.1.2. Chaos in the Proposed System

This sub-portion explores the possible existence of some chaos in system (9) via the introduction of a perturbation term. This section analyzes 2D and 2D vs. time phase diagrams for the proposed system as well. After adding the perturbation, the following is obtained:

$$\begin{array}{c}\hfill \{\begin{array}{c}{\displaystyle \frac{d{\mathcal{S}}_1\left(t\right)}{d\delta }}=S_2\left(t\right),\hfill \\ {\displaystyle \frac{d{\mathcal{S}}_2\left(t\right)}{d\delta }}={\mathcal{G}}_1{\mathcal{S}}_1\left(t\right)+{\mathcal{G}}_2{\mathcal{S}}_1{\left(t\right)}^2+{\mathcal{G}}_3{\mathcal{S}}_1{\left(t\right)}^3+\Omega cos\left(\phi t\right),\hfill \end{array}\end{array}$$

(10)

In Figure 2 and Figure 3, we analyze the impact of the term $\Omega cos\left(\phi t\right)$ on the dynamical evolution of system (10). Here, $\Omega$ stands for the amplitude, while $\phi$ stands for the frequency of the proposed system.

The 2D and 2D vs. time phase diagrams for the analysis of the system are presented with the utilization of the values as $\mathrm{a}=1,{{\rm Y}}_1=2,\mathrm{b}=1,{{\rm Y}}_2=1,\mathrm{c}=1,{{\rm Y}}_3=1,{{\rm Y}}_4=1,\beta =1,\alpha =1$, while the frequency and amplitudes are varied as $\Omega$ and $\phi$ are varied as in Figure 2 [(a),(b)] $\Omega =1,\phi =1$ and in [(c),(d)] $\Omega =0.1,\phi =3.2$. Further, in Figure 3, we consider for [a),(b)] $\Omega =1,\phi =0.9$ and in [(c),(d)] $\Omega =0.5,\phi =5.1$.

After the analysis of the phase portraits, one can see very complex and interesting dynamics. In Figure 2a, we see multi-scroll torus-like dynamics, while, in Figure 2c, complex and unusual dynamics are obtained. Moreover, in Figure 3a, complex multi-scroll behavior can be seen, while Figure 3c shows heart-shaped oscillations. This shows the system’s behaviors’ susceptibility to the perturbations arising in $\phi$, offering significant insights into the manner in which the perturbation term $\Omega cos\left(\phi t\right)$ influences the behavior of the proposed system. These new understandings of the system’s vulnerability to parameter variations enhance our comprehension of the complex relationships among $\phi$ and the overall dynamics of the system. These findings facilitate a broader understanding of the ways in which small changes can affect the trajectories of the proposed dynamical system, eventually paving the way for more accurate and informed forecasts of its behavior in different conditions.

3.1.3. Effects of Parameters on the Chaotic Flow of the System

This part of the manuscript considers different parameter values for the parameters $\beta$ and ${{\rm Y}}_4$, to observe their effects on the dynamics of the proposed system’s evolution. Therefore, Figure 4a–c show the effects of varying $\beta$. The values of $\beta$ considered for Figure 4a–c are $0.71$, $1.99$ and $6.5$, respectively. We observe multi-scroll dynamics in Figure 4a, four-scroll dynamics in Figure 4b and finally oval-shaped multi-scroll attraction in Figure 4c. Furthermore, the values of ${{\rm Y}}_4$ considered in Figure 4d–f are $0.1$, $0.5$ and $6.5$, respectively. Here, we see that there are multi-scroll dynamics of a torus-like nature in Figure 4d, multi-scroll dynamics in Figure 4e and finally O-shaped attraction in Figure 4f.

3.1.4. Sensitivity Analysis

Here, we study the sensitivity of the proposed dynamical model presented in Equation (9). For this, let us consider the following dynamical model:

$$\begin{array}{c}\hfill \left\{\begin{array}{c}{\displaystyle \frac{d{\mathcal{S}}_1\left(\mathsf{t}\right)}{d\mathsf{t}}}={\mathcal{S}}_2,\hfill \\ {\displaystyle \frac{d{\mathcal{S}}_2\left(\mathsf{t}\right)}{d\mathsf{t}}}={\mathcal{S}}_1{\mathcal{S}}^3\left(\mathsf{t}\right)+{\mathcal{S}}_1{\mathcal{S}}_1\left(\mathsf{t}\right).\hfill \end{array}\right.\end{array}$$

(11)

The parameters are considered as follows: $\mathrm{a}=1,{{\rm Y}}_1=2,\mathrm{b}=1,{{\rm Y}}_2=1,\mathrm{c}=1,{{\rm Y}}_3=1,{{\rm Y}}_4=1,\beta =1,\alpha =1$. Meanwhile, the initial values are used as the blue waves to represent the behavior of the system with $({\mathcal{S}}_1\left(0\right),{\mathcal{S}}_2\left(0\right))=(0.1,0)$. Similarly, for the green oscillatory dynamics of the proposed system, the initial values are supposed to be $({\mathcal{S}}_1\left(0\right),{\mathcal{S}}_2\left(0\right))=(0.3,0)$, while, for the red-colored oscillations, the initial values are $({\mathcal{S}}_1\left(0\right),{\mathcal{S}}_2\left(0\right))=(0.6,0)$.

The obtained results are demonstrated in Figure 5. Observations of the figures reveal that small perturbations in the initial values consequently provide huge changes in the dynamics of the proposed system.

4. Mathematical Analysis

This section primarily focuses on the use of our proposed technique to validate its reliability, performance, and effectiveness. This will allow us to obtain the governing model’s soliton solution. We obtain $M=1$ by applying the balance principle found in Equation (8). Thus, using $M=1$, the generic solution provided in Equation (4) becomes

$$\mathcal{V}\left(\omega \right)={\mathcal{P}}_1\mathcal{K}\left(\omega \right)+{\mathcal{P}}_0.$$

(12)

A coefficient with the same power of ${\mathcal{K}\left(\omega \right))}^n$ is equated, where $n=0,1,2,3,...$. After inserting Equation (12) into Equation (8) and solving the system of equations through symbolic computation, we acquire the families and solutions listed below. This process yields a system of algebraic equations.

Family 1:

$$\begin{array}{c}\hfill \left\{{\mathcal{P}}_1\to -{\displaystyle \frac{8i\sqrt{2}\sqrt{q_2}}{5\alpha \sqrt{q_1}}},{\mathcal{P}}_0\to -{\displaystyle \frac{8}{5\alpha }},\beta \to {\displaystyle \frac{12}{5{{\rm Y}}_1^2{q}_1}},\mathrm{b}\to {\displaystyle \frac{(4-5\mathrm{a}){{\rm Y}}_1-5c{{\rm Y}}_3+5{{\rm Y}}_4}{5{{\rm Y}}_2}}\right\}.\end{array}$$

(13)

The following solutions satisfy Family 1 based on the analysis.

$$\begin{array}{c}\hfill {\mathcal{U}}_{1,1}=-{\displaystyle \frac{8}{5\alpha }}-{\displaystyle \frac{8i\sqrt{2}\sqrt{-\frac{q_1}{q_2}}\sqrt{q_2}\mathrm{sech}\left(\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)}{5\alpha \sqrt{q_1}}},q_2<0,\end{array}$$

(14)

$$\begin{array}{c}\hfill {\mathcal{U}}_{1,2}=-{\displaystyle \frac{8}{5\alpha }}-{\displaystyle \frac{8i\sqrt{2}\sqrt{-\frac{q_1}{q_2}}\sqrt{q_2}\mathrm{csch}\left(\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)}{5\alpha \sqrt{q_1}}},q_2<0,\end{array}$$

(15)

$$\begin{array}{c}\hfill {\mathcal{U}}_{1,3}=-{\displaystyle \frac{8}{5\alpha }}-{\displaystyle \frac{32i\sqrt{2}\sqrt{q_2}{r}_1}{5\alpha \left(\left(4r_1^2-q_2\right)sinh\left(\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)\right)}}\\ \hfill +{\displaystyle \frac{32i\sqrt{2}\sqrt{q_2}{r}_1}{+\left(4r_1^2-q_2\right)cosh\left(\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)}},q_2<0,\end{array}$$

(16)

$$\begin{array}{c}\hfill {\mathcal{U}}_{1,4}=-{\displaystyle \frac{8}{5\alpha }}-{\displaystyle \frac{8i\sqrt{-\frac{q_1}{q_2}}\sqrt{q_2}tanh\left(\frac{\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}{\sqrt{2}}\right)}{5\alpha \sqrt{q_1}}},q_1<0,\end{array}$$

(17)

$$\begin{array}{c}\hfill {\mathcal{U}}_{1,5}=-{\displaystyle \frac{8}{5\alpha }}-{\displaystyle \frac{8i\sqrt{-\frac{q_1}{q_2}}\sqrt{q_2}coth\left(\frac{\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}{\sqrt{2}}\right)}{5\alpha \sqrt{q_1}}},q_1<0,\end{array}$$

(18)

$$\begin{array}{c}\hfill {\mathcal{U}}_{1,6}=-{\displaystyle \frac{8}{5\alpha }}-{\displaystyle \frac{8i\sqrt{-\frac{q_1}{q_2}}\sqrt{q_2}\left(tanh\left(\sqrt{2}\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)\right)}{5\alpha \sqrt{q_1}}}\\ \hfill {\displaystyle \frac{+i\mathrm{sech}\left(\sqrt{2}\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)}{5\alpha \sqrt{q_1}}},q_1<0,\end{array}$$

(19)

$$\begin{array}{c}\hfill {\mathcal{U}}_{1,7}=-{\displaystyle \frac{8}{5\alpha }}-{\displaystyle \frac{4i\sqrt{-\frac{q_1}{q_2}}\sqrt{q_2}\left(tanh\left(\frac{\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}{2\sqrt{2}}\right)+icoth\left(\frac{\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}{2\sqrt{2}}\right)\right)}{5\alpha \sqrt{q_1}}},q_1<0,\end{array}$$

(20)

$$\begin{array}{c}\hfill {\mathcal{U}}_{1,8}=-{\displaystyle \frac{8}{5\alpha }}-{\displaystyle \frac{8i\sqrt{-\frac{q_1}{q_2}}\sqrt{q_2}\left(\sqrt{e_1^2+e_2^2}-e_{1}cosh\left(\sqrt{2}\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)\right)}{5\alpha \sqrt{q_1}\left(e_{1}sinh\left(\sqrt{2}\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)+e_2\right)}},q_1<0,\end{array}$$

(21)

$$\begin{array}{c}\hfill {\mathcal{U}}_{1,9}=-{\displaystyle \frac{8}{5\alpha }}-{\displaystyle \frac{8i\sqrt{-\frac{q_1}{q_2}}\sqrt{q_2}cosh\left(\sqrt{2}\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)}{5\alpha \sqrt{q_1}\left(sinh\left(\sqrt{2}\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)+i\right)}},q_1<0,\end{array}$$

(22)

$$\begin{array}{c}\hfill {\mathcal{U}}_{1,10}=-{\displaystyle \frac{8}{5\alpha }}-{\displaystyle \frac{8i\sqrt{2}\sqrt{\frac{-q_1}{q_2}}\sqrt{q_2}sec\left(\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)}{5\alpha \sqrt{q_1}}},q_1<0,\end{array}$$

(23)

$$\begin{array}{c}\hfill {\mathcal{U}}_{1,11}=-{\displaystyle \frac{8}{5\alpha }}-{\displaystyle \frac{8i\sqrt{2}\sqrt{\frac{-q_1}{q_2}}\sqrt{q_2}csc\left(\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)}{5\alpha \sqrt{q_1}}},q_1<0,\end{array}$$

(24)

$$\begin{array}{c}\hfill {\mathcal{U}}_{1,12}=-{\displaystyle \frac{8}{5\alpha }}-{\displaystyle \frac{8i\sqrt{\frac{q_1}{q_2}}\sqrt{q_2}tan\left(\frac{\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}{\sqrt{2}}\right)}{5\alpha \sqrt{q_1}}},\end{array}$$

(25)

$$\begin{array}{c}\hfill {\mathcal{U}}_{1,13}=-{\displaystyle \frac{8}{5\alpha }}-{\displaystyle \frac{8i\sqrt{\frac{q_1}{q_2}}\sqrt{q_2}cot\left(\frac{\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}{\sqrt{2}}\right)}{5\alpha \sqrt{q_1}}},\end{array}$$

(26)

$$\begin{array}{c}\hfill {\mathcal{U}}_{1,14}=-{\displaystyle \frac{8}{5\alpha }}+{\displaystyle \frac{8i\sqrt{\frac{q_1}{q_2}}\sqrt{q_2}\left(tan\left(\sqrt{2}\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)\right)}{5\alpha \sqrt{q_1}}}\\ \hfill \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}-{\displaystyle \frac{sec\left(\sqrt{2}\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)}{5\alpha \sqrt{q_1}}},\end{array}$$

(27)

$$\begin{array}{c}\hfill {\mathcal{U}}_{1,15}=-{\displaystyle \frac{8}{5\alpha }}-{\displaystyle \frac{4i\sqrt{\frac{q_1}{q_2}}\sqrt{q_2}\left(tan\left(\frac{\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}{2\sqrt{2}}\right)-cot\left(\frac{\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}{2\sqrt{2}}\right)\right)}{5\alpha \sqrt{q_1}}},\end{array}$$

(28)

$$\begin{array}{c}\hfill {\mathcal{U}}_{1,16}=-{\displaystyle \frac{8}{5\alpha }}-{\displaystyle \frac{8i\sqrt{\frac{q_1}{q_2}}\sqrt{q_2}\left(\sqrt{e_1^2-e_2^2}-e_{1}cos\left(\sqrt{2}\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)\right)}{5\alpha \sqrt{q_1}\left(e_{1}sin\left(\sqrt{2}\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)+e_2\right)}},\end{array}$$

(29)

$$\begin{array}{c}\hfill {\mathcal{U}}_{1,17}=-{\displaystyle \frac{8}{5\alpha }}-{\displaystyle \frac{8i\sqrt{\frac{q_1}{q_2}}\sqrt{q_2}cot\left(\sqrt{2}\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)}{5\alpha \sqrt{q_1}}},\end{array}$$

(30)

$$\begin{array}{c}\hfill {\mathcal{U}}_{1,18}=-{\displaystyle \frac{8}{5\alpha }}-{\displaystyle \frac{32i\sqrt{2}\sqrt{q_1}\sqrt{q_2}{e}^{\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}}{5\alpha \left(e^{2\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}-4q_1{q}_2\right)}},q_2<0,\end{array}$$

(31)

$$\begin{array}{c}\hfill {\mathcal{U}}_{1,19}=-{\displaystyle \frac{8}{5\alpha }}-{\displaystyle \frac{32i\sqrt{2}\sqrt{q_1}\sqrt{q_2}{e}^{\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}}{5\alpha \left(1-4q_1{q}_2{e}^{2\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}\right)}},q_2<0,\end{array}$$

(32)

$$\begin{array}{c}\hfill {\mathcal{U}}_{1,20}=-{\displaystyle \frac{8}{5\alpha }}-{\displaystyle \frac{8i\sqrt{2}}{5\alpha \sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}},q_1<0,\end{array}$$

(33)

$$\begin{array}{c}\hfill {\mathcal{U}}_{1,21}={\displaystyle \frac{8\sqrt{2}\sqrt{q_2}}{5\alpha \sqrt{q_1}\sqrt{-q_2}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}}-{\displaystyle \frac{8}{5\alpha }},q_2<0.\end{array}$$

(34)

Family 2:

$$\begin{array}{c}\hfill \left\{{\mathcal{P}}_1\to {\displaystyle \frac{8i\sqrt{2}\sqrt{q_2}}{5\alpha \sqrt{q_1}}},{\mathcal{P}}_0\to -{\displaystyle \frac{8}{5\alpha }},\beta \to {\displaystyle \frac{12}{5{{\rm Y}}_1^2{q}_1}},a\to {\displaystyle \frac{4{{\rm Y}}_1-5b{{\rm Y}}_2-5c{{\rm Y}}_3+5{{\rm Y}}_4}{5{{\rm Y}}_1}}\right\}.\end{array}$$

(35)

The following solutions satisfy Family 2 based on the analysis.

$$\begin{array}{c}\hfill {\mathcal{U}}_{2,1}=-{\displaystyle \frac{8}{5\alpha }}+{\displaystyle \frac{8i\sqrt{2}\sqrt{-\frac{q_1}{q_2}}\sqrt{q_2}\mathrm{sech}\left(\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)}{5\alpha \sqrt{q_1}}},\end{array}$$

(36)

$$\begin{array}{c}\hfill {\mathcal{U}}_{2,2}=-{\displaystyle \frac{8}{5\alpha }}+{\displaystyle \frac{8i\sqrt{2}\sqrt{-\frac{q_1}{q_2}}\sqrt{q_2}\mathrm{csch}\left(\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)}{5\alpha \sqrt{q_1}}},\end{array}$$

(37)

$$\begin{array}{cc}\hfill {\mathcal{U}}_{2,3}=& -{\displaystyle \frac{8}{5\alpha }}+{\displaystyle \frac{32i\sqrt{2}\sqrt{q_2}{r}_1}{5\alpha \left(\left(4r_1^2-q_2\right)sinh\left(\sqrt{q_1}\right({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\left)\right)\right)}}\end{array}$$

(38)

$$\begin{array}{cc}& \hspace{1em}\hspace{1em}\hspace{1em}+{\displaystyle \frac{32i\sqrt{2}\sqrt{q_2}{r}_1}{+\left(4r_1^2-q_2\right)cosh\left(\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)}},\end{array}$$

(39)

$$\begin{array}{c}\hfill {\mathcal{U}}_{2,4}=-{\displaystyle \frac{8}{5\alpha }}+{\displaystyle \frac{8i\sqrt{-\frac{q_1}{q_2}}\sqrt{q_2}tanh\left(\frac{\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}{\sqrt{2}}\right)}{5\alpha \sqrt{q_1}}},\end{array}$$

(40)

$$\begin{array}{c}\hfill {\mathcal{U}}_{2,5}=-{\displaystyle \frac{8}{5\alpha }}+{\displaystyle \frac{8i\sqrt{-\frac{q_1}{q_2}}\sqrt{q_2}coth\left(\frac{\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}{\sqrt{2}}\right)}{5\alpha \sqrt{q_1}}},\end{array}$$

(41)

$$\begin{array}{cc}\hfill {\mathcal{U}}_{2,6}=& -{\displaystyle \frac{8}{5\alpha }}+{\displaystyle \frac{8i\sqrt{-\frac{q_1}{q_2}}\sqrt{q_2}\left(tanh\left(\sqrt{2}\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)\right)}{5\alpha \sqrt{q_1}}}\end{array}$$

(42)

$$\begin{array}{cc}& \hspace{1em}\hspace{1em}+{\displaystyle \frac{\left(i\mathrm{sech}\left(\sqrt{2}\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)\right)}{5\alpha \sqrt{q_1}}},\end{array}$$

(43)

$$\begin{array}{c}\hfill {\mathcal{U}}_{2,7}=-{\displaystyle \frac{8}{5\alpha }}+{\displaystyle \frac{4i\sqrt{-\frac{q_1}{q_2}}\sqrt{q_2}\left(tanh\left(\frac{\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}{2\sqrt{2}}\right)+icoth\left(\frac{\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}{2\sqrt{2}}\right)\right)}{5\alpha \sqrt{q_1}}},\end{array}$$

(44)

$$\begin{array}{c}\hfill {\mathcal{U}}_{2,8}=-{\displaystyle \frac{8}{5\alpha }}+{\displaystyle \frac{8i\sqrt{-\frac{q_1}{q_2}}\sqrt{q_2}\left(\sqrt{e_1^2+e_2^2}-e_{1}cosh\left(\sqrt{2}\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)\right)}{5\alpha \sqrt{q_1}\left(e_{1}sinh\left(\sqrt{2}\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)+e_2\right)}},\end{array}$$

(45)

$$\begin{array}{c}\hfill {\mathcal{U}}_{2,9}=-{\displaystyle \frac{8}{5\alpha }}+{\displaystyle \frac{8i\sqrt{-\frac{q_1}{q_2}}\sqrt{q_2}cosh\left(\sqrt{2}\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)}{5\alpha \sqrt{q_1}\left(sinh\left(\sqrt{2}\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)+i\right)}},\end{array}$$

(46)

$$\begin{array}{c}\hfill {\mathcal{U}}_{2,10}=-{\displaystyle \frac{8}{5\alpha }}+{\displaystyle \frac{8i\sqrt{2}\sqrt{-\frac{q_1}{q_2}}\sqrt{q_2}sec\left(\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)}{5\alpha \sqrt{q_1}}},\end{array}$$

(47)

$$\begin{array}{c}\hfill {\mathcal{U}}_{2,11}=-{\displaystyle \frac{8}{5\alpha }}+{\displaystyle \frac{8i\sqrt{2}\sqrt{-\frac{q_1}{q_2}}\sqrt{q_2}csc\left(\sqrt{-q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)}{5\alpha \sqrt{q_1}}},\end{array}$$

(48)

$$\begin{array}{c}\hfill {\mathcal{U}}_{2,12}=-{\displaystyle \frac{8}{5\alpha }}+{\displaystyle \frac{8i\sqrt{\frac{q_1}{q_2}}\sqrt{q_2}tan\left(\frac{\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}{\sqrt{2}}\right)}{5\alpha \sqrt{q_1}}},\end{array}$$

(49)

$$\begin{array}{c}\hfill {\mathcal{U}}_{2,13}=-{\displaystyle \frac{8}{5\alpha }}+{\displaystyle \frac{8i\sqrt{\frac{q_1}{q_2}}\sqrt{q_2}cot\left(\frac{\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}{\sqrt{2}}\right)}{5\alpha \sqrt{q_1}}},\end{array}$$

(50)

$$\begin{array}{cc}\hfill {\mathcal{U}}_{2,14}=& -{\displaystyle \frac{8}{5\alpha }}-{\displaystyle \frac{8i\sqrt{\frac{q_1}{q_2}}\sqrt{q_2}\left(tan\left(\sqrt{2}\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)\right)}{5\alpha \sqrt{q_1}}}\end{array}$$

(51)

$$\begin{array}{cc}& \hspace{1em}\hspace{1em}-{\displaystyle \frac{\left(sec\left(\sqrt{2}\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)\right)}{5\alpha \sqrt{q_1}}},\end{array}$$

(52)

$$\begin{array}{c}\hfill {\mathcal{U}}_{2,15}=-{\displaystyle \frac{8}{5\alpha }}+{\displaystyle \frac{4i\sqrt{\frac{q_1}{q_2}}\sqrt{q_2}\left(tan\left(\frac{\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}{2\sqrt{2}}\right)-cot\left(\frac{\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}{2\sqrt{2}}\right)\right)}{5\alpha \sqrt{q_1}}},\end{array}$$

(53)

$$\begin{array}{c}\hfill {\mathcal{U}}_{2,16}=-{\displaystyle \frac{8}{5\alpha }}+{\displaystyle \frac{8i\sqrt{\frac{q_1}{q_2}}\sqrt{q_2}\left(\sqrt{e_1^2-e_2^2}-e_{1}cos\left(\sqrt{2}\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)\right)}{5\alpha \sqrt{q_1}\left(e_{1}sin\left(\sqrt{2}\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)+e_2\right)}},\end{array}$$

(54)

$$\begin{array}{c}\hfill \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}{\mathcal{U}}_{2,17}=-{\displaystyle \frac{8}{5\alpha }}+{\displaystyle \frac{8i\sqrt{\frac{q_1}{q_2}}\sqrt{q_2}cot\left(\sqrt{2}\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)\right)}{5\alpha \sqrt{q_1}}},\end{array}$$

(55)

$$\begin{array}{c}\hfill {\mathcal{U}}_{2,18}=-{\displaystyle \frac{8}{5\alpha }}+{\displaystyle \frac{32i\sqrt{2}\sqrt{q_1}\sqrt{q_2}{e}^{\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}}{5\alpha \left(e^{2\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}-4q_1{q}_2\right)}},\end{array}$$

(56)

$$\begin{array}{c}\hfill {\mathcal{U}}_{2,19}=-{\displaystyle \frac{8}{5\alpha }}+{\displaystyle \frac{32i\sqrt{2}\sqrt{q_1}\sqrt{q_2}{e}^{\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}}{5\alpha \left(1-4q_1{q}_2{e}^{2\sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}\right)}},\end{array}$$

(57)

$$\begin{array}{c}\hfill {\mathcal{U}}_{2,20}=-{\displaystyle \frac{8}{5\alpha }}+{\displaystyle \frac{8i\sqrt{2}}{5\alpha \sqrt{q_1}\left({{\rm Y}}_{1}x+{{\rm Y}}_{2}y+{{\rm Y}}_{3}z-{{\rm Y}}_{4}t+w\right)}},\end{array}$$

(58)

5. Physical Interpretations of the Solutions

This section offers the physical representations and a comparison to earlier approaches, which frequently produce a limited number of solution types. The addition of rational exponential, trigonometric and hyperbolic solutions in this work offers a larger and more flexible solution set. These novel solutions provide several benefits. By contrasting these solutions with earlier findings, we show that the proposed technique allows for a greater variety of solution behaviors and more successfully handles complex initial values and a range of boundary conditions. Numerous real-world phenomena can be modeled by the various solution types. For example, hyperbolic solutions apply to shock waves and solitary wave patterns in fluid dynamics and optics, whereas trigonometric solutions are very helpful for periodic behavior in wave and oscillation investigations. The extracted solutions as shown in Figure 6, Figure 7, Figure 8, Figure 9 and Figure 10, which were obtained using the proposed analytical technique. The proposed technique extracts a wide range of soliton solutions, including dark, bright, kink and anti-kink types, each of which has useful applications in fluid mechanics and nonlinear optics. By increasing the analytical solution space and providing computational efficiency, this approach enables us to handle intricate, high-dimensional wave models without sacrificing the solution diversity or correctness. The demonstrated results are depicted in the form of 3D plots.

Figure 6 shows the dynamics of the analytical solution ${\mathcal{U}}_{1,6}$ with parameters $q_1=-1,q_2=1,w=1,{{\rm Y}}_1=1,{{\rm Y}}_2=0.1,{{\rm Y}}_3=0.9,{{\rm Y}}_4=0.2,\alpha =1,z=0,y=0$. From the simulation, we observe the kink soliton solution, which shows that the wave profile underwent a sharp transition between two distinct, stable states. Contrary to the regular waves, which disperse as time evolves, kink solitons preserve their shape and propagate at a uniform speed, suggesting the information or energy can travel in a medium without dissipation. Kink solitons are very important in real life, describing several scenarios, ranging from magnetic domain walls in materials science to cosmology and quantum field theory as well as optical fibers.

Figure 7 shows the dynamics of the analytical solution ${\mathcal{U}}_{1,7}$ with parameters $q_1=-1,q_2=1,w=1,{{\rm Y}}_1=1,{{\rm Y}}_2=0.1,{{\rm Y}}_3=0.9,{{\rm Y}}_4=0.2,\alpha =1,z=0,y=0$. From the simulation, a singular wave is observed, which is identified with sharp troughs or peaks that reach extreme values, even becoming infinite at some points, suggesting singularities in the wave profile. These dynamics often correspond to very intense physical phenomena, which the conventional waves fail to describe. In the field of nonlinear optics, singular waves describe pulse formations having extremely concentrated power, which is crucial in advancing high-intensity laser applications.

Figure 8 depicts the dynamics of the analytical solution ${\mathcal{U}}_{1,8}$ with the parameters $q_1=-1,q_2=1,w=1,e_1=1,e_2=1,{{\rm Y}}_1=1,{{\rm Y}}_2=0.1,{{\rm Y}}_3=0.9,{{\rm Y}}_4=0.2,\alpha =1,z=0,y=0$. From the simulation, an anti-kink soliton is observed, which is the reverse of a kink soliton, showing a rapid change between two different states in the opposite direction. These waves are important in those systems that support bidirectional wave propagation. For example, in condensed matter physics, anti-kinks are essential to understand how domain walls in ferroelectric or ferromagnetic materials reverse their direction, a phenomenon that underpins the functionality of devices used for memory, such as magnetic random-access memory. Similarly, in biology, the transmission of electric signals in neurons is modeled with kink and anti-kink solitons, which show how information moves from one place to another through the nervous system.

Figure 9 depicts the dynamics of the analytical solution ${\mathcal{U}}_{1,19}$ with the parameters $q_1=1,q_2=-1,w=1,{{\rm Y}}_1=1,{{\rm Y}}_2=0.1,{{\rm Y}}_3=0.9,{{\rm Y}}_4=0.2,\alpha =1,z=0,y=0$. From the simulation, a bright soliton is observed. The bright soliton is a localized wave with a pronounced peak, suggesting a region with high intensity. These waves result from a precise balance between dispersion and nonlinear effects, which force them to preserve their shape over long distances. Such types of waves find considerable applications in fiber optics, where they are utilized to transmit information with no distortion, even over very long distances. This characteristic gives them a great advantage and makes them integral to the advanced communication systems of the modern world, resulting in high-speed data transfer. Additionally, in hydrodynamics, these can be used to model shallow water waves, which preserve their integrity.

Figure 9 depicts the dynamics of the analytical solution ${\mathcal{U}}_{2,1}$ with the parameters $q_1=1,q_2=-1,w=1,{{\rm Y}}_1=1,{{\rm Y}}_2=0.1,{{\rm Y}}_3=0.9,{{\rm Y}}_4=0.2,\alpha =1,z=0,y=0$. From the simulation, a dark soliton is observed. These waves present depressions or dips in the wave profile, where the intensity decreases at the center as compared to the surroundings. Unlike bright solitons, they behave as “holes” in an otherwise uniform background. Expanding the proposed technique for larger classes of nonlinear equations, particularly those with complex or variable coefficients, is the goal of future research. This method’s usefulness could be increased by creating adaptive modifications that could handle a larger variety of differential equations. The process of solving complicated, high-dimensional equations could be streamlined by creating automated frameworks. This would increase the computing effectiveness and enable the solution of equations in high-dimensional areas without requiring a lot of manual labor. Through further modification, future research can investigate methods to obtain a wider range of solution forms, such as oscillatory or periodic solutions. These solitons are crucial in optical communications, where variations in the light intensity denote data bits. They facilitate the encoding and transmission of information. In fluids, they can model the dynamics of long internal waves, contributing to our understanding of complex oceanic and atmospheric interactions.

6. Conclusions

In this study, we used the modified Sardar sub-equation (MSSE) method to find multiple soliton solutions for the unidirectional wave model. Numerous significant solutions, such as dark, bright, singular, kink, and anti-kink solitons, were obtained using this technique. MATLAB (https://www.mathworks.com/products/matlab.html) was used to create 3D charts that showed the solutions. Additionally, we used sensitivity, bifurcation, and chaos studies to examine the model’s dynamic structure. These results provide valuable insights into the nonlinear characteristics of the model and establish a foundation for future studies in soliton dynamics and nonlinear phenomena in related systems. Potential avenues for future research include practical implementations of the derived soliton solutions and further analysis of their robustness to perturbations and parameter variations. A deeper understanding of bifurcations, chaos, and soliton dynamics contributes valuable insights for various scientific and mathematical fields. Future research should look into integrating the Hirota bilinear method and the MSSE method with other analytical approaches, as this could improve the accuracy of the solutions for intricate wave interactions in high-dimensional systems. This approach may also reveal new soliton forms and intricate behaviors that are crucial for domains such as fluid dynamics and nonlinear optics when applied to more generalized, multidimensional nonlinear wave equations.