Dynamics of Rössler Prototype-4 System: Analytical and Numerical Investigation

Abstract

In this paper, the dynamics of a 3D autonomous dissipative nonlinear system of ODEs-Rössler prototype-4 system, was investigated. Using Lyapunov-Andronov theory, we obtain a new analytical formula for the first Lyapunov’s (focal) value at the boundary of stability of the corresponding equilibrium state. On the other hand, the global analysis reveals that the system may exhibit the phenomena of Shilnikov chaos. Further, it is shown via analytical calculations that the considered system can be presented in the form of a linear oscillator with one nonlinear automatic regulator. Finally, it is found that for some new combinations of parameters, the system demonstrates chaotic behavior and transition from chaos to regular behavior is realized through inverse period-doubling bifurcations.

1. Introduction

In the last thirty years, many authors have been aroused by the search for the mathematically simplest systems of various species that can exhibit chaos [1,2,3,4,5,6,7]. A good example is the book of J. Sprott [4] in which he discovered: (1) some new systems that are simpler than those previously know; (2) these new systems are otherwise more “elegant” on the part of the number of system parameters-their values, special symmetry and economy of notation.

The theory of dynamical systems which describes the qualitative changes in phase spaces at a variation of one or several parameters is called bifurcation theory [7,8,9,10,11]. Generally, all bifurcations can be separated (classed) into two main kinds: local and global. Although it is difficult to make a precise distinction between them, the dynamics of local bifurcations is determined by information contained in the terms of Taylor series or Poincarѐ map at a point, and the dynamics of global bifurcations requires information about the vector field along an entire (non-periodic) orbit [8,9,10]. With respect to dynamical systems theory [9,10,12], the Hopf bifurcation theorem [13], and some other elements of the local bifurcation theory (e.g., the first Lyapunov value—$L_1\left({\lambda }_0\right)$) are basic analytical instruments to investigate the transition between different dynamical states in a neighborhood of equilibrium states (fixed points). The (in) stability of these transitions may have essential consequences on the dynamics of the system. Practically, the complete investigation of complex phenomena in a dynamical system is impossible without a thorough study of both bifurcations-local and global.

An important element towards the understanding of the global behavior of a dynamical system of nonlinear differential equations is the analysis of the existence of homoclinic/heteroclinic trajectories (cycles) [14,15,16,17]. Note that a classical homoclinic cycle ${\mathsf{\Gamma }}_0$ is a loop that consists of a saddle equilibrium state, i.e., this as a trajectory which is bi-asymptotical to a saddle periodic orbit as time $t\to \pm \infty$. According to Peixoto’s theorem [18,19], homoclinic bifurcations are structurally unstable, i.e., they can be destroyed by small perturbations. Therefore, they are more difficult to identify than a local bifurcation, because knowledge of the global properties of the phase space trajectories is required. According to [11,20], the bifurcations in the dynamical systems can be divided into three types: (1) bifurcations not originating from the Morse-Smale class of systems; (2) bifurcations in a class of systems with non-trivial hyperbolicity; and (3) bifurcations associated with the transition from Morse-Smale to systems with non-trivial hyperbolicity. The last bifurcations are of particular interest because they can explain the mechanisms of emergence of chaotic properties in deterministic dynamical systems.

