Coexisting Hidden and self-excited attractors in an economic system of integer or fractional order

According to Shilnikov criteria, chaos emergence requires at least one unstable equilibrium [Shilnikov(1965)]. So, if for some values of the key parameter the equilibria are unstable and one tries to generate attractors with initial conditions near these equilibria, one obtains chaotic attractors and could deduce that the system evolves chaotically. On the other side, this conclusion could be incomplete, or even false, since for initial conditions taken from an attraction basin which does not intersect with any small neighborhoods of unstable equilibria, the underlying trajectories could lead to some other attractor that might be a regular motion, but not a chaotic motion. Recently this problem has been managed by introducing the term of hidden attractor (see e.g. [Leonov & Kuznetsov(2013)], one of the first works on this subject, or [Kuznetsov(2015)]). If the attraction basin intersects with any open neighborhoods of an equilibrium the attractor is called self-excited, otherwise, it is called hidden. Usually, the self-excited attractors are numerically derived from unstable equilibria, while hidden attractors are difficult to be localized, because their attraction basins have no relation with small neighborhoods of any equilibria. Hidden attractors can be found in a nonlinear system with one or more stable equilibria [Molaie et al.(2013), Wang et al.(2017), Cang et al.(2019), Qigui et al.(2019)], with a line of equilibria (infinite equilibria) [Jafari & Sprott(2013)], with both unstable equilibria and stable equilibria [Danca(2017)], with coexistence of various attractors [Zhou et al.(2018)], or even without equilibria [Wei(2011)]. Note that there is a lack of correlation between hidden attractors and equilibria. Moreover, the precise localization of hidden attractors seems to be an intractable issue and, to our best knowledge, there is still no general analytical way, but only numerically with luck, to find hidden attractors (for the particular case of Chua’s circuit, see [Leonov et al.(2011)]).

To understand better the importance of hidden attractors consider the crash of aircraft YF-22 Boeing in April 1992. The analysis of aircrafts and launchers control systems with saturation regards the linear stability and also hidden oscillations which might occur [Andrievsky et al.(2013)]. The difficulties of rigorous analysis and design of nonlinear control systems with saturation related to this case are presented in [Lauvdal et al.(1997)].

Therefore, the identification of hidden attractors in a given system can help drive the system along the desired attractor, which can be a self-excited attractor or a hidden attractor.

The fractional order (FO) calculus is as old as the integer-order (IO) one, but its application was exclusive in mathematics. In the latter years FO systems where proved to describe, generally more accurately and in compact expressions, the behavior of real dynamical systems, compared to the IO models. This happens taking into account the nonlocal characteristic of “infinite memory”, i.e. the evolution of the system depends at each moment on the entire previous history. Basic aspects of the theory of FO can be found in [Oldham & Spanier(2006)]. While the definition of FO derivative for continuous-time real functions has been formulated by Liouville, Grunwald, Letnikov, and Riemann, in the late 19th century, the first definition of a fractional difference operator was proposed in 1974 [Diaz & Olser(1974)]. In this paper the FO derivative is considered in the sense of Caputo, because it allows the choosing initial conditions for the IO systems. One of the most utilized numerical methods to integrate continuous systems of FO is the predictor-corrector Adams-Bashforth-Moulton predictor-corrector method for fractional-order differential equations [Diethelm et al. (2002)].

There are many continuous but also discrete chaotic systems of IO, such as the supply and demand systems which are among the oldest and simplest economic discrete models, having complicated dynamics [Lorenz(1993)]. It is proved that the economical models exhibit generally unstable steady states and also fluctuations if the income distribution varies sufficiently and if shareholders save more than workers (see e.g. [Miloslav(2001)] and references therein).

In this paper a continuous economic system is considered, and proved numerically that it presents coexisting attractors for both IO and FO variants.

In [Miloslav(2001)] the question of whether using and/or un-using of foreign investments can change the qualitative properties is addressed for a growth path of a proposed dynamical economic system with three endogenous variables and only one non-linear term, modeled by 3-dimensional autonomous differential equations, while in [Pribylova(2009)] the system is further analyzed via bifurcations routes including supercritical and subcritical Hopf bifurcation, and generalized Hopf bifurcation as well. It is proved that the existing cycle exhibits period-doubling bifurcation as a route to chaos.

