Chaos Prediction in Fractional Delayed Energy-Based Models of Capital Accumulation

. This paper presents the nonlinear dynamic analysis of energy-based models arisen from the applied systems characterized by the energy transport in the presence of fractional order derivative and time delay. The studied model is the fractional version of Bianca-Ferrara-Dalgaard-Strulik (BFDS) model of economy which is viewed as a transport network for energy in which the law of motion of capital occurs. By considering the time delay as bifurcation parameter, a proof to investigate the existence of Hopf bifurcation and the phase lock solutions using the Poincare-Linstedt and the harmonic balance methods is given. At deﬁnite values of time delay, period-doubling bifurcations followed up by the consequences of chaotic states are detected. Simulation results assure that the BFDS model can generate new (hyper) chaotic attractors beyond half order derivatives through the eﬀect of the time delay on that system. In accordance with the literatures related to the problem of chaos, the concluded results conﬁrm the proposed theorem by El-Borhamy in which the time delay possesses the ability to change the dynamic state of nonlinear systems from regular to chaotic within the fractional order derivative domain.


Introduction
Nowadays, the mathematical modeling using the fractional calculus, leading to the fractional differential equations (FDEs) or systems, becomes a new and hot challenge topic due to their frequent appearance in various physical and engineering applications, cf.[1,2].Some of the semianalytical methods are developed for solving fractional problems such as the adomian decomposition method, the homotopy perturbation method, the homotopy analysis method, and the perturbation techniques, cf.[3][4][5][6].Another powerful method which can also give an explicit form for the solution is the variational iteration method (VIM), cf.[7,8].Consequently, considerable attention has been given to the qualitative and numerical analysis of systems with fractional order of the physical interest, because it is difficult to find exact solutions, so that a necessary object for several investigations of any scientific study in physical engineering and economical disciplines in this stream are focused on the stability, bifurcation, and chaos analysis, cf.[9][10][11][12].e fractional order systems are characterized by two important mathematical features.First is the existence of an unlimited memory due to the fractional order while the integer order systems have a limited memory.
e importance of this characterization roughly generalizes the mathematical treatment of the studied systems.Secondly, the existence of regular oscillations in fractional order systems requires very restricted conditions on the used fractional operator or the sustained periodic motion is not easy to be detected as well, cf.[13,14].In contrast, integer order systems possess the existence of regular oscillations which are periodic solutions or isolated periodic solutions (limit cycles).
Retarded functional differential equations form a class of equations characterized by the existence of the time delay such that their solutions are influenced by the history.Such retarded equations generate infinite-dimensional dynamical systems, and they can be classified into various types of mathematical problems with delay (there are point delays, distributed delays, state-dependent delays, and integrals within or taken over the delay), cf.[15,16].
It is generally accepted that the presence of delays in systems is a potent source of nonstationary phenomena such as periodic oscillations, instabilities, and chaos, cf.[17,18].
is can manifest itself as the loss of stability of an otherwise stable-steady state if the delay exceeds a certain threshold related to the dominant time-scale of a system.However, also there exists an evidence for that the time delay enhanced the stability, and short delays were shown to stabilize otherwise unstable dynamical systems.Also, it has been suggested that delays can damp out oscillations; this was shown with models for coupled oscillators under the condition that the delays in mutual interactions exceed a threshold value, cf.[19].Indeed, in fractional systems with a history, the existence of time delay can support the system to create a sustained periodic output with regular modes.
e study of stability and routes to chaos has become crucial for delayed fractional systems which has been taken care by many researchers, cf.[20][21][22].Particularly, the study of bifurcation and chaos is an interesting and frustrating phenomenon and it is studied and developed intensively since founding the first chaotic attractor by Lorenz in [23].In general, any (hyper) chaotic system is a nonlinear deterministic system that displays a complex and unpredictable behavior.Several methods are used to measure the existence of chaotic attractors, for instance, calculation of Lyapunov exponents, bifurcation diagrams, phase plane diagrams, and Poincare's maps.Somehow, Lyapunov exponents are considered an excellent tool to decide whether the attractor is (hyper) chaotic or not, cf.[24,25].Meaningly, if the system is chaotic, then at least one of the Lyapunov exponents is positive and in the case of a hyperchaotic system, two or more Lyapunov exponents are positive.e contribution of this work is centered on to generate a new chaotic attractor from delayed fractional model of capital accumulation for the considered BFDS model in [26].
e study will be based on the influence of the fractional order and the time delay parameters on the chaotic dynamic state.Basic properties of the new chaotic model will be analyzed by means of bifurcation diagrams, phase plane diagrams, and Lyapunov exponents.
e present work is organized as follows: Section 2 presents the fractional version of BFDS model.Section 3 handles the linear stability analysis of the fractional model and the characteristic function.Section 4 introduces Hopf bifurcation analysis for the system.Section 5 presents existence of the limit cycle using the perturbation analysis of the system.Section 6 presents the influence of time delay on the chaotic state of the system and the fractional order derivative as well.In the last, the conclusion is given.