It is well-known that chaos cannot occur in dissipative two-dimensional systems (one degree of freedom) of ordinary differential equations (ODEs). Chaos requires systems with at least one and a half degrees of freedom. Such systems have so-called strange attractors (the trajectory winds around forever, never repeating) with non-integer dimension. In a three-dimensional continuous dissipative dynamical systems, the only possible spectra, and the attractors, are as follows: a strange attractor, a two-torus, a limit cycle, a fixed point [21,22,23]. Mathematical representation of a spatial order and chaos are saddle equilibria, saddle periodic movements or complex saddle invariant set. According to [24], around a saddle-focus equilibrium a systematic characterization of homoclinicity can be provided. It is well-known that, if the Shilnikov condition is satisfied, i.e., the saddle-focus index $\delta =\left|\mathrm{Re}\left(\raisebox{1ex}{${\chi }_2$}\!\left/ \!\raisebox{-1ex}{${\chi }_1$}\right.\right)\right|<1$ (where ${\chi }_1$ (or ${\chi }_s$) and ${\chi }_2$ (or ${\chi }_u$) are the leading stable and unstable eigenvalues (manifolds)), then an infinite number of non-periodic trajectories coexist in the vicinity of a homoclinic trajectory (orbit) ${\mathsf{\Gamma }}_0$ bi-asymptotic to the saddle-focus. For more details see [19,24,25]. In other words, if exactly one of the leading eigenvalues is complex (i.e., ${\chi }_s\in C$ and ${\chi }_u\in R$), then the homoclinic orbit is called unifocal and according to the Shilnikov condition close to it: (1) there exist horseshoes in every neighborhood of ${\mathsf{\Gamma }}_0$-when $\delta >1$; and (2) there is a neighborhood of ${\mathsf{\Gamma }}_0$ which contains no periodic orbits-when $\delta <1$. Here we note that Shilnikov condition $\delta \ne 1$ is also called non-resonance condition [15,16,26]. If we look at further conditions on the eigenvalues we can determine a second genericity requirement: $\delta \ne 2$. If $1<\delta <2$, then the linear flow at the equilibrium state(s) is contracting and the horseshoes are attracting, while for $\delta >2$ the linear flow at the origin is expanding [16].

The generic homoclinic orbit with complex principal eigenvalues (bifocal homoclinic orbit) has been studied by numerous researchers [25,27,28,29]. In a theorem they prove that horseshoes exist in every neighborhood of a generic bifocal homoclinic orbit.

For a long period of time the Lorenz and Rössler systems [30,31] were regarded as the simplest chaotic autonomous dissipative systems of ODEs. In this paper we consider the so-called Rössler prototype-4 system given by

$$\begin{array}{l}\dot{x}=-y-z,\\ \dot{y}=x,\\ \dot{z}=-bz+a\left(y-y^2\right),\end{array}$$

(1)

where $\left(x,y,z\right)\in R^3$ are the state variables [1]. Originally, Rössler in [1] has proposed a few abstract three-component chemical models (systems) which have chaotic solutions. Latter, Olsen and Degn [32,33] showed that these systems (including system (1)) are an example of chaos in experimental (real) chemical (enzyme) reaction systems. It is seen that system (1) has only two parameters, six terms and a single quadratic nonlinearity. To the best of our knowledge, the system was only numerically studied in [1,3,4,34]. It was shown that for $a=b=0.5$ the system (1) has chaotic behavior (see for more details Figure A3 in Appendix C), and for $a\in \left(0,1\right], b\in \left(0,1\right]$ has different dynamics-unbounded, periodic and chaotic solutions. In recent years, some interesting results about stability of periodic orbit and possible modifications of Rössler prototype-4 system (1) have emerged. In Garcia et al. [35] existence and stability of periodic orbits of system (1) are proved. The work described in Sprott and Linz [36] is another example where for some values of the parameters the dynamics of the system (1) is quasiperiodic and a trajectory lies on an invariant torus. In [37,38] the authors: (1) studied the condition and type of Hopf bifurcation in system (1) based on normal form; (2) investigated the amplitude control of limit cycle; (3) modified system (1) based on Chua’s diode nonlinearity; and (4) studied multistability and multiscroll generation of system (1). However, a detailed analysis of the transitory phenomena leading to irregular (chaotic) dynamics of system (1) is missing, which prompted our interest.

In this paper, using analytical and numerical tools, we investigate the dynamical behavior of system (1). The plan is as follows: in Section 2 and Section 3 we present analytical and numerical results concerning the system (1) for different values of bifurcation (control) parameters $a$ and $b$. In Section 4 we summarize our results.

2. Analytical Investigation

2.1. Local Analysis and First Lyapunov Value-$L_1$

