#### Abstract

The Lorenz-Stenflo mathematical model describes a complex dynamical behavior related to atmospheric acoustic-gravity waves. In this study, qualitative analysis of the four-dimensional hyperchaotic Lorenz-Stenflo system via the Caputo fractional derivative is implemented. By using the Matignon stability criterion, the local stability analysis of the system showed that all the equilibrium points of the system are locally unstable. Calculation of the Lyapunov exponents along with the relevant bifurcation diagrams with respect to different fractional orders exposed the hyperchaotic dynamical behavior for the system. Bifurcation diagrams for all the four parameters in the system also showed the hyperchaotic nature of the Lorenz-Stenflo system. Different phase attractors of the system corresponding to different fractional derivatives and parameters are presented to specify the dynamical nature of the system. The Lorenz-Stenflo system showed sensitivity to initial conditions. The master and slave systems showed a strong correlation among themselves, as verified by graphs of time series solutions of the two systems.

#### 1. Introduction

In recent years, fractional operators have been applied to develop mathematical models for which we can investigate different dynamical systems in some areas such as mathematical biology, epidemiology, and engineering [1].

Different researchers introduced several concepts of fractional derivatives. Some of them are the Caputo fractional-order derivative [2–4], the Riemann Liouville fractional-order derivative [4], the Caputo-Fabrizio fractional-order derivative [5], and the Atangana-Baleanu fractional derivatives [6].

In the mathematical models with integer derivatives, the derivative orders give the instantaneous rate of change of the function. On the other hand, in the case of mathematical models with fractional derivatives, the parameters of dynamical systems represent the memory index of variation of the function [7].

Thus, the advantage of using the memory index of a dynamical system with fractional derivative made fractional derivatives advantageously applied in the fields of chaotic dynamics, epidemiological modeling, and several other fields. In particular, one can see these applications in the modeling of COVID-19 based on real information from Pakistan [8], designing the SEIR model of COVID-19 [9], the modeling of the memristor-based hyperchaotic circuit via nonsingular operator [10, 11], the thermostat model via Bernstein polynomials [12], the modeling of coronavirus by the Caputo operator [13], the optimal control of nonsingular tumor-immune surveillance [14], the SEIRA model [15], designing a bank data with fractal-fractional operators [16], the nonlinear modeling of the Navier system [17], analyzing an Atangana–Baleanu SEAIR model [18], compartmental disease modeling [19], and the references therein.

Lorenz developed the famous three-dimensional mathematical model to study the dynamics and behaviors of the atmosphere in 1963 [20–22]. Later on, in 1996, Lennart Stenflo modified the Lorenz system, adding a fourth parameter to consider the evolution of finite-amplitude acoustic-gravity waves in a rotating atmosphere, and developed a new four-dimensional Lorenz-Stenflo mathematical model. Since then, many studies have been conducted on the complex dynamical behaviors of the system by applying both integer and fractional derivatives. For further information, see [21] and the references therein.

Some of the studies conducted on the complex dynamics of Lorenz-Stenflo systems are reviewed as follows: Zhang et al. [22] investigated the qualitative properties of the higher integer-order Lorenz-Stenflo Chaotic system which appeared in mathematical physics. The authors proved that the higher-order Lorenz-Stenflo system is globally stable. A globally attractive set of Lorenz-Stenflo systems indicating the evolution of a finite amplitude of acoustic gravity waves contained in the rotating atmosphere is investigated by Zhang et al. [20]. Adaptive control synchronization and circuit implementation of the ordinary Lorenz-Stenflo system with an integer order was established by Yang and Wu [23]. An analysis on the dynamics of the Lorenz-Stenflo system via fractional operators along with sets of parameters is done through the spectra of the Lyapunov exponent type, bifurcation graphs, and 0-1 test.

The complex dynamics of Lorenz-Stenflo dynamical systems are analyzed with Lyapunov exponents, bifurcation diagrams, and the 0-1 test by using the Adomian decomposition numerical scheme and fractional-order representation of the model in the sense of the Riemann-Liouville integral [24]. A robust chaos suppression control is designed and applied via sliding mode control to the integer-order Lorenz-Stenflo system constrained to uncertainties and nonlinearities [25].

