Study of Chaotic and Regular Modes of the Fractional Dynamic System of Selkov

. The paper investigates the dynamic modes of the Sel’kov fractional self-oscillating system in order to simulate the interaction of cracks. The spectra of the maximum Lyapunov exponents, constructed depending on the parameters of the dynamic system, are used as a research tool. The maximum Lyapunov exponents were constructed according to the Benettin-Wolf algorithm. It is shown that the existence of chaotic regimes is possible. In particular, the spectrum of the maximum Lyapunov exponents of the order of the fractional derivative con-tains positive values, which indicates the presence of a chaotic regime. Phase trajectories were also constructed to conﬁrm these results. It was also conﬁrmed that the orders of fractional derivatives are responsible for dissipation in the system under consideration.


Introduction
In [1], the authors proposed an interesting approach to describing the interaction of cracks in an elastic-brittle medium, which is based on the use of the Selkov's nonlinear dynamic system, which is investigated within the framework of the theory of biological systems [2]. The Selkov dynamic system [1], according to the authors, describes well the interaction of two types of cracks: the first type is seed cracks with a lower energy, which, when the critical level of their concentration is reached, transforms into the second type of cracks with a higher energy. Cracks of the second type are a source of microseismic phenomena (vibrations) and after the release of their energy, they partially pass into seed cracks. Further, this self-oscillatory process is repeated.
In the work of the authors [3], the generalized dynamic system of Selkov was considered, taking into account heredity (memory), which is studied within the framework of hereditary mechanics [4,5]. This model takes into account the property of a dynamic system to remember the impact on it for some time. Accounting for the heredity of the system is determined using derived fractional orders in the sense of Gerasimov-Caputo. In [3], using the numerical method of Adams-Buschfort-Multon, the calculated curves were constructed and the rest points of the fractional Sel'kov system were investigated.
In this paper, we consider the dynamic modes of the Sel'kov fractional system by constructing the maximum Lyapunov exponents.

Basic concepts and definitions
Let's introduce some concepts and definitions. Definition 1. A Sel'kov fractional dynamical system is a system of the form: where x(t) is a function that determines the concentration of seed cracks of the first type; y(t) is a function that determines the concentration of cracks of the second type that generate microseisms, t ∈ [0; T ] is a coordinate that is responsible for the current process time, T > 0 is a constant, modeling time; x 0 , y 0 , v, a, b -given positive constants; operators of fractional differentiation are understood in the sense of Gerasimov-Caputo of orders 0 < α i < 1 and are defined: properties that can be found in the monographs [6][7][8].
Consider a more general fractional non-autonomous dynamical system: where the functions f 1 (x(t), y(t), t) and f 2 (x(t), y(t), t) possess the necessary smoothness properties for the existence of the first partial derivatives. Definition 2. The system of equations in variations for the fractional dynamical system (2) is the following system of the form: Remark 1. For the Sel'kov fractional dynamical system (1), the variational equations, according to (3), will have the form: The system of equations in variations (3) for the original dynamical system (1) is key in constructing the maximum Lyapunov exponents. These equations describe the evolution of an infinitesimal perturbation of the phase trajectory of the original dynamical system.

Algorithm for constructing the maximum Lyapunov exponents
Consider an algorithm for constructing the maximum Lyapunov exponents for the Selkov fractional dynamical system (1), which consists of the following stages: 1. We choose the starting point -the initial vector x 0 and together with it we will track the K perturbed trajectories. In our case K = 2.
2. We solve numerically, for example, by the Adams-Bashfort-Multon method [3,9], the original system (1) together with three sets of perturbed equations or equations in variations (4). The total number of equations to be solved can be determined by the formula n(n + 1), where n is the number of equations in the original system. As the initial vectors for the variational equations, it is necessary to choose a set of vectorsx 0 0 ,ỹ 0 0 , which are orthogonal and normalized to unity.
3. After time T , the trajectory will go to the point -vector x 1 , the perturbation vectors x 1 ,ỹ 1 are renormalized using the Gram-Schmidt method according to the formulas: Here the notation ·, · is the scalar product of vectors.
4. Then we continue counting from the point x 1 and the vectors of perturbationsx 0 1 ,ỹ 0 1 . After the next time interval T , we obtain a new set of perturbation vectorsx 2 ,ỹ 2 , which is subjected to orthogonalization and renormalization. 5. Stages 2 -4 are repeated M times and in the course of calculations, the sums: ln ỹ i are calculated, in which the disturbance vectors appear before renormalization, but after orthogonalization. 6. The maximum Lyapunov exponents are calculated by the formula: Remark 2. The Gram-Schmidt orthogonalization procedure is necessary in order to exclude the dominance of the component of the maximum Lyapunov exponent at large times when calculating vectors along the phase trajectory. Otherwise, the task will be illconditioned.
Remark 3. It is important to study the spectrum of the maximum Lyapunov exponents, which is constructed depending on the values of the control parameter of interest to us. The following can be used as control parameters for system (1): v, a, b, α 1 , α 2 .
Therefore, we will investigate the following dependences of the maximum Lyapunov exponent: Λ max (v), Λ max (a), Λ max (b), Λ max (α 1 ), Λ max (α 2 ). If the spectra of the maximum indicators have positive values, then there is a chaotic regime, and negative ones -a regular regime.

Research results
Of greatest interest is the spectrum of the maximum Lyapunov exponents Λ max (α 1 ), since there are positive values of the spectrum (figure 1). Figure 1. Spectrum of maximum Lyapunov exponents Λ max (α 1 ).
In figure 1 shows that there is a range of variation of the alpha parameter, at which it takes positive values. The spectrum was plotted for the following parameter values: v = 0.5, a = 0.1, b = 1, t ∈ [0; 100], x 0 = 0.2, y 0 = 0.3 with a step τ = 0.05. Also, figure 1 shows the phase trajectories for different values of α 1 . At α 1 = 0.5, the phase trajectory has the shape of a twisting spiral (stable focus) (figure 1a), for α 1 = 0.9 the phase trajectory is a limit cycle (figure 1b), and for α 1 = 0.1 the open phase trajectory has a complex (chaotic shape) ( figure 1c).
Also, the phase trajectories in figure 1 show that the order of the fractional derivative α 1 is responsible for the energy dissipation in the system, i.e. as the value of alpha decreases, the limit cycle goes into a stable focus. Similarly, it can be shown that for order α 2 these conclusions remain valid.
The spectra of the remaining parameters are shown in figure 2. We see that all the spectra of the maximum Lyapunov exponents take negative values, therefore, we have no chaotic regimes. The spectra were plotted according to the following parameter values (t ∈ [0; 100], x 0 = 0.2, y 0 = 0.3 for all cases): • Λ max (a): v = 0.5, α 1 = α 2 = 1, a ∈ (0; 10], with a step τ = 0.1;

Conclusion
In this work, we investigated the spectra of the maximum Lyapunov exponents of the parameters of the Sel'kov fractional system (1). The Λ max (α 1 ) spectrum has positive values and, therefore, a chaotic regime is observed. The rest of the spectra of Lyapunov exponents have negative values, which indicates regular regimes. A more detailed picture of dynamic modes can be provided by building a map of dynamic modes, which will require large computational resources.
We can say that the orders of fractional derivatives introduce dissipation into the process under consideration (figure 1). In figure 1 it can be seen that with a decrease in the value of the parameter, the limit cycle is destroyed, and the phase trajectory becomes open. This confirms the results of studies of various fractional oscillatory systems [10][11][12][13].
Further continuation of the work can be based on a generalization of the orders of fractional derivatives in system (1) as functions of time: α 1 (t), α 2 (t).