In this section, we are going to consider the system (1) which presents an autonomous dynamical model. Clearly, the equilibrium (steady state) points of system (1) are: $O_1\left(0,0,0\right)$ and $O_2\left(0,1+\frac{b}{a},-1-\frac{b}{a}\right).$

The divergence of the flow (1) is:

$$D_3\left(t\right)=\frac{\partial \dot{x}}{\partial x}+\frac{\partial \dot{y}}{\partial y}+\frac{\partial \dot{z}}{\partial z}=D_3\left(t_0\right)\times e^{-bt}=-b,$$

(2)

where $D_3\left(t_0\right)$ is a volume element. The system (1) is dissipative and has an attractor, when $D_3\left(t\right)<0$, i.e., $b>0$.

The system (1) can be presented in the form of a (non)-linear oscillator with one automatic regulator. In general form we have

$$\ddot{y}=-\frac{\partial V}{\partial y},\hspace{1em}\dot{z}=-\epsilon \left[z-g\left(y\right)\right],$$

(3)

where $\epsilon$ is a small parameter, $g\left(y\right)$ is a nonlinear polynomial function and $V$ is the potential in the form $V=\frac{{\alpha }_0{y}^4}{2}-{\displaystyle \sum _{i=1}^2\left({\alpha }_i-{\beta }_{i}z\right)\frac{y^i}{i}}$ [39]. For system (1) we have: ${\alpha }_0={\alpha }_1={\beta }_2=0$, ${\alpha }_2=-1$ and ${\beta }_1=1$. Hence, for the system (1) we can write

$$\begin{array}{l}\ddot{y}=-\frac{\partial V}{\partial y}=-y-z,\\ \dot{z}=-b\left[z-kg\left(y\right)\right]=-bz+a\left(y-y^2\right),\end{array}$$

(4)

where $k=\raisebox{1ex}{$a$}\!\left/ \!\raisebox{-1ex}{$b$}\right.$, and the potential of the system is $V=\frac{y^2}{2}+yz$. The energy of the system (1) is $E=V+\frac{{\dot{y}}^2}{2}$. Here we note that Lorenz system also can be presented in the form (3) as ${\alpha }_0=1, {\alpha }_1={\alpha }_2={\beta }_1=0$, ${\beta }_2=1$ and $g\left(y\right)=y^2$ [39].

To determine analytically first Lyapunov value $\left(L_1\right)$ and Routh-Hurwitz conditions for stability of fixed points, we must accomplish some transformations of system (1). Hence, we obtain

$$\stackrel{⃛}{y}=-b\ddot{y}-\dot{y}-\left(a+b\right)y+a{y}^2.$$

(5)

Let us denote

$$w=y-y_0.$$

(6)

After substitution of (6) into (5), Equation (5) takes the form

$$\stackrel{⃛}{w}=-b\ddot{w}-\dot{w}+cw+a{w}^2,$$

(7)

where $c=a\left(2y_0-1\right)-b$ and $y_0$ is the equilibrium states $O_1$$\left(\overline{y}=0\right)$ or $O_2\left(\overline{y}=1+\frac{b}{a}\right)$.

On the other hand, let

$$\ddot{w}=y_1, \dot{w}=y_2, w=y_3.$$

(8)

Substituting (8) in (7), we obtain a system of three first-order differential equations in the form

$$\begin{array}{l}{\dot{y}}_1=-b{y}_1-y_2+c{y}_3+a{y}_3^2,\\ {\dot{y}}_2=y_1,\\ {\dot{y}}_3=y_2.\end{array}$$

(9)

The Routh-Hurwitz conditions for stability of $O_1$ and $O_2$ can be written in the form [12]:

$$\begin{array}{l}p=b>0,\\ q=1>0,\\ r=-c=-a\left(2y_0-1\right)+b>0,\\ R=pq-r=b+c=a\left(2y_0-1\right)>0.\end{array}$$

(10)

The notations $p, q, r$ and $R$ are taken from [12]. The characteristic equation of the system (9) (which is equivalent to system (1)) has the form