Fractional Delayed BFDS Model
By building on the contribution by Bianca-Ferrara and Dalgaard-Strulik in [26][27][28] for the economy model, this work will develop the mathematical modeling for an economy growth which is viewed as a transport network for energy by including the fractional order derivative.In general, the energy conservation principles of dynamical systems have been proposed for the modeling of the law of motion of certain grown phenomena in physical sciences, cf.[28][29][30][31][32][33], so that the law of motion of capital occurs with a fractional time rate of energy and exhibiting whenever the system behavior is dependent at least in part on its history.e principle of energy conservation implies where μ is the energy requirement to operate and maintain the generic capital good, ] is the energy costs to create a new capital good, k(t) denotes capital stock assuming that the population grows at a constant rate of n, and D q t k(t) represents Riemann-Liouville fractional derivative of k(t).Specifically, the generalization of the BFDS model by including the energy conservation equation with a fractional derivative and time delay takes care of the previous occurring dynamics.
Due to Dalgaard and Strulik mathematical model in [28], the energy is derived as where 0 < a < 1 is a real constant, proportional to the dimension and efficiency of the network, where D is the dimension of the network, η ∈ [1, D] depends on the efficiency of the network, and ε > 0 is a real constant in the sense that it is independent of capital per worker.In specific words, the base ground of the scaling coefficient, a, should fall in [(1/2), (D/(D + 1))] and for a threedimensional network, the scaling exponent (a) should fall in the interval [(1/2), (3/4)], depending on the efficiency of the system.But, in this work, the domain of the scaling coefficient (a) has been extended within the interval (0, 1).Consequently, the law of motion for capital is described by the following general nonlinear fractional delay differential equation: for some initial functions k(t) � ϕ(t)∀t ∈ [− τ, 0], where the ratio (μ/]) captures the physical phenomenon of capital depreciation.e common definitions concerning the fractional derivative are usually being Gr € unwald-Letnikov fractional derivative, Riemann-Liouville fractional derivative, and Caputo's fractional derivative which are mostly used to express the fractionalized dynamical systems, cf.[1,2].For a function k(t) defined on an interval [t o , T], the Riemann− Liouville fractional integral of k(t) of order q > 0 is defined by and Riemann-Liouville fractional derivative of k(t) of order q > 0 is defined by where ℓ − 1 < q < ℓ(i.e., ℓ is the integer order) and Γ(q) is the gamma function.e definition of Caputo's fractional derivative of order q > 0 is e relation between Caputo and Riemann-Liouville fractional derivatives is illustrated by the following expression: Moreover, Caputo's fractional derivative of k(t) � c (constant) is always 0, whereas in the case of zero value of However, by setting t o ⟶ − ∞ in both definitions and requiring reasonable behavior of k(t) and its derivatives, we end up with the same formula for both of them as the following: e fractional derivatives of sine or cosine functions in the Riemann-Liouville sense have no periodic behavior and they are where E q 1 ,q 2 (t) is the two-parameter Mittag-Leffler function, When the time tends to infinity, it is corresponding to Weyl's definition of fractional derivative, and then equation (11) converges to the following periodic functions: Also, by setting t o ⟶ − ∞ in equation ( 6), it is corresponding to Liouville's definition of fractional derivative, and it also gives the following periodic results: − ∞ D q t exp(ωt) � ω q exp(ωt). (14c)