There is only limited literature devoted to studying the complex dynamics of the Lorenz-Stenflo systems by using integer-order derivatives and fractional-order derivatives. Moreover, to the best of the authors’ knowledge, no study is conducted on the qualitative analysis of the Lorenz-Stenflo mathematical model in the sense of the Caputo fractional derivative and by using a numerical scheme developed by Garrappa [26], which is the main focus of this study.

Among the several notions of fractional derivatives presented above, the Caputo fractional derivative gives the opportunity of including the classical initial conditions in a mathematical model and the Caputo fractional derivative of a constant is zero, which is not the case for instance in the Riemann-Liouville fractional derivative. Moreover, it has a MATLAB code that can obtain phase portraits and time series solutions of chaotic or hyperchaotic systems using a numerical scheme named the predictor-corrector method reported by Garrappa and applied in several studies on chaotic systems [27–29].

The manuscript is organized as follows: Section 2 involves the statement of the problem. In this section, the fractional representation of the Lorenz-Stenflo system (LS) is designed via the Caputo derivatives after recalling some critical preliminary definitions of fractional operators. Section 3 is devoted to analyzing the local stability of the mentioned LS system. In Section 4, the numerical scheme for the fractional-order LS system is developed based on a research reported by Garrappa [26]. The broader part of the manuscript is devoted to Lyapunov exponents, bifurcation, and chaos of the LS system in Section 5, followed by sensitivity analysis to initial conditions in Section 6. In Section 7, the development of the master and slave system to create a strong relationship between the systems using a coupling function is considered. Finally, a conclusion and a reference list are provided.

#### 2. Statement of the Problem

The Lorenz-Stenflo (LS) system of differential equations is given by [30, 31].

The variable is the intensity of motion of the fluid; and denote the horizontal and vertical direction temperature variations of the atmosphere. Parameters , , and are all positive real numbers and are the Prandtl number, is the rotation number, is the Rayleigh number, and is a geometric parameter. Parameters and depend on the material and geometrical properties of the layer of the fluid. Parameter is proportional to the difference in temperature.

We used fractional operators to find the hidden properties of the nature of solution of the LS system that are not observable via integer-order derivatives.

##### 2.1. The Formulation of the LS System via the Caputo Fractional Derivative

Here, we shortly give several definitions of fractional operators pertinent to our study.

*Definition 1 (see [2, 3]). *The Riemann-Liouville fractional integral for a continuous function is defined by

*Definition 2 (see [2, 3]). *The Riemann-Liouville fractional derivative for a continuous function is defined by

*Definition 3 (see [2, 3]). *The Caputo fractional derivative for a continuous function is defined by

This section develops the fractional order representation of the LS system (1). The Caputo fractional derivative is used because it can incorporate customary initial conditions in the model, unlike the Riemann-Liouville fractional derivative. Furthermore, the Caputo fractional derivative of a constant is zero, which is not the case in the Riemann-Liouville fractional derivative. Accordingly, the Caputo fractional representation of the LS system (1) is given by where with the initial value conditions and

The fractional-order derivative and the parameter values are and

The specific objectives of this study are to analyze the local stability of system (5), to investigate the chaotic behavior of system (5) via Lyapunov exponents, and to plot the bifurcation diagrams, attractors, and time series trajectories of the system by variation of fractional orders and four parameters of system (5). In addition, sensitivity to initial conditions and synchronization of the LS fractional system are considered.

#### 3. Stability Analysis of the LS System (5)

In this section of the manuscript, analysis of the local stability of the mentioned system (5) is conducted. There are several methods for performing local stability analysis of systems. Some of these methods are the Laplace transform techniques and the Matignon criterion. For the local stability analysis of the LS system (5), we used the Matignon method, because it is most commonly used in literature [30, 31]. The Matignon condition for fractional stability is given by where is the Jacobian matrix of the system, denotes the class of all eigenvalues of , and is the order of the LS system (5). The LS system (5) is said to be locally asymptotically stable whenever all the eigenvalues of satisfy the Matignon criterion. Firstly, we considered the equilibrium points of (5), given by , and constant values of the parameters given by , , , and .

