*Article* **About Some Possible Implementations of the Fractional Calculus**

**María Pilar Velasco 1,\*,†, David Usero 2,†, Salvador Jiménez 1,†, Luis Vázquez 2,†, José Luis Vázquez-Poletti 2,† and Mina Mortazavi 3,†**


Received: 6 May 2020; Accepted: 26 May 2020; Published: 2 June 2020

**Abstract:** We present a partial panoramic view of possible contexts and applications of the fractional calculus. In this context, we show some different applications of fractional calculus to different models in ordinary differential equation (ODE) and partial differential equation (PDE) formulations ranging from the basic equations of mechanics to diffusion and Dirac equations.

**Keywords:** fractional calculus; fractional differential equations; nonlocal effects

#### **1. From Elementary Mathematical Analysis to Fractional Derivatives**

A first reason to justify the use of fractional operators is the need to introduce memory terms into differential models in a natural form. In this sense, we can consider the classical Calculus Fundamental Theorem with the introduction of a convolution kernel *F* associated to a function *g* leads to the natural form of introducing a memory term by changing the convolution kernel of the integral:

$$F(\mathbf{x}) = F(a) + \int\_{a}^{\mathbf{x}} K(\mathbf{g}(\mathbf{x}) - \mathbf{g}(\mathbf{s})) \cdot f(\mathbf{s}) d\mathbf{s} \tag{1}$$

Now *F* is the generalized primitive of the function *f* and the convolution kernel *K* is a memory term that could be different in each specific problem, changing the definition of fractional operator consequently. This generalization of the integral can be considered as a base to construct possible definitions for fractional integrals.

From these considerations, the fractional calculus emerges in the mathematical world as the study of integral and derivative operators of non-integer orders on domains of real or complex functions. Several definitions of fractional derivatives *Dα* have been developed progressively with the objective to generalize the concept of ordinary derivative *D*, such that for *α* = 1 the ordinary operator can be recovered [1–3].

Fractional calculus is a powerful mathematical tool that allows to create intermediate-order parameters equations and offers modeling scenarios where fundamental mathematical questions converge and appropriate numerical algorithms can de developed. A lot of fractional operators have been defined in the literature; however, not all of them can be used in each real-world application. In this context, we appreciate very much the enthusiastic and clarifier paper of D. Baleanu and A. Fernandez [4] and M.D. Ortigueira and J.A.T. Machado [5–7].

#### *1.1. From Factorial to the Gamma Function*

The first definitions of fractional operators are related to the use of the Gamma function is a function that generalizes the definition of the factorial to non-positive numbers. Its definition is:

$$
\Gamma(z) = \int\_0^\infty s^{z-1} e^{-s} ds\tag{2}
$$

for any complex number *z* with positive real part.

Using integration by parts in Equation (2), a fundamental property of the Gamma function is obtained:

$$
\Gamma(z) = (z - 1)\Gamma(z - 1),
\tag{3}
$$

which allows to give the Gamma function of a positive integer number as

$$
\Gamma(n) = (n-1)! \,. \tag{4}
$$

In this context, the Gamma function is a generalization of the concept of factorial.

In 1738, Euler introduces the first generalization of ordinary derivative, verifying that the fractional derivation made sense for the potential function *<sup>x</sup>a*. And in 1819, Lacroix starts from the *m*-order derivative of the function *<sup>x</sup>n*, with *m* and *n* positive integer numbers

$$\frac{d^m}{dx^m} \mathbf{x}^n = \frac{n!}{(n-m)!} \mathbf{x}^{n-m},\tag{5}$$

to determine the 1/2 order derivative of the function *<sup>x</sup>a*, using the generalization of the factorial function by the Gamma function:

$$\frac{d^{1/2}}{dx^{1/2}}\mathfrak{x}^a = \frac{\Gamma(a+1)}{\Gamma\left(a+\frac{1}{2}\right)}\mathfrak{x}^{a-\frac{1}{2}}\,\_{\prime} \tag{6}$$

such that for *a* = 1:

$$\frac{d^{1/2}}{dx^{1/2}}\mathbf{x} = \frac{\sqrt{\overline{\mathbf{x}}}}{\Gamma\left(\frac{3}{2}\right)} = \frac{2\sqrt{\overline{\mathbf{x}}}}{\sqrt{\pi}}.\tag{7}$$

 This is the result that will be obtained with the called Riemann–Lioville fractional derivative.