$${\chi }^3+p{\chi }^2+q\chi +r=0.$$

(11)

According to [9,12], the boundaries of the stability region for three-dimensional systems, are two surfaces: ${\mathsf{\Psi }}_1\left(R=0, p>0, q>0\right)$ and ${\mathsf{\Psi }}_2\left(r=0, p>0, q>0\right)$. On the first surface $R=0$, the characteristic equation (11) has a pair of purely imaginary roots as in this case an Andronov-Hopf (AH) bifurcation takes place, and at least one zero root on the surface $r=0$. It is well-known that Andronov-Hopf bifurcation can be: (i) supercritical (soft loss of stability) and (ii) subcritical (hard loss of stability) [9,12,13,23,40]. In the bifurcation points (where $R=0$), a positive first Lyapunov value represents a subcritical (hard or irreversible) bifurcation and determines that the system possesses unstable solutions but may fold back and exhibit unstable periodic solutions (oscillations) (unstable limit cycle) coexisting with stable steady state. In this case the boundary of stability is called “dangerous”. Inversely, a negative value for $L_1$ predicts a supercritical (soft or reversible) Andronov-Hopf bifurcation. Thus, the loss of stability takes place when stable periodic solutions (self-oscillations/stable limit cycle) emerge from a transition through the bifurcation point. Now, the boundary of stability is called ‘safe’.

It is seen that for first and second fixed points we have: $r^{\left(1\right)}=a+b, R^{\left(1\right)}=-a; r^{\left(2\right)}=-a-b; R^{\left(2\right)}=a+2b$. According to [12,23,40,41], we obtain the formula describing the first Lyapunov value for our system and calculate its value in the calculated bifurcation points (for more details, see Appendix A). Thus, for $L_1$, we write:

$$L_1\left({\lambda }_0\right)=\frac{\pi a^2}{2b\left(b^2+4\right){\Delta }_0^2}\left[bc\left(b^2+6\right)-3b^2-8\right],$$

(12)

where ${\Delta }_0=1-bc$, and ${\lambda }_0$ is defined as a value of $a$ and $b$ for which the relation $R=0$ takes place. It is clear that: (i) for $a=0$, $b>0$, then $R^{\left(1\right)}=0$ and $L_1=0$, i.e., the first fixed point $O_1$ is structurally unstable; (ii) for $a=-2b, b>0$, then $R^{\left(2\right)}=0$ and $L_1=-\frac{2\pi b}{\left(b^2+4\right){\left(1+b^2\right)}^2}\left(b^4+9b^2+8\right)<0$. Hence, according to [9,10,12], a soft stability loss takes place, i.e., in the case of transition through the AH boundary ($R=0$) from positive values to negative ones, a stable limit cycle (self-oscillations) emerges. Conversely, in the case of a transition from negative values to positive ones the stable limit cycle disappears, i.e., the self-oscillations cease. In the following Section 3, we demonstrate numerically this type of behavior.

2.2. Special Cases

2.2.1. Conservative Case

It is well-known that conservative systems play an important role in many mathematical and mechanical problems [11,42,43]. The equations of motion of such systems can be obtained from their Hamiltonian. A one-dimension system, with position coordinate $y$, momentum $x$ and Hamiltonian $H\left(x,y\right)$ has equations of motion

$$\dot{y}=\frac{\partial H\left(x,y\right)}{\partial x},\hspace{1em}\dot{x}=-\frac{\partial H\left(x,y\right)}{\partial y}.$$

(13)

Here $H\left(x,y\right)$ is a first integral because

$$\dot{x}\frac{\partial H}{\partial x}+\dot{y}\frac{\partial H}{\partial y}=\frac{d}{dt}H\equiv 0,$$

(14)

i.e., $H$ remains constant along the trajectories. In other words, $H$ is a conserved quantity (constant of the motion). It is important to note that the existence of a first integral is useful because of the connection between its level curves and the trajectories of the system.

For $a=b=0$, the system (1) reduces to a planar $\left(R^2\right)$ Hamiltonian system with Hamiltonian (first integral)