The matrix at is given by

The eigenvalues of the matrix are and , and the corresponding arguments of the eigenvalues are 0 for and for the remaining eigenvalues. It is easy to conclude that is locally unstable, because the absolute value of the arguments does not satisfy the Matignon criteria (7). The presence of an eigenvalue with a positive real part, , is necessary for the SL system to exhibit a double-scroll attractor [31].

The general equilibrium points of the system in terms of the parameters are given by

Substituting the parameter values , , , and , the equilibrium points are given by and .

The Jacobian matrix corresponding to the equilibrium point is given by

The eigenvalues of are and The arguments of the eigenvalues are and . It can then be inferred that the Matignon criterion of local stability is satisfied for , since , On the other hand, the Matignon criterion is not satisfied for since only for That is, the equilibrium point is unstable for

The Jacobian matrix corresponding to the equilibrium point is given by

The eigenvalues and their corresponding arguments of Jacobian matrix are identical to those of Jacobian matrix . Hence, the conclusion in relation to the local stability of equilibrium point is the same as that in equilibrium point .

In summary, the equilibria of the fractional-order system (5) are unstable for the fractional order of and

#### 4. Numerical Solution

This section presents the numerical method used to obtain the phase portraits and time series solutions of the fractional-order LS system (5). To solve fractional-order systems, there are several numerical and analytical methods such as the predictor-corrector method, Homotopy method, and Adomian decomposition method. In this study, we use the predicator-corrector method, because it has a MATLAB code that can be used to obtain phase portraits and time series solutions of chaotic or hyperchaotic systems.

The numerical method in [26] is applied to our system as follows: by starting from the first equation of (5) and applying the Reimann-Liouville fractional integral given by (2), we obtain the following system

Applying the predictor-corrector method to (12), the integral terms are approximated and it becomes

where , and are defined by (5) and

Moreover, the predictors are given as follows:

By Garrappa [26], the discretization method applied in this study is stable and convergent. Thus, the advantages of this discretization method are that it is both stable and convergent, and of course, it has a MATLAB code that can be used for plotting the attractors of the LS system (5). Several other advantages of the method are described in [27].

As mentioned above, in this research, the predictor-corrector technique is utilized to obtain the attractors and time series trajectories of the LS system (5). Moreover, the Lyapunov exponents are obtained using an algorithm by Danca et al. [32]. The code is developed to determine all Lyapunov exponents of a class fractional-order system modeled by the Caputo derivative. The predictor-corrector Adams-Bashforth-Moulton numerical method is the underlying numerical method used in this code.

#### 5. Lyapunov Exponents, Bifurcation Diagrams, and Hyperchaotic Behavior of the LS System

This section is devoted to investigating the chaotic or hyperchaotic nature of the LS system (5). The magnitude of the chaos is quantified via the Lyapunov exponents (LE). Furthermore, bifurcation diagrams caused by the variation of the parameters of (5) are portrayed using the values of the parameters and the initial value condition presented in (5).

##### 5.1. Lyapunov Exponents for the Fractional Order

Based on the Danca algorithm mentioned above, some of the Lyapunov exponents corresponding to different fractional orders of the LS system (5) are shown in Table 1. The simulation is made to run for 300 s.

Based on the LE values in Table 1, it is followed that the dynamical system (5) possibly exhibits hyperchaotic behavior for , because in each column of Table 1, there are two positive LEs. The positive Lyapunov exponent ensures the sensitive dependence on the choice of initial conditions (local instability of the system in the state space), which is only one property of a chaotic system. Crudely, for a given dynamical system to be chaotic, it must have the following properties: sensitivity to initial values, topological transitivity, and also dense periodic orbits.

System (5) is dissipative for as the sum of the LEs in every column of Table 1 is negative. The dissipativity of the system ensures the existence of an attractor. From the experimental result in Table 1, it is observed that for values of , the system is not dissipative. For instance, the sum of LEs for the fractional order is 0.0211. Therefore, to get a dissipative hyperchaotic system in this study, we considered the fractional order .