This concept can be generalized to any order and the following relation between the ordinary and fractional case with the Riemann–Liouville fractional derivative is obtained:

$$\frac{d^m}{dx^n} \mathbf{x}^m = \frac{m!}{(m-n)!} \mathbf{x}^{m-n} \quad \Rightarrow \quad \frac{d^n}{dx^n} \mathbf{x}^\mu = \frac{\Gamma(\mu+1)}{\Gamma(\mu-n+1)} \mathbf{x}^{\mu-n}.\tag{8}$$

Later, some provisional definitions of fractional operators were introduced by Fourier, Abel, Liouville and Riemman, without much success. Until, in 1870, N. Ya. Sonine started from the Cauchy formula for repeated integration:

$$(\,\_a\mathbf{I}\_x^n f\mathbf{)}(\mathbf{x}) = \int\_a^\mathbf{x} d\mathbf{x}\_1 \int\_a^{\mathbf{x}\_1} d\mathbf{x}\_2 \cdot \cdots \int\_a^{\mathbf{x}\_{n-1}} f(\mathbf{t}) dt = \frac{1}{(n-1)!} \int\_a^\mathbf{x} (\mathbf{x} - \mathbf{t})^{n-1} f(\mathbf{t}) dt\tag{9}$$

and, using the generalization of the factorial function by the Gamma function, he obtained the actual definition of the fractional integral of Riemann–Liouville:

$$(\,^{a}I\_{x}^{a}f)(\mathbf{x}) = \frac{1}{\Gamma(a)} \int\_{a}^{\mathbf{x}} (\mathbf{x} - t)^{a-1} f(t) dt, \qquad \Re(a) > 0 \tag{10}$$

although in 1884 Laurent formulated it definitively. *1.2. Some Definitions of Fractional Integrals and Derivatives*

An important definition of fractional integral and derivative corresponds to Riemann–Liouville:

• Left-side Riemann–Liouville Fractional Integral of order *α* > 0:

$$\,\_aI\_a^a\phi(\mathbf{x}) = \frac{1}{\Gamma(a)}\int\_a^\mathbf{x} (\mathbf{x} - t)^{a-1} \phi(t) dt, \qquad \mathbf{x} > a. \tag{11}$$

• Right-side Riemann–Liouville Fractional Integral of order *α* < 0:

$$\, \_\times I\_b^a \phi(\mathbf{x}) = \frac{1}{\Gamma(a)} \int\_\mathbf{x}^b (\mathbf{x} - t)^{a-1} \phi(t) dt, \qquad \mathbf{x} < b. \tag{12}$$

• Left-side Riemann–Liouville Fractional Derivative of order *α* > 0:

$$\,\_{a}D\_{x}^{a}\phi(\mathbf{x}) = \frac{1}{\Gamma(n-a)} \left(\frac{\partial}{\partial \mathbf{x}}\right)^{n} \int\_{a}^{x} (\mathbf{x}-\mathbf{t})^{n-(n-1)} \phi(\mathbf{t}) d\mathbf{t} = D^{n}(\,\_{a}I\_{x}^{n-a}f)(\mathbf{x}), \qquad \mathbf{x} > a. \tag{13}$$

• Right-side Riemann–Liouville Fractional Derivative of order *α* < 0:

$$\partial\_x D\_b^a \phi(\mathbf{x}) = \frac{1}{\Gamma(n-a)} \left(-\frac{\partial}{\partial \mathbf{x}}\right)^n \int\_{\mathbf{x}}^b (\mathbf{x}-\mathbf{t})^{n-(n-1)} \phi(\mathbf{t}) d\mathbf{t} = (-D)^n ({}\_x l\_b^{n-a} f)(\mathbf{x}), \qquad \mathbf{x} < b. \tag{14}$$

In all cases *n* ∈ N, such that 0 ≤ *n* − 1 < *α* < *n*.

These operators recover the classical operators for the parameter *α* = 1 and the algebra of these operators is different to the classical operators:

• Let *f* ∈ *Lp*(*<sup>a</sup>*, *b*) (1 ≤ *p* ≤ ∞) and *Re*(*α*), *Re*(*β*) > 0. Then:

$$(({}\_a I\_x^a {}\_a I\_x^{\beta} f)(\mathbf{x}) = ({}\_a I\_x^{a+\beta} f)(\mathbf{x}) \tag{15}$$

in [*a*, *b*].

• Let *f* ∈ *<sup>L</sup>*1(*<sup>a</sup>*, *b*), *α*, *β* > 0, such that *n* − 1 < *α* ≤ *n*, *m* − 1 < *β* ≤ *m* (*<sup>n</sup>*, *m* ∈ N) and *α* + *β* < *n*, *fm*−*<sup>α</sup>* = *a Im*−*<sup>α</sup> x f* ∈ *ACm*([*<sup>a</sup>*, *b*]). Then:

$$({}\_{a}D\_{x}^{a}D\_{x}^{\mathcal{S}}f)({}\_{a}) = ({}\_{a}D\_{x}^{a+\beta}f)({}\_{a}) - \sum\_{j=1}^{m}({}\_{a}D\_{x}^{\mathcal{S}-j}f)({}\_{a} + )\frac{({}\_{a} - a)^{-j-a}}{\Gamma(1-j-a)}.\tag{16}$$

• Let *f* ∈ *<sup>L</sup>*1(*<sup>a</sup>*, *b*), *α* ≥ *β* > 0. Then:

$$(\_{a}D\_{\mathbf{x}^{a}}^{a}I\_{\mathbf{x}}^{\beta}f)(\mathbf{x}) \ = \ \_{a}(\_{a}I\_{\mathbf{x}}^{\beta-a}f)(\mathbf{x}), \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \$$

$$({}\_aD\_{\mathbf{x}^a}^{\mathbf{a}}I\_{\mathbf{x}}^{\beta}f)(\mathbf{x}) \quad = \quad ({}\_aD\_{\mathbf{x}}^{\mathbf{a}-\beta}f)(\mathbf{x}), \quad {}\_a\mathbf{a}\geq\beta \tag{18}$$

Related to the integrals of Riemann–Liouville, the definition of Caputo fractional derivative appears:

• Left-side Caputo Fractional Derivative of order *α* > 0:

$${}^{C}\_{a}D\_{x}^{a}\phi(\mathbf{x}) = {}\_{a}D\_{x}^{a}\left(\phi(\mathbf{x}) - \sum\_{k=0}^{n-1} \frac{\phi^{(k)}(a)}{k!}(\mathbf{x} - a)^{k}\right)$$

$$= \frac{1}{\Gamma(n-\alpha)} \int\_{a}^{\mathbf{x}} \frac{\phi^{(n)}(t)}{(\mathbf{x}-t)^{a+1-n}}dt, \qquad \mathbf{x} > a. \tag{19}$$

• Right-side Caputo Fractional Derivative of order *α* < 0:

$$\, \_x^C D\_b^a \phi(\mathbf{x}) = \frac{(-1)^n}{\Gamma(n-a)} \int\_x^b \frac{\phi^{(n)}(t)}{(\mathbf{x}-t)^{a+1-n}} dt, \qquad \mathbf{x} < b. \tag{20}$$

We have, as before, *n* ∈ N such that 0 ≤ *n* − 1 < *α* < *n*, and now the *n* + 1 derivatives of function *φ* must be continuous and bounded in [*a*, *b*].

The following identity established the relation between the Riemann–Liouville and Caputo fractional derivatives, for *f* a suitable function (for instance, *f* n-derivable):

$$(({}\_{a}D\_{x}^{a}f)({\bf x}) = ({}\_{a}^{C}D\_{x}^{a}f)({\bf x}) + \sum\_{j=0}^{n-1} \frac{f^{(j)}({\bf a})}{\Gamma(1+j-a)}({\bf x}-{\bf a})^{j-a}.\tag{21}$$

The extension of Riemann–Liouville fractional operators to infinity intervals leads to the Liouville fractional operators:

• Left-side Liouville Fractional Integral of order *α* > 0:

$$\int\_{+}^{L} I\_{x}^{\mathbf{x}} \phi(\mathbf{x}) = \frac{1}{\Gamma(a)} \int\_{-\infty}^{\mathbf{x}} (\mathbf{x} - t)^{a - 1} \phi(t) dt, \qquad \mathbf{x} \in \mathbb{R}. \tag{22}$$

• Right-side Liouville Fractional Integral of order *α* < 0:

$$\int\_{\mathcal{X}} \mathbf{I}\_{-}^{\alpha} \boldsymbol{\phi}(\mathbf{x}) = \frac{1}{\Gamma(\alpha)} \int\_{\mathcal{X}} ^{\infty} (\mathbf{x} - t)^{\alpha - 1} \boldsymbol{\phi}(t) dt, \qquad \mathbf{x} \in \mathbb{R}. \tag{23}$$

• Left-side Liouville Fractional Derivative of order *α* > 0:

$$\frac{1}{\tau^L\_+} D\_x^a \phi(\mathbf{x}) = \frac{1}{\Gamma(n-a)} \left(\frac{\partial}{\partial \mathbf{x}}\right)^n \int\_{-\infty}^{\mathbf{x}} (\mathbf{x} - t)^{a - (n-1)} \phi(\mathbf{t}) dt = D^n(+l\_x^{n-a} f)(\mathbf{x}), \qquad \mathbf{x} \in \mathbb{R}. \tag{24}$$

• Right-side Liouville Fractional Derivative of order *α* < 0:

$$\mathbf{u}\_x^L D\_-^a \phi(\mathbf{x}) = \frac{1}{\Gamma(n-a)} \left(-\frac{\partial}{\partial \mathbf{x}}\right)^n \int\_x^\infty (\mathbf{x} - t)^{a - (n-1)} \phi(t) dt = (-D)^n (\mathbf{x} I\_-^{n-a} f)(\mathbf{x}), \qquad \mathbf{x} \in \mathbb{R}. \tag{25}$$

In all cases *n* ∈ N, such that 0 ≤ *n* − 1 < *α* < *n*.

#### *1.3. Mittag-Leffler Functions*

Special functions related to the eigenfunctions of fractional operators are the Mittag-Leffler functions. They appear in the solution of many fractional differential equations. The Mittag-Leffler functions are generalizations of the exponential function and they was introduced by the mathematician G.M. Mittag-Leffler in 1903:

$$\begin{aligned} E\_a(t) &= \sum\_{k=0}^{\infty} \frac{t^k}{\Gamma(ak+1)} \qquad (a > 0, a \in \mathbb{R}),\\ E\_{a,\beta}(t) &= \sum\_{k=0}^{\infty} \frac{t^k}{\Gamma(ak+\beta)} \qquad (a, \beta > 0, a, \beta \in \mathbb{R}).\end{aligned} \tag{26}$$

Here, we have some elementary properties of the Mittag-Leffer function. For some values of the parameters *α*, *β*, Mittag-Leffer functions return known classical functions, for example:

$$E\_1(t) = e^t,\tag{27}$$

$$E\_2(t) = \cosh(\sqrt{t}).\tag{28}$$

The relevance of the Mittag-Lefler functions is their behavior as generalized exponential functions associated to the Riemann–Liouville and Caputo fractional derivatives:

$$
\lambda\_0 D\_t^a t^{a-1} E\_{a,a}(\lambda t^a) = \lambda t^{a-1} E\_{a,a}(\lambda t^a), \tag{29}
$$

and

$$
\lambda\_0^\mathbb{C} D\_t^a E\_\mathfrak{a}(\lambda t^a) = \lambda E\_\mathfrak{a}(\lambda t^a). \tag{30}
$$

On the other hand, a very useful extension of the Mittag-Leffer function is related to extend the fractional derivatives defined by using a Mittag-Leffler kernel which is non-local and non-singular ([8,9]).

$${}^{ABR}D\_{a+}^{a}f(t) = \frac{B(a)}{1-a} \frac{d}{dt} \int\_{a}^{t} f(\mathbf{x}) E\_{\mathbf{a}}\left(\frac{-a}{1-a}(t-\mathbf{x})^{a}\right) d\mathbf{x},\tag{31}$$

$${}^{ABC}D\_{a+}^{a}f(t) = \frac{B(a)}{1-a} \int\_{a}^{t} f'(\mathbf{x}) E\_{a}\left(\frac{-a}{1-a}(t-\mathbf{x})^{a}\right) d\mathbf{x} \tag{32}$$

valid for 0 < *α* < 1, with *<sup>B</sup>*(*α*) being a normalisation function. This new definition has many applications at the same time that satisfies the extensions of the product rule and chain rule.

Futhermore, these functions have many uses, for instance, they allow to address fractal kinetics from basic functions ([10]) or they can be used as a simplification tool combined to exponent/powerlaw mathematical congruence ([11]).

#### *1.4. Some Ideas for Numerical Integration*

It is simple to extend some classical methods from integer to non-integer orders. For instance, the most basic first order explicit method for numerical integration of ordinary differential equations with a given initial value, known as Euler methods.

In classical dynamical models, symplectic integration schemes preserve the flow of the hamiltonian while other classical integrators as Runge–Kutta schemes do not necessarily conserve it. This is immediately translated into the conservation of the first integral of motion and long term stability of the scheme. In fractional dynamics, energy is not conserved. Despite this, using a symplectic scheme in fractional mechanics ensures that the observed instabilities will be certainly due to the fractional operators. Long-term stability is inherited in this fractional mapping as is shown below.

Let us consider the following initial value problem with Caputo fractional derivative:

$$\,\_0^C D\_t^a y(t) = y^a, \qquad t \in [0, T], \quad 0 < a < 1 \tag{33}$$

$$y(0) = y\_0 \tag{34}$$

where *yα* is a Lipschitz function with constant *L*.

By using a truncated Taylor series for the Caputo fractional derivative, with step *tn* = *nh* with *h* = *TN* y *n* = 0, 1, 2, . . . , *N*, we obtain:

$$y(t+h) = y(t) + h^a \frac{^C\_0 D\_t^a y(t)}{\Gamma(a)}\tag{35}$$

and then the fractional Euler method for this initial value problem is:

$$y\_{n+1} = y\_n + \frac{h^\alpha}{\Gamma(\alpha)} f(t\_n, y\_n) \tag{36}$$

with convergence order *h<sup>α</sup>*.

Other more complex generalizations of classical methods are possible. For instance, fractional Hamiltonian–Jacobi methods. The equations of motions for a one-dimensional Hamiltonian system *H* = 12*p*2 + *<sup>V</sup>*(*x*) with unit mass, are defined as

$$\dot{\mathfrak{x}} = \frac{\partial H}{\partial p} = p, \qquad \dot{p} = -\frac{\partial H}{\partial \mathfrak{x}} = -V'(\mathfrak{x}), \tag{37}$$

that is associated to the second order equation *x*¨ + *<sup>V</sup>*(*x*) = 0. This system can be generalized by using Caputo time-fractional derivatives

$$\begin{cases} \,^{\mathbb{C}}\_{0}D\_{t}^{\mathbf{x}}\mathbf{x} = \mathbf{p},\\ \,^{\mathbb{C}}\_{0}D\_{t}^{\mathbf{x}}\mathbf{p} = -V'(\mathbf{x}),\end{cases}\tag{38}$$

where 0 < *α* ≤ 1.

And the system Equation (38) with initial condition *x*(0) = *x*0, *p*(0) = *p*0 is equivalent to

$$\begin{cases} \mathbf{x}(t) = \mathbf{x}(0) + \mathbf{o}I\_t^\mathbf{a} p(t) \\ p(t) = p(0) - \mathbf{o}I\_t^\mathbf{a} V'(\mathbf{x}(t)) \end{cases} \tag{39}$$

For the numeric solution of system Equation (39) we have developed a map (see [12])

$$\begin{cases} p\_n = p\_0 - \frac{(\Delta t)^a}{\Gamma(a+1)} \sum\_{k=0}^{n-1} V'(\mathbf{x}\_k)[(n-k)^a - (n-k-1)^a] \\ \mathbf{x}\_n = \mathbf{x}\_0 + \frac{(\Delta t)^a}{\Gamma(a+1)} \sum\_{k=0}^{n-1} p\_{k+1} [(n-k)^a - (n-k-1)^a] \end{cases} \tag{40}$$

When *α* = 1, this is equivalent to a second order symplectic integrator *pn* = *pn*−<sup>1</sup> − Δ*t <sup>V</sup>*(*xn*−<sup>1</sup>), *xn* = *xn*−1 + Δ*t pn*, and mapping Equation (38) provides an orbit (*xn*, *pn*) approaching the exact orbit at *t* = *n*Δ*t* when Δ*t* → 0. Futhermore, the term *pk*+<sup>1</sup> can be replaced by *pk* in order to return the second order Euler scheme which is not symplectic as Δ*t* → 0.

The orbit at step *n* depends on all the previous states up to the initial one due to the memory kernel of the fractional integral and then we have an infinite dimensional mapping. The computational complexity of the orbit up to (*xn*, *pn*) is of order *n*<sup>2</sup> whereas it is of order *n* for *α* = 1.

This map has been tested taking *α* = 1 with standard models and using different potentials which solutions are known. In particular for the harmonic oscillator and initial conditions *x*0 = 1, *p*0 = 0 and Δ*t* = 0.01, solution has been compared with *x*(*t*) = cos(*t*), *p*(*t*) = sin(*t*) with error smaller than 0.5% after 10,000 steps.

The map Equation (40) has been used to simulate numerical solutions to significant non-linear fractional generalized Hamiltonian problems with a potential *<sup>V</sup>*(*x*), where the explicit solution cannot be found. For instance, standard academic cases like free particle motion (*V* = 0) and a uniformly accelerated particle (*V* = *kx*) [13], or the simple oscillator *V* = 12*ω*2*x*2, the double well potential *V* = 14 *x*4 − 12 *x*2 (see Figure 1) and the pendulum *V* = cos(*x*) [12]. In this pioneer numerical work, it was observed that the fractional derivative introduces a damping effect which can be either algebraic or exponential depending on the time scales of the system. It is considered in a more general context in [10] or in [11].

Riemann–Liouville time fractional problem could be integrated in this way, changing Riemann– Liouville derivative by Caputo's and applying a similar mapping.

**Figure 1.** Phase portrait for the nonlinear oscillator with double well potential *V* = 1 4 *x*4 − 1 2 *x*2 (corresponding to the time interval 0 ≤ *t* ≤ 20).

Non-homogeneous systems can be studied through the generalizations of this map. For instance, by the introduction of an external force *f*(*t*) to simulate a forced-damped oscillator [14]:

$$\begin{cases} \,^{\mathbb{C}}\_{0}D\_{t}^{a}\mathbf{x} = p, \\ \,^{\mathbb{C}}\_{0}D\_{t}^{a}p = -\omega\_{0}^{2}\mathbf{x} + f(t), \end{cases} \tag{41}$$

and this change means the introduction of an extra force term *fn* = *f*(*n*Δ*t*) in the first equation of (40).

In particular, the system evolves to a limit cycle for an harmonic forcing ( *f*(*t*) = *A*0 cos(*ω<sup>t</sup>*)), similar to the classic case with the forced-damped oscillator. Varying the forcing frequency *ω* a resonance motion is reproduced, with amplitude

$$A\_{\rm res} = \frac{A\_0}{2\omega\_0} \left| \frac{1}{(i\omega)^a - i\omega} \right| \tag{42}$$

(see Figure 2).

**Figure 2.** Plot of the amplitude of the limit cycle of *x*(*t*) for different forcing frequencies *ω* versus the theoretical amplitude Equation (42) with *α* = 0.95.

#### **2. Variational Problems and Euler—Lagrange Equations**

#### *2.1. Nonlocal, Fractional Calculus of Variations*

Classical mechanics can be viewed as founded on Hamilton's principle of stationary action. It is an elegant theory based on a very simple axiom. Some authors, appealed both by this and by the potential effectiveness of fractional modeling, have developed a fractional mechanics.

The guiding idea is to keep the same axiom (Hamilton's principle) and allow to have fractional derivatives as variables. Apparently, the founding papers are due to Riewe [15,16]. He uses fractional derivatives and builds the corresponding Euler–Lagrange equations in a systematic way. He also gives the corresponding Hamilton equations. As an application, he provides a formulation that includes a dissipative force. His presentation is not without some limitation.

Agrawall and other authors have generalized from them a Lagrangian and a Hamiltonian formalism [17–19] that also includes constraints.

From these references we see that right and left fractional derivatives appear in the Lagrangian and in the fractional Euler–Lagrange equations. The implications that this poses to the causality have not been dealt with.

We propose a slightly more general formalism that allows systems with different nonlocal terms in the Lagrangian, including fractional integrals and fractional derivatives of Caputo and Riemann–Liouville kind.