The economic system considered in this paper, with a foreign financing, has the following form

where the variables $x_{1}$ represents the savings of households, $x_{2}$ the Gross Domestic Product (GDP), $x_{3}$ the foreign capital inflow, and the variables $a$ the variation of the marginal propensity to savings, $b$ the ratio of capitalized profit, $c$ the potential GDP, $d=1/\nu$, with $\nu$ the capital/output ratio, $e$ the capital inflow/savings ratio and $f$ the debt refund/output ratio [Lorenz(1993)]. As explained in [Lorenz(1993)], if $x_{2}<c$, the activities are not constrained and high profits are derived from the sectors where the markets are not yet saturated. The single non-linearity represents the quality of a capitalization of the profits. If $x_{2}>c$, then the inflation is possible. A lack of new investment opportunities modifies it by savings.

Note that the system can also be viewed as an extension of the van der Pol s system on the plane $(x_{1},x_{2})$, i.e. when economically the system has no foreign investment ($x_{3}=0$)

which, as a two-dimensional autonomous system, cannot display chaotic behavior.

In this paper one considers the 3-dimensional system (1) with $a$ the bifurcation parameter and $b=0.01$, $c=1$, $d=0.031847$, $e=0.19$, $f=0.25$ (there are other interesting possible parameters choices¹¹1Parameter $d$ has been chosen as in [Miloslav(2001)], with six decimals, because rounding to fewer decimals affects significantly the results.)

While in [Miloslav(2001)] coexisting stable cycles and equilibria and also the existence of a chaotic attractor are studied analytically, in this paper it is shown numerically that for $a<f$ the system presents richer dynamics, including coexisting hidden and self-excited chaotic attractors and stable cycles. This behavior could be extremely interesting from an economical point of view.

Because the right hand side of the system (1) is an odd function, the system presents symmetries in the bifurcation planes, and also in the phase space.

The equilibria, which are collinear (see Fig. 4), have the following expression

Note that $X_{1,2}^{*}$ with the considered 4 decimals, are approximations of the exact equilibria $X_{1,2}^{*}$.

The stability of equilibria and the occurrence of the Hopf bifurcation for $a\geq 0.25$ are studied in [Miloslav(2001)].

In this paper, the bifurcation parameter $a$ is considered, $0\leq a\leq 0.2$, where all equilibria are unstable and the system presents most interesting dynamics.

The following result, presented in [Miloslav(2001), Pribylova(2009)], is slightly improved and re-proved numerically here, to be useful for finding hidden attractors later.

Fig. 1 (c) presents a chaotic attractor for $a=0.001$. As can be seen, before reaching the chaotic attractor, the trajectory connects the saddles $X_{0}^{*}$ and $X_{1}^{*}$, i.e. it is a heteroclinic connection.

Hereafter attractors are denoted as follows:

• Self-Excited Cycles: $SEC$;

• Self-Excited Chaotic Attractors: $SECH$;

• Hidden Chaotic Attractors: $H$;

• Hidden Cycles: $HC$.

The numerical method utilized for the integration of the system (3) of IO is the Matlab variant of RK method, $ode45$.

The numerical trajectories for different attractors are represented overplotted in the phase space, as time series and the local Lyapunov exponents are determined. Generally, transients are removed. Also, excluding the Perron effect (when a negative largest Lyapunov exponent does not, in general, indicate stability, and also when a positive largest Lyapunov exponent does not, in general, indicate chaos [Perron(1930), Leonov & Kuznetsov(2007)]), the positiveness of the maximal Lyapunov exponent (MLE) is considered as an evidence for chaos. In this paper the zero MLE is considered as having a maximum value of order $10^{-5}-10^{-4}$. If $MLE$ is zero, then, after transients removed, the trajectory is considered as regular (stable cycle), while if $MLE>0$ the trajectory is considered chaotic.