The Kaplan-Yorke dimension denoted by of system (5) for can be calculated. For instance, for the fractional order the Kaplan-Yorke dimension is

The Kaplan-Yorke dimension corresponding to is

The Kaplan-Yorke dimension corresponding to is

The Kaplan-Yorke dimension corresponding to is

From the above s, it can be said that the Kaplan-Yorke dimension is decreased from 3.9975 to 3.6224, as the fractional order is increased from 0.857 to 1. Note that the Kaplan-York dimension determines the upper bound of the Hausdorff dimensions. In this study, the Hausdorff dimensions of the attractors corresponding to different fractional orders are nonintegers and system (5) has strange attractors with fractal structures. Moreover, in Table 1 and the Kaplan-Yorke dimensions calculated above, one can find that the significance of the hyperchaotic attractors decreases as the fractional order increases from 0.857 to 1. The loss of significance is described by decreasing the sum of the positive LEs in Table 1 and also by decreasing the dimensions as the fractional order increases from 0.857 to 1. The loss of significance is also observable from the bifurcation diagram of the fractional orders shown in Figure 1 in the next section.

##### 5.2. Bifurcation for the Fractional Order

The bifurcation diagram of the fractional order is obtained by varying its value in the interval with a time step of . The bifurcation diagram due to the variations of the order is illustrated in Figure 1. According to this figure, for values of , the given system exhibits oscillation with stability and the system excites a chaotic behavior for fractional values .

The simulation results corresponding to different fractional orders are portrayed in Figures 2–4 to verify the conclusion made in Figure 1.

**(a)**

**(b)**

**(c)**

**(a)**

**(b)**

**(c)**

**(d)**

**(e)**

It can be seen in Figures 2(a)–2(c) that the system exhibits oscillation with stability for the fractional order , which is in agreement with the conclusion in Figure 1. The time series solution of system (5) is shown for in Figure 3, where the solution trajectories have the property of oscillation and stability as claimed from the bifurcation diagram shown in Figure 1. Moreover, Figures 4 and 5 characterize that system (5) exhibits a hyperchaotic behavior at confirming what is predicted by the bifurcation diagram in Figure 1.

##### 5.3. Bifurcation for Parameter

To get the bifurcation diagram of parameter of the hyperchaotic system (5), the fractional order used is and , , , and with a time step of and the initial value condition is . The simulation is made to run for 100 s. The bifurcation graph due to the variation of is illustrated in Figure 6.

It can be inferred in Figure 6 that system (5) exhibits a hyperchaotic behavior for the parameters considered in the simulation for . The phase portraits of system (5) shown in Figures 7 and 8 for and verify the conclusion from the bifurcation diagram of Figure 6.

**(a)**

**(b)**

**(c)**

**(a)**

**(b)**

**(c)**

##### 5.4. Bifurcation for Parameter

For the bifurcation diagram of parameter of the hyperchaotic system (5), the fractional order used is , , and and the bifurcation parameter is made to vary in the interval [3.5, 4.5] with a time step of , with initial value condition , and the simulation is made to run for 100 s. The bifurcation diagram due to the variation of parameter is illustrated in Figure 9.

It is observable in Figure 5 that system (5) shows hyperchaotic dynamics for parameter in [3.5, 4.5]. Moreover, the phase portraits of system (5) for shown in Figure 10 verify the conclusion drawn in Figure 9.

**(a)**

**(b)**

##### 5.5. Bifurcation for Parameter

To obtain the bifurcation diagram due to parameter of the hyperchaotic system (5), the fractional order used is . The values of the remaining parameters are and , and the bifurcation parameter is made to vary in the interval [0.5, 1] with a time step of . The initial value condition is and the simulation is made to run for 100 s. The bifurcation graph due to the variation of is illustrated in Figure 11.

It is clear in Figure 11 that system (5) shows hyperchaotic dynamics for parameter of the interval [0.5, 1], the given values of the parameters, and fractional orders. Moreover, the attractors of system (5) for values of parameter in [0.5, 1] are approximately similar to the figures shown in Figure 10.

##### 5.6. Bifurcation for Parameter