$$H\left(x,y\right)=\frac{x^2}{2}+\frac{y^2}{2}-Cy=h,$$

(15)

where $C$ is a real constant. This integrable system has a fixed point $O\left(0,-C\right)$ of center type, i.e., the trajectories in a neighborhood of fixed point are periodic (if $C=0$), or of saddle type (if $C\ne 0)$. Note here that, the following center forms exist: linear type center, nilpotent center and degenerate center. If the fixed point is of a focus or a center type, we say that it is monodromic singular point [44]. When $C=1$, the equilibrium state $O$ is a saddle at the level $H=\frac{3}{2}$. It is easy to see that, in this case, the general solution reeds

$$\begin{array}{l}x\left(t\right)=c_1\cos \left(t\right)-c_2\sin \left(t\right),\\ y\left(t\right)=c_2\cos \left(t\right)+c_1\sin \left(t\right)-c_3,\\ z\left(t\right)=c_3,\end{array}$$

(16)

where $c_1, c_2,$ and $c_3$ are real constants.

2.2.2. Dissipative Case

Let us view the situation when $a=0$ and $b>0$. Then for system (1) we have

$$\begin{array}{l}\dot{x}=-y-z,\\ \dot{y}=x,\\ \dot{z}=-bz,\end{array}$$

(17)

i.e., we obtain a linear $3D$ autonomous system. It is seen that the last equation into (17) contains only the variable $z$. Thus, after a direct calculation, we derive from the third equation of system (17) that for $z\ne 0$,

$$\left|z\right|=c_3{e}^{-bt},$$

(18)

where $\left|-bt\right|<+\infty$. On the other hand, the differential equation

$$\frac{dx}{dy}=\frac{-y-z}{x},\hspace{1em}x\ne 0$$

(19)

has solutions on the plane

$$\frac{x^2}{2}+\frac{y^2}{2}+y{c}_3{e}^{-bt}=h_1,$$

(20)

where $h_1$ is a real constant. Note that the module of $z$ is strictly bounded below by $b>0$ and $t\to +\infty$, i.e., $z\left(t\right)$ satisfies $\left|z\right|\to 0$.

A noteworthy, easily derivable analytical solution is this:

$$\begin{array}{l}x\left(t\right)=c_1\cos \left(t\right)-c_2\sin \left(t\right)-c_3\frac{1}{1+b^2}\left(b\cos \left(t\right)+\sin \left(t\right)-b\exp \left(-bt\right)\right),\\ y\left(t\right)=c_2\cos \left(t\right)+c_1\sin \left(t\right)-c_3\frac{1}{1+b^2}\left(b\sin \left(t\right)-\cos \left(t\right)+b\exp \left(-bt\right)\right),\\ z\left(t\right)=c_3\exp \left(-bt\right),\end{array}$$

(21)

where similar to the previous case $c_1, c_2,$ and $c_3$ are real constants.

2.3. Global Analysis

These bifurcations for which it is insufficient to consider small local surroundings near the equilibrium states or the limit cycles are called global. They lead to qualitative change in the stable and the unstable manifolds of the states of equilibrium and/or the limit cycles. One of these situations is the occurrence of a loop of the separatrices of one and the same saddle (steady state) [14,45]. For $a\in \left[-1.51,1\right]$ and $b\in \left(0.29,0.8\right]$ we obtain [34] that the two fixed points $O_1$ and $O_2$ at bifurcation parameters $a$ and/or $b$ are of saddle-focus type: (i) unstable focus-negative real eigenvalue and complex eigenvalues with positive real part (see Figure 1 left panel); (ii) stable focus-positive real eigenvalue and complex eigenvalues with negative real part (see Figure 1 right panel). Note that for some subintervals the fixed point $O_2$ can be of stable focus type. These two fixed points can be included in homoclinic/heteroclinic structures of Shilnikov type, where their stable and unstable invariant manifolds ($W^s$ and $W^u$), are meeting each other in a most intricate manner. The structure of phase space near homoclinic loop essentially depends on the saddle index $\delta$, i.e., when $\delta >1$, a simple dynamics takes place, and when $\delta <1$, a complex dynamics can be seen.