The bifurcation diagram of the maximum component $x_{1}$, for $a\in[0,0.2)$, is presented in Fig. 2 (a).

As can be seen, for all values of $a$, the underlying attractors are symmetric. Due to the symmetry, the coexistence of attractors features this system. Thus, depending on the initial condition $x_{0}$, there exists one or even two pairs of symmetric attractors (blue-red and black-yellow in the bifurcation diagram). Beside the standard period-doubling cascade which leads to chaos, the zoomed region in Fig. 2 (b) reveals an interesting coexisting window for $a\in[a_{1},a_{2}]$, with $a_{1}=0.0485$ and $a_{2}=0.0524$, where exterior and interior attractors crises, thin periodic windows interleaved between large coexisting chaotic-periodic windows, and also several miniature bifurcation diagrams within a large chaotic window can be seen. Moreover, at $a=a_{2}$, there exist both exterior and interior crisis.

However the most interesting characteristic of the window presented in Fig. 2 (b) represents the coexistence of self-excited chaotic attractors with hidden chaotic attractors.

Because there exists no general analytical way to find the attraction basin of hidden attractors of a given nonlinear system, in the case when the considered system admits unstable equilibria, an acceptable way is to search randomly initial conditions in neighborhoods of all these equilibria. This can be done either in planar (usually horizontal) sections through each equilibrium or/and analyzing three-dimensional (usually spherical) neighborhoods centered on equilibria.

For the considered system (3), one considers lattices $B_{0,1}=[-5,5]\times[-5,5]$ of $400\times 400$ points $(x_{1},x_{2})$ and zoomed regions, as horizontal sections through $X_{0}^{*}$ (with $x_{3}=0$) and $X_{1}^{*}$ (with $x_{3}$ the third coordinate of $X_{1}^{*}$), respectively, and also spherical neighborhoods of unstable equilibria.

Because of the symmetry, the case of the equilibrium $X_{2}^{*}$ is similar to $X_{1}^{*}$. Therefore only $X_{0}^{*}$ and $X_{1}^{*}$ are analyzed.

Probably the most interesting characteristic of this system is the coexistence of attractors in the space of the fractional-order systems in the sense that for some fractional order $q$ and fixed parameter $a$, there exist two symmetric pairs of coexisting chaotic and stable cycles attractors.

The commensurate FO variant of the system is

where $D_{*}^{q}$ denotes the Caputo differential operator of order $q\in(0,1)$ with starting point 0 (see e.g. [Podlubny(1999), Diethelm et al. (2002)])

where $\lceil\cdot\rceil$ denotes the ceiling function that rounds up to the next integer and $D^{\lceil q\rceil}$ represents the standard differential operator of order $\lceil q\rceil\in\mathbb{N}$. Due to the use of Caputo’s operator, the initial conditions can be taken as for the integer-order case.

Consider next $a=0.05$. The bifurcation diagram with respect the fractional order $q$ is presented in Fig. 9 (a). As can be seen, at the right of the bifurcation diagram, there exists a FO coexisting window where two symmetric pairs of attractors coexist.

Note that, as for continuous integer-order systems and contrary to FO difference systems, where the chaos strongness decreases for increasing $q$ (see e.g. [Danca(2020)]), for this system of FO, the chaos strongness increases with $q$ tending to 1.

Because the inherent time history of the ABM method, a deep numerical analysis to detect hidden attractors is tedious and time consuming. Therefore, the quality of the obtained attractors is determined by several empirical tests.

Now, the cycles have $MLE=5e-4$ while for the chaotic attractors $MLE=0.0056$.

Therefore, for a given parameter value $a$, the FO system (4) admits at least one FO coexisting window in the space of fractional order $q$, with two symmetric pairs of coexisting attractors. Other potential FO-coexisting windows, not examined here, are at the bifurcation points, e.q. at about $a=0.985$.

This important characteristic indicates that for some given fractional order $q$ and parameter $a$, there could be two different systems.

Choosing a suitable level of precision for computations is a real challenging task. The numerical integrator for (3) in the IO case is the matlab ode45 with implicit absolute and relative error tolerances. For the particular case of $a=0.0509$ (Subsection 3.3), errors have been reduced by using $options=odeset(^{\prime}RelTol^{\prime},1e-8)$.