Linear Stability Analysis
3.1.e Linearized Form.Equilibria or the stationary points, of course, coincide with the corresponding points for the zero fractional rate and zero delay, i.e., D q t k(t) � 0, at τ � 0(corresponding to the classical case when q ⟶ 1, D q t k(t) ⟶ k ′ (t) � 0).Hence, there exists a unique positive steady state k * , satisfying the relation εk * a− 1 � μ + n].After setting the following translation k(t) � x(t) + k * , one can obtain Hence, Taylor's expansion at the zero equilibrium reads and its linearized form is where k p � (a − 1)(μ/]) < 0. Equation ( 17) can be reduced to where  k p � k p + n.

Particular Solutions of the Linearized Form.
e prototype linear rate equation without delay and fractional derivative for a rate constant k p will have the following decaying solution: where x o is an arbitrary constant.
In fact, on the other hand, the same equation with a constant delay and integer order derivative will have multiple solutions ranging from monotonic decay to damped oscillations, stable periodic oscillations, and undamped Complexity growing oscillations based on the range of time lag values.For instance, a particular solution for the abovementioned linear rate equation with constant delay and integer order derivative at the zero rate of population can be obtained as Some methods are developed to obtain the solution for the delay differential equations (DDEs) for a given set of initial values, where the function x(t) is defined over the interval (− τ, 0) such as the method of steps and the extended one-step techniques, cf.[34,35].
To follow-up, the same equation, without a constant delay but with a fractional derivative in Caputo's sense, will have the following decaying solution: e general solution of a linear fractional delay problem (n � 0), seek in Liouville's sense, is

e Characteristic Function and Eigenvalues.
Using the fundamental tools of fractional calculus, the characteristic function of the linear fractional time delay case is It is easy to see that, en, one can find that there is at least one λ < 0, satisfying Λ(λ) � 0. Using such value of λ, the solution  x(t) � x o e − λt becomes a positive solution with so that as long as x(t) remains positive.But x(t) could become zero without first reaching  x(t), hence under the condition By considering equation (23), en, using Lagrange's expansion, let z be that root of which has the value F(z) � F(z o ) when ξ � 0. If F(z) and ϕ(z) are analytic functions inside and on a circular contour (C) containing z � z o , then the roots satisfy Hence, if we consider that then so that one obtains en, en, the eigenvalues can be represented by Regarding the zero rate of population case (n � 0) in equation (29), one can get the following eigenvalues using the following expansion: If 0 < τ ≪ 1, we have the following linearized part of equation (37): 4 Complexity en, the root of the characteristic function is approximated by For the case of τ q k p � − 1, the eigenvalues of equation (39) become In equation ( 39), in the case of integer derivative (q � 1), it can be reduced to If we consider in equation ( 29) that q � 1, then the eigenvalue is given by if the following condition holds: By considering that λ � α + iω and to obtain the transition values (α � 0), we assume that λ � iω c , and one obtains which is a root of equation (23) and ω c > 0 is a real number.en, one gets Separate the real and imaginary parts to obtain From squaring and adding the above two equations, one obtains us, one can get the following constraint: which means that the characteristic equation has no purely imaginary roots.Consequently, we can state the following theorem.
, where m ∈ Z, then the zero solution of the linearized fractional time delay equation is globally asymptotically stable at zero rate of population.