It is well-known that the existence of a homoclinic/heteroclinic cycle is one of the common manners of the formation or disappearance of a limit cycle, when the limit cycle emanates from, or approaches the heteroclinic cycle as a singular limit, respectively [14,45,46]. There are known situations in which a unique limit cycle is born, and certain criteria can be used to determine if this limit cycle must be stable or unstable. As we mentioned above, the system (1) has steady states of saddle-focus type. Therefore, in our case the Shilnikov theory for analyzing of bifurcations of homoclinic/heteroclinic trajectories and existence of complex dynamics (chaos) in flow in $R^3$ can be used [24]. In this connection, in the Appendix B, we explain shortly the general results obtained by Shilnikov. On the other hand, historically, the four-dimensional case was considered in 1967 [47], and the general case-three years later in [25].

For the system (1), if: (1) $a=-1.5, b=0.6$, (2) $a=-1, b=0.35$ and (3) $a=b=0.5$ the equilibriums (fixed points) and their eigenvalues are given by:

Case 1.

$O_1\left(0,0,0\right),$ then $\left({\chi }_1,{\chi }_2,{\chi }_3\right)=\left(0.5508,-0.5754\pm 1.1414i\right)$, $0.5<\delta =0.957<1$;

$O_2\left(0,0.6,-0.6\right),$ then $\left({\chi }_1,{\chi }_2,{\chi }_3\right)=\left(-0.7855,0.0928\pm 1.0664i\right)$, $\delta =0.1181<0.5$.

Case 2.

$O_1\left(0,0,0\right),$ then $\left({\chi }_1,{\chi }_2,{\chi }_3\right)=\left(0.4694,-0.4097\pm 1.1031i\right)$, $\delta =1.146>1$;

$O_2\left(0,0.65,-0.65\right),$ then $\left({\chi }_1,{\chi }_2,{\chi }_3\right)=\left(-0.5754,0.1127\pm 1.0569i\right)$, $\delta =0.196<0.5$.

Case 3.

$O_1\left(0,0,0\right),$ then $\left({\chi }_1,{\chi }_2,{\chi }_3\right)=\left(-0.8038,0.1519\pm 1.105i\right)$, $\delta =0.189<0.5$;

$O_2\left(0,2,-2\right),$ then $\left({\chi }_1,{\chi }_2,{\chi }_3\right)=\left(0.6015,-0.5507\pm 1.1659i\right)$, $\delta =1.092>1$.

Since the complex exponents ${\chi }_{2,3}$ for the saddle-foci $O_2$ (in Cases 1 and 2) and $O_1$ (in Case 3) are nearest to the imaginary axis, the homoclinic loop implies the emergence of infinitely many periodic orbits of saddle type. Moreover, since the second saddle value ${\sigma }_1$ (see Appendix B) is negative, then near the homoclinic loop may exist stable periodic orbits along with saddle ones. Note that these stable orbits are with long periods and have weak attraction basins, i.e., they are practically invisible in numerical experiments.

Interestingly, when the equilibrium states $O_1$ (in Cases 1 and 2) and $O_2$ (in Case 3) have exponents ($\dim W^u=1$, $\dim W^s=2$), the second saddle value ${\sigma }_1$ is also negative and therefore in a small neighborhood of the homoclinic loop there are stable periodic orbits. Moreover, the condition $0.5<\delta <1$ in this case is necessary to exclude the transition to complex dynamics via codimension-two bifurcations of homoclinic loop(s) [9,14,48]. Note here that these bifurcations of homoclinic loops are essential for bifurcation behavior in the Lorenz system [30].

Finally, we will remark the following fact: at the non-trivial equilibrium $O_2$, the homoclinic bifurcations always started/finished from the super-critical Andronov-Hopf bifurcation (the first Lyapunov value is negative).