To obtain the bifurcation diagram due to parameter of the hyperchaotic system (5), the fractional derivative, the initial condition, the time, and the values of the remaining parameters are the same as those used in Section 5.5. The bifurcation diagram due to the variation of parameter is illustrated in Figure 12.

It is clear in Figure 12 that system (5) shows the hyperchaotic dynamics for parameter in [9.5,10.5], the given parameters, and the fractional orders. Moreover, the attractors of system (5) for different values of are shown in Figure 13, verifying the hyperchaotic nature of system (5).

**(a)**

**(b)**

**(c)**

**(d)**

#### 6. Sensitivity to Initial Conditions

This section examines the effect of different initial values on the dynamics of the hyperchaotic system (5). Since sensitivity to initial conditions is one of the important properties of chaotic or hyperchaotic systems, it is necessary to examine this property for the LS system (5) using different initial conditions. In this study, to examine sensitivity to initial conditions, we used the parameter values indicated in (5) and . The initial condition considered is . The different values of are shown in Figure 14.

**(a)**

**(b)**

It can be observed in Figure 14 that all the trajectories coincided for about the first 5.1 seconds and they were divided into two for about 10.1 seconds. At about 15.1 seconds, they divided into three, diverging from each other. The values of of two trajectories that overlapped for the first 15.1 seconds are closer to each other than the value of of the third trajectory.

#### 7. Synchronization of the Lorenz-Stenflo System

In this part of the study, we developed a master-slave system of the dynamic system (5). Firstly, two identical copies of system (5) are associated via a coupling function. Then, we performed a simulation of the coupled systems using different initial conditions. Finally, the error dynamics showed asymptotic stability and the phase portraits of the master and slave system showed a strong correlation.

Let the master hyperchaotic model be given by and the slave system be given by where the coupling function is is the Jacobian matrix of the system, is a Hermitian matrix of appropriate size, and . Moreover, the coupling function is given by

where a Hermitian matrix is chosen to be

It can be shown that all the real parts of the eigenvalues of are negative. Thus, the Matignon criterion is satisfied; is a Hermitian matrix.

Now, considering the coupling function (23), the slave system becomes

Defining the error as the error dynamics is shown in

To verify if the master and slave systems in equations (21) and (25) are correlated, the time series trajectories of the two systems, including the graphs of the dynamics of the error, are illustrated in Figures 15–19. The values of the parameters used are together with the fractional order Moreover, the initial value conditions used for the master and the slave system are and , respectively.

It can be observed in Figure 15 that there is a robust correlation between the master and slave systems because the error graphs converge to zero approximately in the first 2.5 seconds, as shown in Figure 15.

#### 8. Conclusion

This research study implemented a qualitative analysis on the four-dimensional Lorenz-Stenflo mathematical model in the sense of the Caputo operator. The Matignon local stability criterion showed that the equilibrium points of the LS system are locally unstable, implying that the LS system led to hyperchaotic dynamical behavior. The scheme developed by Garrappa approximated the numerical solution, and the corresponding MATLAB code was used to obtain all the figures in the study. The authors believe that the qualitative analysis made in this manuscript, the displayed figures, and the synchronization systems developed have revealed complex dynamical behaviors of the Lorenz-Stenflo system that are not obtained earlier by other studies on the system. It is also worth mentioning that this study have been done by including other concepts of fractional derivatives and comparing the results to the results obtained due to the Caputo fractional derivative, to be treated in the future by the authors or interested researchers in the area.

#### Data Availability

No data were used to support this study.

#### Conflicts of Interest

The authors declare that they have no competing interests.

#### Authors’ Contributions

Chernet Tuge Deressa did the actualization, methodology, formal analysis, validation, investigation, and initial draft. Sina Etemad did the actualization, validation, methodology, formal analysis, investigation, and initial draft. Mohammed K. A. Kaabar did the actualization, methodology, formal analysis, validation, investigation, initial draft, and supervision of the original draft and editing. Shahram Rezapour did the actualization, validation, methodology, formal analysis, investigation, and initial draft. All authors read and approved the final version.

#### Acknowledgments

The second and fourth authors would like to thank Azarbaijan Shahid Madani University.