Local Stability Analysis of the Algebraic Retarded
Equation.
e algebraic retarded equation is obtained when q ⟶ 0 as a nondifferentiable version for equation (4) as follows: By considering that the equilibrium solution is k e , it reads so that the equilibrium solution, represented by equation (50), is locally asymptotically stable, if all the roots λ satisfy the following condition: where

Hopf Bifurcation Analysis
Now, we are going to derive an expression for the critical time delay that separates stable and unstable regions bifurcated from k * , the unique positive equilibrium for the mathematical model.Moreover, we are going to prove a theorem which characterizes the nature of the equilibrium point k * .Let λ j � α j (τ) + iω j (τ) in general, so at α j (τ c ) � 0 and ω j (τ c ) � ω c by substituting in the characteristic equation, then en, one obtains us, at n � 0, one easily obtains Complexity erefore, to characterize the existence of the transition curves, we state the following theorem.Theorem 2. Let k * be the unique positive equilibrium point for the mathematical model.
en, there exists a positive number τ c such that the equilibrium k * is asymptotically stable for τ ∈ [0, τ c ) and unstable for τ > τ c .
Proof.Let λ j � α j (τ) + iω j (τ) denote a root of equation ( 23) near τ � τ c , satisfying α j (τ c ) � 0 and ω j (τ) � ω c .Differentiating the characteristic equation (equation ( 23)) with respect to τ, obtaining the change of sign of the real part and at q � 1, this gives us, from equation ( 56), if we have that (dα j (τ))/dτ| τ�τ c > 0, this implies that all the roots crossing the imaginary axis at iω c from left to right as τ increase and this result yields in the loss of stability.en, the conclusion holds.

Construction of Periodic Solution
In this section, we investigate the existence of periodic solution without periodic input; otherwise, it is generated by the existence of time delay.Necessary conditions for generating periodic solutions, and the stability of periodic solution of equation ( 4) at τ c are given by using the method based on Poincare-Lindstedt method and harmonic balance method, cf.[5,36,37].
Following the procedures in [17], we start by considering Taylor's expansion of equation ( 4) up to the third order at the zero equilibrium using Liouville's sense of fractional derivative as obtained by equation (16).en, we begin by introducing a small parameter ε via the scaling en, for values of delay τ which lie close to τ c , expansion series is introduced as follows: Next, we stretch time by replacing the independent variable t by s, as follows: s � ωt, and d q dt q � ω q d q ds q . (61) After substitution, we get    6 Complexity x(s) , and x(s − ωτ) in a power series in ε, By recognizing that ω 1 � 0 and τ 1 � 0, then we obtain Finally, substituting and collecting terms, we obtain Order of ε 1 : We take the solution of the x 0 (s) equation as follows: x 0 (s) � A 0 cos(s). (66) en, we obtain an expression for x 1 (s): Next, step substitute the expressions of x 0 and x 1 into x 2 equation after long simplifications; the coefficients of the resonant terms, sin s and cos s, are equated to zero. is results in A 2 o as follows: By including the zero rate of population n � 0, one can obtain where Complexity where us, A 0 is real so A 2 0 > 0, this means that τ 2 should have the same sign as (M/N).If q � 1, this leads to the classical model with In order for A 2 0 to be positive, τ 2 is required to be more than zero, which means that the limit cycle is born out of unstable fixed point.Since the stability of the limit cycle should be the opposite of the stability of fixed point from which it is born, our conclusion is that the limit cycle is stable and we have supercritical Hopf bifurcation, so that the following theorem is concluded.Theorem 3.For q ∈ (0, 1], a stable limit cycle is born out of unstable fixed point with radius A 0 .
Numerically, these results are verified by Figures 3-6 to confirm this conclusion at the values a � 0.001, ε � 1, μ � 1, and ] � 2. e economy oscillator settles down to a steady oscillation with uniform periods at each fractional order derivative. is is clearly appeared in the phase plane for t ∈ [900, 1000] by avoiding the transients which have completely disappeared at the scale shown.