3. Numerical Investigation

In the considered case here, some of the known analytical techniques and results are not applicable [16,49,50], and we are forced to use numerical calculations and specific features of system (1). Hence, in this section we find some new results for its bifurcation behavior and route to chaos. In order to compare the analytical predictions with numerical results, the governing equations of system (1) were solved numerically using MATLAB (The MathWorks, Inc., Natick, MA, USA). Consider the numerical technique (based on the initial value problem) used here, we note that it does not allow to detect and numerically trace periodic orbits and predict their bifurcations. In this case, a possibility is to use the systematic numerical approach presented in [51].

Figure 2 shows bifurcation diagrams for the system (1): (a) and (c) values of $y$ coordinate versus $b$; and (b) and (d) values of $z$ coordinate are plotted against $b$ regarded as a continuously varying control (bifurcation) parameter, when the parameter $a=-1.51$ is fixed. It is seen that for $b\in \left[0.495,0.524\right]$ the system is chaotic. On a further increase of the bifurcation parameter $b$, the system (1) exhibits inverse period-doubling bifurcations leading to a periodic motion.

A more detailed investigation of the system (1) in the region $b\in \left[0.52,0.6\right]$ is shown in Figure 2c,d, i.e., when the bifurcation parameter $b$ varies in narrower ranges. Practically, we observe in this interval a cascade of inverse period-doubling bifurcations when the fixed points $O_1$ and $O_2$ are of saddle-focus type, i.e., $O_1$ possesses a 2-dimensional stable manifold and a 1-dimensional unstable manifold, and $O_2$ possesses a 2-dimensional unstable manifold and a 1-dimensional stable manifold. It is interesting to note that $O_2$ becomes stable focus ($\dim W^s=3, \dim W^u=0$) after $b=0.75$ and remains stable till the end of the interval $b\in \left(0.75,0.8\right]$. In other words, as a result of the evidences shown in Figure 2a,b, we can conclude that a stable limit cycle (small amplitude attenuation oscillations) with period one occurs when $O_2$ becomes stable focus-see for details Figure 3. Denote here that for $b\in \left[0.77,0.8\right]$ the system (1) has only stable solutions-the self-oscillations cease. These numerical results (calculations) are in exact accordance with the corresponding analytical results obtained for first Lyapunov value in previous Section 2.1.

In Figure 4 the bifurcation diagrams of the system (1) when $a=-1$ are shown. It can be seen that at $b\in \left[0.315,0.32\right]$ chaotic solution occurs. In analogy with the previous case, the system passes from chaos to regular motion after inverse period-doubling bifurcations. In this case the fixed points $O_1$ and $O_2$ are also of saddle-focus type, i.e., $O_1$ possesses a 2-dimensional stable manifold and a 1-dimensional unstable manifold, and $O_2$ possesses a 2-dimensional unstable manifold and a 1-dimensional stable manifold. Here $O_2$ becomes stable focus after $b=0.5$ and remains stable till the end of the interval $b\in \left(0.75,0.65\right]$. It is interesting to note that here the system (1) has regular solutions at the beginning and in the end of the interval for the control parameter $b$. Comparing Figure 2 and Figure 4 we conclude that in the case, when $a=-1.51$ the chaotic zone is longer than those obtained in Figure 4 for $a=-1$. It is interesting to note that inverse period-doubling bifurcations can be seen in many experimental systems, for example triple physical pendulum [52]. Thus, it is possible to detect the hysteresis phenomena, i.e., it is expected that bifurcation diagrams presented in Figure 4 will be different with increase and decrease of the parameter $b$.

Figure 5 shows bifurcation diagrams for the system (1): (a) and (c) values of $y$ coordinate versus $a$; and (b) and (d) values of $z$ coordinate are plotted against $a$ regarded as a continuously varying control (bifurcation) parameter, when the parameter $b=0.5$ is fixed. Similar to the previous case for $a=-1.51$ (see Figure 2), initially the system (1) is chaotic. On a further increase of the bifurcation parameter $a$, the system (1) exhibits inverse period-doubling bifurcations leading to a periodic motion-stable limit cycle according to analytical results obtained for first Lyapunov value in Section 2.1. The white zones, seen in Figure 5a,b, correspond to very fast inverse period-doubling bifurcations-see Figure 5c,d.

