Results of numerical simulation.
The chapter proposes a mathematical model for a wide class of hereditary oscillators, which is a Cauchy problem in the local formulation. As an initial model equation, an integrodifferential equation of Voltaire type was introduced, which was reduced by means of a special choice of difference kernels to a differential equation with nonlocal derivatives of fractional-order variables. An explicit finite-difference scheme is proposed, and questions of its stability and convergence are investigated. A computer study of the proposed numerical algorithm on various test examples of the hereditary oscillators Airy, Duffing, and others was carried out. Oscillograms and phase trajectories are plotted and constructed.
- mathematical model
- cauchy problem
- derivative of fractional order
- finite-difference scheme
- phase trajectory
In the paper of the Italian mathematician Vito Volterra , the notion of heredity (memory), a property of a dynamical system characterized by nonlocality in time, is introduced, which consists in the dependence of its current state on a finite number of previous states. In another paper , Volterra investigated the hereditary oscillator—a vibration system with memory, which was written in the form of an integrodifferential equation with a difference kernel, a function of memory. Further, for such an oscillator, Volterra derived the law of total energy, in which an additional term appeared, responsible for the dissipation of energy in the vibrational system. This fact was confirmed in subsequent works.
In papers [3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21], fractal oscillators were considered, which represent the class of hereditary oscillators with a power-law function of memory. The peculiarity of such oscillators is that their mathematical description can be reduced to differential equations with nonlocal derivatives of fractional constant orders, which are investigated within the framework of the theory of fractional calculus [22, 23, 24].
In papers [6, 7, 8, 9, 11, 12, 13, 14, 18, 19, 20, 21], models of fractal linear oscillators were investigated in the sense of the Gerasimov-Caputo derivative and in papers [4, 10, 16]—in the sense of the Riemann-Liouville derivative. Analytical solutions of model equations in terms of a special function of Mittag-Leffler-type and generalized Wright-type function, oscillograms, and phase trajectories are constructed. It is shown that in the regime of free oscillations, the presence of memory effects in the system leads to attenuation of oscillations as a result of energy dissipation, and with allowance for external periodic action, it is possible to stabilize the amplitude of the oscillations, with the phase trajectories reaching the limit cycle and also the resonance effect.
In papers of the author [25, 26], the fractal parametric resonance (the fractal Mathieu oscillator) was investigated, and the Strutt-Ince diagrams of parametric resonance existence areas were constructed. It is shown that these regions strongly depend on the orders of the value fractional derivatives entering into the initial equation.
In a monograph by the Slovak mathematician Ivo Petras , the fractal nonlinear oscillator models whose differential equations contained fractional derivatives in the sense of Riemann-Liouville were considered and analyzed using numerical methods and considered the stability of the rest point of oscillatory systems. However, the stability and convergence of numerical methods have not been considered.
A further continuation of the investigation of hereditary oscillators is associated with the introduction of the derivatives of fractional variable orders in the model equations. This is due to the fact that the orders of fractional derivatives are related to the properties of the medium in which this or that process takes place and changes with time under the influence of external influence. Therefore, papers [27, 28, 29, 30] proposed that the models of fractal nonlinear oscillators (Duffing, Van der Pol, Van der Pol-Duffing, FitzHugh-Nagumo) were proposed and investigated using explicit finite-difference schemes, whose equations contain both the derivatives of the constants and variable fractional orders of the Gerasimov-Caputo and Riemann-Liouville types. With the help of computer experiments, the convergence of finite-difference schemes was shown, and estimates of the computational accuracy of the method were obtained; oscillograms and phase trajectories were constructed. However, the questions of stability and convergence were not formulated in the form of corresponding theorems.
In [31, 32], a new class of fractal oscillators was proposed and investigated; active fractal oscillators (AFOs)—nonlinear oscillators with external influences, which include the fractional Riemann-Liouville integral—were investigated. Such oscillators are constructed on the basis of the scheme of the radioelectronic аutogenerator with a fractional feedback circuit. Authors use the method of equivalent linearization to investigate AFOs and come to the conclusion that the self-oscillator is isochronous.
From the analysis of the above publications on the study of hereditary oscillator, we can conclude that the main tool for their study is numerical methods, for example, finite-difference schemes. In most cases, the authors leave without considering the questions of stability and convergence of finite-difference schemes and, even if they touch, then without formulating the corresponding theorems and proofs. Therefore, the goal of the present paper is to construct a finite-difference scheme for a wide class of hereditary (fractal) linear and nonlinear oscillators, prove its stability and convergence, formulate results in the form of corresponding theorems, and study finite-difference schemes on specific test examples.
2. Formulation of the problem
Consider the following model integrodifferential equation for the function , where :
where , —given constant.
Eq. (1) describes a wide class of hereditary, depending on the form of the right-hand side (function ) of linear or nonlinear oscillators.
Definition 1. Functions and —difference kernels in Eq. (1)—will be called memory functions, since they define the notion of heredity (memory), which was introduced in the work of the Italian mathematician Vito Volterra .
Definition 2. A nonlinear function on the right-hand side of Eq. (1) satisfies a Lipschitz condition with respect to a variable :
Eq. (1) describes a broad class of hereditary nonlinear oscillators, depending on the form of the function on its right-hand side, and the parameter has the meaning of the coefficient of friction.
Note that in Definition 1, the memory functions and can be chosen arbitrarily, depending on the conditions of the particular problem. We will choose these functions power law, since power laws are often found in various fields of knowledge . We choose the memory functions and in the form
where , —Euler gamma function.
We give the following definitions.
Definition 3. Derivatives of fractional variables of orders and Gerasimov-Caputo type: we call the following operators of fractional differentiation:
We note that in the case when the functions and in the relations (4) and (5) are constants, we arrive at the definitions of the fractional derivative in the sense of Gerasimov-Caputo [34, 35], and in the case when these constants and , the operators of fractional differentiation (4) become classical derivatives of the second and first orders.
Taking into account Definition 3, the model Eq. (1) can be rewritten in a more compact form:
For Eq. (5), the initial conditions in the local formulation are valid:
where and —given constants. As a result, we arrive at the Cauchy problems (5) and (6), which we will investigate.
3. Explicit finite-difference scheme
We construct an explicit finite-difference scheme for the Cauchy problems (5) and (6). We divide the time interval into equal parts with a step . We introduce the grid function , where .
then the formulas will refer to the formula (7)
Here, to shorten the record, .
Taking into account relation (7), the Cauchy problems (5) and (6) in the difference formulation will have the form
Here, to shorten the record, . We write the problem (8) explicitly:
We note that the weight coefficients and have properties, which we formulate in the form of the following lemmas.
Lemma 1. For any weights, coefficients and , as well as coefficients and , have the following properties:
Proof. The first property follows the definition of weight coefficients:
The second property is proven in the following way. We introduce two functions:
where . These functions are decreasing. Really derived from these functions
Therefore, the second property holds. The third property follows also the properties of the gamma function. The lemma is proven.
Let and —approximations of differential operators of Gerasimov-Caputo types and . Then, we have the following lemma.
Lemma 2. Approximations and operators of the Gerasimov-Caputo types and satisfy the following estimates:
where and —constants that are independent of the parameter .
Proof. Using the first property of Lemma 1 and Definition 3, we obtain
Similarly, we can obtain the second estimate from Eq. (10).
The lemma is proven.
Investigation. According to Lemma 2, it can be shown that the explicit finite-difference scheme (9) has an error . This fact will be used in computer experiments in determining the computational accuracy of the numerical method.
Lemma 3. The sums in the finite-difference scheme (9) have the following representations:
Proof. The representation (11) follows the properties of the sum. Indeed, by opening the sums in Eq. (9) and grouping the terms properly, we arrive at the representation (11).
Using Lemma 3, the finite-difference scheme (9) can be rewritten as
Or in matrix form
where the matrix , :
Theorem 1. An explicit finite-difference scheme (9) converges with the first order if the following condition is satisfied:
Proof. Let be the exact solution of system (8) and the error vector . Then, system (8), with allowance for Lemma 2, can be written as follows:
We note that for any , the inequality holds . Consider the norm for the matrix :, we obtain
According to Lemma 1, we note that the inequality holds . Suppose that the condition is satisfied , then the second diagonal element satisfies the inequality in the matrix , and the remaining diagonal elements are equal to the inequality ; these elements with , in view of properties 2 and 3, tend to be zero.
The remaining elements of the matrix (16) also possess these properties. We also note that the matrix is a matrix with a diagonal predominance for small values .
Therefore, the sum of the elements of the second row in the matrix satisfies the condition . Further, by virtue of properties 2 and 3 of Lemma 1, it is obvious that the sum of the remaining terms also satisfies these conditions. Therefore, the following estimate is valid for the norm: .
Note that for the values of the parameter the norm , however, the condition number is violated and the diagonal transformation is violated; therefore, it is necessary to decrease the step .
Further, for any constant independent of , and the error rate, the following estimate holds:
We introduce the notation in Eq. (17): . Then, we obtain the following estimate:
Substituting into Eq. (18) , we obtain
From the second initial condition (6) it follows: and .
Now according to our assumption , which leads us to the relation
Condition (19) begins to work at such values , when many of conditionalities arise , and for sufficiently small values , it suffices that the step satisfies the inequality . Therefore, we arrive at the relation (14). The theorem is proven.
We note that in  the authors used the classical Lax theorem, which holds for local finite-difference schemes, to prove the convergence of the scheme. For nonlocal finite-difference schemes, the convergence must be proven independently.
We consider the stability of an explicit finite-difference scheme (4). Suppose that and two are different solutions of the matrix Eq. (12) with initial conditions and .
Theorem 2. An explicit finite-difference scheme (9) is conditionally stable if condition (14) is satisfied and the estimate holds for any , where does not depend on the step .
Proof. We introduce the notation: Then, Eq. (12) can be written in the form . Here, as it was said in.
According to Theorem 1, we have the following estimate:
Therefore, the following estimate holds
With , we obtain and .
The last inequality follows the second condition of problem (6). Therefore, if there is a perturbation, then it does not lead to a large increase in the error of the numerical solution. However, for large values , many of conditionalities arise, and therefore it is necessary to decrease the step ; according to Eq. (19), for small values , the estimate is valid . Then, the system is stable if condition (14) is satisfied. The theorem is proven.
4. Results of modeling
Consider the work of the explicit finite-difference scheme (9) on specific examples. We show that the scheme (9) has the first order of accuracy. Since in the general case, the exact solution of the Cauchy problems (5) and (6) cannot be written in analytical form, we will use the double conversion method. For this, we introduce two parameters: —absolute error between the numerical solution in step and the numerical solution in step . Then, the order of computational accuracy can be estimated by the formula
We note that in the case when the fractional parameters in the scheme (9) do not change and have the following values of and , we arrive at the classical local explicit finite-difference scheme with the second order of accuracy.
The numerical algorithm (9) was implemented in Maple software.
Example 1. Suppose that the right-hand side in Eq. (1) has the form
We choose the initial condition (6) for simplicity by homogeneous:
We note that the Airy oscillator is used in optics in the simulation of Airy laser beams.
In this case, the explicit finite-difference scheme (9) has a more specific form:
For the explicit finite-difference scheme (20), we choose the following values of the control parameters: . , , , and . And during the simulation, we will change the number of nodes in the calculation grid.
Note that the values of the selected parameters for Example 1 satisfy the conditions of Theorems 1 and 2, which is indirectly confirmed by the results of modeling for different values of the nodes of the computational grid (Table 1).
From Table 1 we can notice that when the number of calculated nodes in the grid doubles in nodes , the maximum error in absolute value decreases twice, and the order of computational accuracy tends to unite.
This confirms that the explicit finite-difference scheme (9) and in particular the scheme (20) for Example 1 have the first order of accuracy, and since condition (14) is satisfied, then convergence with the same order.
In Figure 1 the oscillogram (Figure 1a) and the phase trajectory (Figure 1b) are shown for Example 1 at the parameter value It can be noted that with time the amplitude of the oscillations is established and as a result the phase trajectory reaches the limit cycle. Another situation arises in the case of free oscillations (Figure 2).
The amplitude of the oscillations decays (Figure 2a), and the phase trajectory twists into a spiral (Figure 2b). The dissipation of energy in this case occurs as a result of the presence of friction with a coefficient and also the “memory” effect, which gives an additional term in the ratio for the total energy of the oscillatory system (Figure 3).
This fact is confirmed by the results of . Consider the following example of a nonlinear hereditary oscillator.
Example 2. Let that in Eq. (1) the right-hand side has the form.
and we choose the initial conditions (6) to be homogeneous:
The explicit finite-difference scheme (9) for this case has the form
For the explicit finite-difference scheme (21), we take the values of the control parameters as follows: , , , and .
Note from Table 2 that for Example 2, with an increase in the number of design nodes , the maximum error in absolute value decreases and the order of computational accuracy tends to unite. This indicates that the explicit finite-difference scheme (21) has the first order of accuracy.
Let’s perform numerical simulation according to the scheme (21) with the values of the following parameters, , , and , and leave the remaining parameters unchanged. Let us construct an oscillogram and a phase trajectory (Figure 3).
The oscillogram (Figure 4a) has a constant amplitude of a more complex shape at its minima and maxima, which is reflected in the phase trajectory (Figure 4b). The phase trajectory enters a complex two-loop limit cycle. The presence of such loops, apparently, is associated with the effects of memory in the oscillatory system .
Figure 5 shows the case of free oscillations for Example 2. It is seen that the presence of friction and memory effects in the oscillatory system intensify energy dissipation, which leads to damping of the oscillations (Figure 5a) and a phase trajectory—a twisting spiral (Figure 5b). Indeed, if there is no friction in the oscillatory system, we obtain an oscillogram and a phase trajectory as in Figure 6.
Example 3. Suppose that in Eq. (1) the right-hand side has the form
where is the spring travel speed; is the surface adhesion energy; is the frequency of free oscillations; and is the coefficients of the expansion of the Fourier series.
Eq. (5) with the right-hand side of Eq. (22) describes the hereditary stick–slip effect . The stick–slip effect is encountered in tribology problems, for example, when the movement of a load on a spring along a surface is investigated. Due to adhesion, the load adheres to the surface, and due to the tension of the spring, it breaks and slides along it, and its oscillations occur [41, 42]. The stick–slip effect can also be incorporated into the mechanical model of an earthquake in the subduction zone of lithospheres plates .
In  it was said that in order to obtain a reliable solution it suffices to take the first seven coefficients in the expansion of the function (22). The values of these coefficients are taken from  , , , , , , and . Values of control parameters are , , , , , , and .
Figure 7 shows the calculated displacement curves, displacement velocities, and phase trajectory. Figure 7a shows the oscillogram for Example 3. It can be seen that during the separation, the load experiences oscillations and the rate of such oscillations in the potential well attenuates rather slowly (Figure 7b). This effect is the eradication of the process. The phase trajectory in Figure 7c shows that the potential wells are stable focuses.
A mathematical model characterizing a wide class of hereditary oscillators is proposed and studied. The model is a differential Cauchy problem with derivatives of fractional-order variables of the Gerasimov-Caputo types (5) and (6). Using the theory of finite-difference schemes, a nonlocal explicit finite-difference scheme (9) was constructed with the first order of accuracy. Questions of its stability and convergence, which are formulated in the form of corresponding theorems, were studied.
The main result of the paper can be formulated as follows: an explicit finite-difference scheme is conditionally stable and converges if criterion (14) is satisfied. With the help of computational examples, it was shown that the scheme (9) has the first order of accuracy. It is confirmed that in the case of free oscillations, the presence of friction and heredity increases dissipation of energy, which leads to attenuation of oscillations.
One of the continuations of the investigation of the Cauchy problems (5) and (6) is a generalization of it:
Another continuation of the research is related to the introduction of other memory functions into the model Eq. (1), which leads to different model equations with different derivatives of fractional orders, and also the Cauchy problems (5) and (6) can be written in terms of the local fractional derivative [44–46].
The question of the stability of the rest points of dynamical systems described by the Cauchy problems (5) and (6) is also interesting, by analogy with the papers [47, 48].
The work was carried out according to the state task within the framework of research work Vitus Bering Kamchatka State University on the topic “Application of fractional calculus in the theory of oscillatory processes” No.AAAA-A17-117031050058-9 and with the support of the grant of the President of the Russian Federation MK-1152.2018.1.