Routes to Chaos
Dynamic trajectories can be drawn to check the existence of aperiodic or chaotic motion in BFDS model system.More or less most of the analytical methods cannot provide enough information to specify the perfect existence of chaotic motion.erefore, methods such as bifurcation diagrams, Lyapunov exponent diagrams, and Poincare's maps or phase trajectory diagrams can identify the behavior of the system motion through given parameters. is study will consider the most sensitive parameter produced chaos, which is the time delay, cf [17].In most cases, studying the effect of the time delay to create chaos in delayed systems is indispensable and more important due to its inherence with them.
In our analysis, started by using the bifurcation diagrams, we obtain valuable insights into system's nonlinear analysis, since with a range of such parameter variation.Some what, the system dynamic is mostly to be observed under the influence of that parameter.Bifurcation diagrams are drawn under the variation of time delay with constant interval.e bifurcation can be easily detected by graph the relationship between the solution (k(t)) and the time delay under the step size Δτ � 0.001 for the shown figures in Figure 7.It can be vividly observed that its existence affects the characteristic of the dynamic trajectories.Roughly speaking, at different values of fractional orders q, the system has period-doubling behavior according to the fixed parameters, n � 1, ε � 15, μ � 0.95, and ] � 0.29, for a � 0.125.e figures started with q � 0.5, q � 0.75, and q � 1, respectively.It is obvious from each figure that, with further increasing of the time delay, the model response gradually started by a stable response up to a certain threshold and after that it enters a chaotic zone by the route of perioddoubling bifurcation supported by the value of fractional derivative. is indicates that it might be predicted a strange attractor, representing the chaotic and irregular motion in the model as shown in the figures: Figures 8-10 at a � 0.125, q � 0.25, and τ �1, 2.5, and 5, respectively.Figures 11-13 at a � 0.125, q � 0.5, and τ �1, 2.5, and 5, respectively.Figures 14-16 at a � 0.125, q � 0.75, and τ �1, 1.5, and 2, respectively.Figures 17-19 at a � 0.125, q � 1, and τ �1, 1.3, and 1.45, respectively.
Indeed, these results assure that the prediction of the chaos is not far away to be obtained within some specific values of the model parameters taking into the account that the branch points have to be changed with further increasing of the parameter a.Moreover, this model with q ∈ (0, 1] can generate hyperchaotic behavior due to the two positive values of Lyapunov exponents as shown in Figure 20, so that the obtained figures indicate the generation of new strange attractor for fractional order equation with delay based on the three measurements (bifurcation diagram, phase trajectories, and Lyapunov exponents).ese results confirm the results in [17] by the responsibility of the time delay to create chaotic motion in delayed system with the support of the fractional order derivative.12 Complexity

Conclusion
is work has studied the dynamic change of a single nonlinear equation represented by the BFDS model with delay supported by the fractional order derivative within the interval (0, 1]. e dynamic of the system fluctuates due to the time delay parameter to generate stable periodic, aperiodic, or chaotic behavior which has been evaluated based on three measurements: bifurcation diagrams, phase diagrams, and Lyapunov exponents.Simulation results showed that the BFDS model can generate period-doubling bifurcation, leading to various new chaotic attractors based on the time delay parameter and the fractional order derivative.Comparatively, the obtained figures indicate that the generation of new hyperchaotic strange attractor is mostly predictable for every value of the fractional derivative of the model.Consequently, based on the results inferred from the BFDS model, it has been confirmed that the time delay is able to change the dynamic state of nonlinear systems from regular to chaotic which is in agreement with El-Borhamy's theorem in [17] within the fractional order derivative domain q ∈ (0, 1]. Figure 20: Lyapunov exponents at a � 0.125 and τ � 15 for q � 0.5, q � 0.75, and q � 1, respectively. Complexity

Figures 1 and 2
Figures 1 and 2 illustrate the validity of eorem 1 at certain value k p and q.

Figure 3 :Figure 6 :Figure 4 :Figure 5 :
Figure3: e birth of Hopf-bifurcated limit cycle from the fixed point k * for the studied case when q � 0.25.