In Figure 6, the Andronov-Hopf bifurcation curve (boundary of stability loss $R=0$) of the fixed point $O_2$ is shown for different values of the bifurcation parameters $a$ and $b$. Following [9,10,12,23] and results obtained in previous Section 2.1 for first Lyapunov value, we conclude that this boundary of stability is “safe”, and a soft stability loss occurs. Therefore, the transformation into/from chaos is related to a soft loss of stability, i.e., the fixed point $O_2$ is a stable weak focus on the stability boundary.

In Figure 7 we show the bifurcation diagrams of the system (1) when $a\in \left[0.1,0.51\right]$ and $b=0.5$. In this case the system (1) demonstrates a period-doubling route to chaos-a steady state loses its stability and simultaneously a new orbit of doubled period is created. Here the fixed points $O_1$ and $O_2$ are of saddle-focus type in all interval for bifurcation parameter $a$, i.e., $O_1$ possesses a 2-dimensional unstable manifold and a 1-dimensional stable manifold, and $O_2$ possesses a 2-dimensional stable manifold and a 1-dimensional unstable manifold. It is interesting to note that when $a=0.51$ and $b=0.5$, the system has pseudo-chaotic (order in chaos) behavior for very short intervals of time-see for details Figure A4.

Finally, discussing the results shown in Figure 5c,d, Figure 7c,d, it is seen that an apparent sudden collapse/appearance in the chaotic attractor size occurs at a value of control parameter $a\approx -1.464$ (or $a\approx 0.465$). According to [53,54,55,56], such a sudden qualitative change in a chaotic attractor is called interior crisis. It is in a good agreement with our analytical results obtained in Section 2, the Theorem explained in Appendix B and the numerical results in the Appendix C—see Figure A1 and Figure A2.

4. Conclusions

The paper presents a study of the dynamic behavior of the so-called Rössler prototype-4 system, using analytical and numerical tools. The governing differential equations of system (1) were solved numerically using MATLAB (The MathWorks, Inc., Natick, MA, USA). For all simulations the initial conditions were $x\left(0\right)=y\left(0\right)=z\left(0\right)=0.1$. The analytical and numerical results lead to the following comments: (1) the system (1) has two fixed points of saddle-focus type, and therefore homoclinic/heteroclinic structures of Shilnikov type take place; (2) for values of the coefficients $a$ and $b$ different from these in [4], the system (1) has chaotic solutions; (3) the original system (1) can be presented in the form of a linear oscillator with one nonlinear automatic regulator; (4) the route chaos starts/finishes from a soft (reversible) loss of stability; (5) for $a=0.51$ and $b=0.5$, the system has pseudo-chaotic (order in chaos) behavior for very short intervals.

Generalizing the numerical results shown in Figure 2, Figure 3, Figure 4, Figure 5, Figure 6 and Figure 7 and those in Appendix C for Rössler prototype-4 system (1), we can conclude that: (i) for values of bifurcation parameter $b\in \left[0.495,0.524\right]$ (for $a=-1.51$), the system is in a chaotic regime. In result of inverse period-doubling bifurcations, for values of $b>0.54$, the system has only regular solutions with different period; (ii) for values $b\in \left[0.315,0.32\right]$ (when $a=-1$), the system (1) is in a chaotic state, as for $b\in \left(0.32,0.65\right]$ period-doubling bifurcations take place. For $b=0.5$, a stable limit cycle emerges whose period is one; iii) for values of bifurcation parameter $a\in \left[-1.51,-1.4\right]$ or $a\in \left[0.41,0.51\right]$ (for $b=0.5$), the system (1) has chaotic/regular solutions. In these intervals for $a$, very fast inverse period-doubling bifurcations and so-called interior crisis take place.