The maximal time used interval is $T_{max}=5000-8000$, depending on the lengths of discarded transients. Due to inherent numerical errors, larger values for $T_{max}$ could lead to invalid results.

To integrate numerically the system (4) of FO, the ABM method [Diethelm et al. (2002)] is utilized.

The algorithms used to determine the finite time MLEs are the Benettin-Wolf algorithm for the IO case, and the algorithm presented in [Danca & Kuznetsov(2018)] for the FO case.

The attractors data are presented in Table 1. Beside the initial conditions presented in Table 1, points $(\pm 1,\pm 1,\pm 1)$ and $(\pm 1,\pm 1,\mp 1)$ can be used as initial conditions.

Note that, because of the relatively large number of decimals of some system coefficients and also due the inherent numerical errors, the 4-digits coordinates of $X_{1,2}^{*}$ (as approximations of the exact $X_{1,2}^{*}$) could also be taken as initial conditions, for cases with neighborhoods of $X_{1,2}^{*}$.

In this paper several type of coexisting attractors of an economic system of IO and FO are presented. Considering one of the several parameters as bifurcation parameter, the bifurcation diagram revealed windows where hidden and self-excited attractors coexist.

First, it is shown numerically that equilibria are unstable for both IO and FO cases, after which the exploration within neighborhoods of these points is investigated. The experiments, realized with the matlab routine ode45 for the IO case and with the ABM method for the FO case, show that the considered attractors fit the definition of hidden or self-excited attractors.

The characteristics of chaotic or regular attractors are revealed, after transients removed, by phase portraits, time series and maximal Lyapunov exponents.

For the bifurcation diagram in Fig. 1, the four used initial conditions are $(\pm 1,\pm 1,\pm 1)$ and $\pm 1,\pm 1,\mp 1)$, respectively. These initial conditions are proved to be valid, besides those mentioned in the text, for all considered cases (Table 1). Nevertheless, to note that generally this is a way to draw bifurcation diagrams are made for fixed initial conditions, it is not a rigorous method, since the attraction basins be modified while the bifurcation parameter varies. Thus some searched attraction basins, such as hidden attractors, might modify once the parameter is varied and, in this way, hidden attractors can “hide” from this searching method. Indeed, one can consider as initial approximations of unstable equilibria (if any), whose expressions depend on parameter, which will lead certainly to all self-excited attractors (if any). But, the chance to identify in this case hidden attractors diminishes significantly. Moreover, sudden changes in the bifurcations branches may indicate bifurcations when some branches may become unstable, but also they might be due to the mentioned dependence of the attraction basins with the parameter modification, when the bifurcation diagram is incomplete. For example, for this system of IO, for values of $a$ close to $a=0.05$, with relative small perturbations of $a$, of order of $1e-3$ or even of $1e-4$ (case $a=0.05$ and $a=0.0509$), the attraction basins change their properties related to hiddenness and self-excited.

A coexisting window in the space of fractional order indicates that analyzing a FO system by following only the bifurcation in the parameter space, without some analysis of the bifurcation diagram in the fractional order space, could be insufficient to identify the entire spectrum of system dynamics.

An important problem, for this system, is a deeper analysis of his FO variant.

The coexistence of attractors in the case of this system is influenced by the presence of the saddles on the separatrix of attraction basins. For all these cases, when the unstable manifolds of the saddle could separate the attraction basins (see the cases in Fig. 4 (c), Figs. 6 (c) and (e) and Fig. 7 (c)), this problem could be analyzed further as a possible sufficient criterion for non-existence of hidden attractors.

Because the variable $x_{3}$ represents the foreign capital inflow [Miloslav(2001)], the presence of hidden attractors could represent an important phenomenon from the economic point of view. In fact under some circumstances such as initial conditions, the dynamics of the system could cover some (hidden) behavior. Also, counterintuitively, due to the existence of hidden attractors, increasing (but also decreasing) the $x_{3}$ variable is not necessarily related with the system economical stability.