# Fractional derivative truncation approximation for real-time applications

Jean-François Duhé, Stéphane Victor, Pierre Melchior, Youssef Abdelmounen, François Roubertie

<figure><img src="image_1.png" /></figure>

## To cite this version:

Jean-François Duhé, Stéphane Victor, Pierre Melchior, Youssef Abdelmounen, François Roubertie. Fractional derivative truncation approximation for real-time applications. Communications in Non-linear Science and Numerical Simulation, 2023, 119, pp.107096. 10.1016/j.cnsns.2023.107096. hal-04317777

HAL Id: hal-04317777

https://hal.science/hal-04317777v1

Submitted on 8 Jan 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution - NonCommercial 4.0 International License
Fractional derivative truncation approximation for real-time applications

Jean-François Duhé<sup>a,c</sup>, Stéphane Victor<sup>a,*</sup>, Pierre Melchior<sup>a</sup>, Youssef Abdelmounen<sup>b,c</sup>, François Roubertie<sup>b,c</sup>

**Abstract**

Fractional derivatives are non local operators as they are well-suited for modeling long-memory phenomena such as heat transfers or fluid dynamics. However, such non-locality implies a constant knowledge of the full past of the functions to differentiate. In the context of real-time (online) system identification, such a global property may limit the implementation as calculations can become slower as time grows. Also, instead of using the full data length to compute the fractional order derivatives, a truncated version of fractional derivative is proposed. However, such a truncation windowing will bring inaccuracy on the computation of the fractional derivative. This study deals with the relationship between the signal frequency content, the fractional derivative truncated approximation as well as the fractional system relaxation in order to get a good approximation of the truncated fractional derivative and consequently, to help providing real-time exploitable system identification algorithms that uses such truncated window fractional derivatives.

**Keywords:** fractional derivative, short memory principle, truncation window, fractional order system, Mittag-Leffler function, fractional relaxation, LMRPEM, system identification

# 1. Introduction

Fractional calculus has been a mere mathematical concept for many centuries. However, during the current and past centuries its properties have proven to be useful for diffusive system modeling (Oustaloup (2014)), as well as robust control design (Shah and Agashe (2016); Oustaloup et al. (1993)). Difusive phenomena have proven to be analytically well expressed and well modeled through fractional derivatives. Thermal impedances at middle and high frequencies has also been proven to behave like a half-order integrator (Malti et al. (2009); Battaglia et al. (2001)). Its properties have proven to be useful in medical scenarios, as is the case of lung impedance modeling (Ionescu et al. (2011); Ionescu et al. (2014); Victor et al. (2020); Duhé et al. (2022)), cardiac tissue (Magin and Ovadia (2006)) and muscle relaxation (Sommacal et al. (2006); Melchior et al. (2019)). Therefore, there is an increasing interest regarding fractional order models and their properties.

Even though fractional differentiation definition is not unique (Kilbas et al. (2006); Garrappa et al. (2019); Samko et al. (1993)), all of the proposed definitions share the same long-memory property to model diffusive phenomena. The past history of a given function *f*(*t*) may still have a notorious influence on its fractional derivative after a long period of time. Fractional

order system dynamics can be interpreted as an infinite distribution of time constants (Oustaloup (1995, 2014)). Unfortunately, this implies the existence of infinitely slow time constants. A main drawback regarding this property is the non local nature of fractional differentiation: the whole knowledge of the function past is required in order to well compute its fractional derivative at any given time. As time increases, a signal history becomes progressively longer, thus leading to slower computation times for computing a fractional derivative. Furthermore, computation time may be critical for real-time (online) applications and a continuous increase in the computation time can lead to infeasible algorithms. As a consequence, an appropriate moving window in time is sought so that the (truncated) approximation of the fractional operator may be accurate enough to still make online implementable algorithms.

The fractional derivative may be elegantly truncated in a really simple way through Podlubny's short memory principle (Podlubny (1999)). The principle allows to use a truncated time window of a function past and still get an accurate estimation of the fractional derivative. Even though Podlubny's principle provides an accuracy guarantee, it gives a pessimistic limit that could still be too slow for implementable real-time scenarios. In system identification, experiments usually imply an input signal exciting a defined bandwidth (see Ljung (1999)), in order to be sufficiently persistent (rich in frequency). Therefore, the spectrum of the input signal may influence the required time window of the past signal and its associated truncation error. One of the aims of this paper is to analyze a signal spectrum and the differentiation order influence on the truncation error in order the design a correct time window so that the fractional deriva-

<footnote><sup>a</sup>Univ. Bordeaux, CNRS, IMS UMR 5218, Bordeaux INP/ENSEIRB-MATMECA,
351 Cours de la Libération, 33405 Talence CEDEX, France</footnote>

<footnote><sup>b</sup>IHU Liryc, Electrophysiology and Heart Modeling Institute, Fondation Bordeaux Université, 33000 Bordeaux, France</footnote>

<footnote><sup>c</sup>INSERM, Centre de recherche Cardio-Thoracique de Bordeaux, 33 000 Bordeaux, France</footnote>

<footnote>*Corresponding author
Email addresses: jean-francois.duhe@ims-bordeaux.fr
(Jean-François Duhé), stephane.victor@ims-bordeaux.fr (Stéphane Victor), pierre.melchior@ims-bordeaux.fr (Pierre Melchior), dr.a.youssef@free.fr (Youssef Abdelmounen), rouberf@me.com
(François Roubertie)</footnote>
tive approximation on a time window gets sufficiently precise enough. This will be performed through a simple scenario.

Besides the fractional differentiation long-memory property, another issue may arise when dealing with real-time fractional order system identification: output estimation at a given time. The fractional order model non-locality and the presence of noise in the recorded data may lead to prefer an estimation of the output based upon an auxiliary model such as the instrumental variable method. The use of a truncated auxiliary model can lead to an initial condition problem for simulation. It has been proven that fractional order systems can be interpreted as an infinite state representation (Trigeassou and Maamri (2011, 2020)). Formal methods to estimate the response of a fractional order system with initial conditions are available. However, their solutions are too complex to be easily computable for real-time scenarios. Therefore, an analysis based on the transient response will be presented for a simple fractional order model of the first kind.

The earliest proposed techniques for fractional order system identification come from the 90s (Le Lay (1998); Oustaloup et al. (1996); Malti et al. (2008a)). These first methods only estimated the coefficients of fractional order transfer functions. However, they were also limited to the Grünwald-Letnikov discrete-time definition of the fractional derivative. By considering continuous-time system identification, methods relying on state variable filters, least squares and instrumental variables have been introduced (Malti et al. (2008b); Victor et al. (2013, 2022c)) as well as validation through experimental applications (Malti et al. (2009); Victor et al. (2013, 2016)). Recently, experiment design for optimal identification has also been studied for fractional order systems (Malti et al. (2022)).

Real-time system identification is based upon recursive estimation methods, such as recursive least-squares or recursive prediction error method (Ljung (1999)). The analysis of recursive identification methods has also been extended to least squares with state variable filters and instrumental variable techniques (Padilla (2017); De Wit (1986)). All such methods have also been adapted and analyzed for continuous-time systems (Garnier and Wang (2008)). Their extensions to fractional order systems have been proposed in Victor et al. (2022b) where the Long-Memory Recursive Prediction Error Method (LMRPEM) algorithm has proven to be a useful tool for accurately estimating all fractional order model parameters (coefficients and differentiation orders). An effective combination of a well-chosen truncation method and the LMRPEM-2 algorithm could lead to an entirely real-time implementable identification technique for fractional order models. Such an algorithm works successfully if the sampling time is not too small so that the computation time at each iteration is far smaller than the sampling time, and also when the experiment does not last too long. Consequently, if the sampling time gets smaller and if the experiments last very long such is the case for open-heart surgeries that can last for several hours, the computation of fractional derivatives should be shortened: therefore, the contribution of this paper is to propose to use a correctly designed time window to get an accurate enough approximation of the truncated fractional derivatives.

The paper is organized as follows: Section 2 introduces the basics of fractional calculus and fractional order transfer functions, Section 3 presents the Short Memory Principle and the frequency study on the truncation error, Section 4 analyzes the fractional order transient response for a simple scenario. Section 5 presents the LMRPEM algorithm and an application of a truncated estimation. Conclusions and perspectives are presented in Section 6.

## 2. Mathematical background

### 2.1. Fractional derivatives

As previously stated, many different definitions for fractional order differentiation exist.

**Definition 2.1** (Riemann and Liouville fractional derivative). The most well-known fractional derivative definition is the one of Riemann and Liouville (Samko et al. (1993)) defined as:

$$D^\nu f(t) = \frac{1}{\Gamma(n-\nu)} \frac{d^n}{dt^n} \int_0^t (t-\tau)^{n-\nu-1} f(\tau) d\tau, \quad (1)$$

where $n-1 \le \nu < n$ with $n$ being an integer and Euler's gamma function is defined as:

$$\Gamma(x) = \int_{0}^{\infty} e^{-t} t^{x-1} dt \quad (2)$$

for $x \in \mathbb{R} \setminus \mathbb{Z}^-$.

If the initial conditions are set to zero, Riemann-Liouville's definition is proven to be equal (Podlubny, 1999, Chap. 7) to the series definition given by Grünwald (1867) and Letnikov (1868).

**Definition 2.2** (Grünwald-Letnikov fractional derivative definition). Another largely used fractional derivative definition is the one of Grünwald-Letnikov defined by:

$$p^{\nu} f(t) = \lim_{h \to 0} \frac{1}{h^{\nu}} \sum_{j=0}^{\lfloor \frac{t}{h} \rfloor} (-1)^j c_j f(t - jh), \quad (3)$$

where $p = \frac{d}{dt}$ stands for the differential operator, $\lfloor \cdot \rfloor$ the floor operator and $c_j$ the generalized Newton's binomial coefficients:

$$c_j = \binom{\nu}{j} = \frac{\Gamma(\nu + 1)}{\Gamma(j + 1)\Gamma(\nu - j + 1)}. \quad (4)$$

Coefficients in equation (4) exhibit a slow tendency towards zero and provide the long-memory behavior for this definition.

If parameter $h$ is replaced by the sampling time $T_s$, definition 2.2 can be seen as a weighted sum of the function past:

$$p^{\nu} f(t) \approx \frac{1}{T_s^{\nu}} \sum_{j=0}^{\lfloor \frac{T_s}{T_s} \rfloor} (-1)^j c_j f(t - jT_s). \quad (5)$$

This definition allows an easy implementation and has led to its significant popularity. Definition (5) will be used for all simulation scenarios in the following.

For null initial conditions, the Laplace transform of the fractional derivative leads to a simple expression:

$$\mathcal{L}\{p^{\nu} f(t)\} = s^{\nu}F(s). \quad (6)$$
## 2.2. Fractional order dynamical model

A single-input single-output (SISO) fractional order model relates its output $y(t)$ to its input $u(t)$ through a fractional order differential equation:

$$y(t) + a_1 p^{\alpha_1} y(t) + \dots + a_{m_A} p^{\alpha_{m_A}} = \\ b_0 p^{\beta_0} u(t) + b_1 p^{\beta_1} u(t) + \dots + b_{m_B} p^{\beta_{m_B}} u(t) \quad (7)$$

where $(a_i, b_j) \in \mathbb{R}^2$ are real numbers and $\alpha_1 < \alpha_2 < \dots < \alpha_{m_A}$ and $\beta_0 < \beta_1 < \dots < \beta_{m_B}$ are non-integer positive numbers.

Applying the Laplace transform to relation (7) provides a fractional order transfer function model:

$$G(s) = \frac{B(s)}{A(s)} = \frac{\sum_{i=0}^{m_B} b_i s^{\beta_i}}{1 + \sum_{j=1}^{m_A} a_j s^{\alpha_j}} \quad (8)$$

**Definition 2.3 (Commensurate order).** A fractional order model is commensurate if all of its differentiation orders are integer multiples of a same order $\mu$, called the commensurate order:

$$G(s) = \frac{B(s)}{A(s)} = \frac{\sum_{i=0}^{m_B} b_i s^{j_i \mu}}{1 + \sum_{j=1}^{m_A} a_j s^{j_j \mu}} \quad (9)$$

## 2.3. Stability of fractional order models

Stability of fractional order systems has been analyzed in different contexts. The most well-known stability theorem was established in Matignon and d'Andréa-Novel (1996) and allows to check the stability of commensurate order systems through the location of its $s^\nu$-poles. The original theorem was established for commensurate orders, $0 < \nu < 1$, but was extended for orders between 1 and 2 in Moze and Sabatier (2005). Orders with commensurate orders beyond 2 can be proven to be unstable (Malti et al. (2011)). Further extensions have been developed in order to check stability of non-commensurate systems (Rivero et al. (2013)). More recently, stability of multivariable and nonlinear fractional systems has also been studied in Lenka (2019).

**Theorem 2.4 (Stability theorem of Matignon).** Let S be a $\mu$-commensurate transfer function, G. $G(s) = \frac{Q_\mu(s)}{P_\mu(s)}$ is bounded-input bounded-output (BIBO) stable if and only if:

$$0 < \mu < 2 \quad (10)$$

and, for every pole $s_k$ ($P_\mu(s_k) = 0$):

$$|\arg(s_k)| > \mu \frac{\pi}{2}. \quad (11)$$

<figure><img src="image_3.png" /></figure>

# 3. Truncated fractional derivative

## 3.1. Short Memory Principle

The Short Memory Principle is based upon Riemann-Liouville fractional derivative definition 2.1. The principle consists in only taking the past of the function in an interval $[t - L, t]$, where L is the time window of the memory.

By expressing the fractional derivative as the sum of a long and a short memory, one gets:

$$D^\nu f(t) = \frac{1}{\Gamma(n-\nu)} \frac{d^n}{dt^n} \left[ \int_0^{t-L} (t-\tau)^{n-\nu-1} f(\tau) d\tau + \int_{t-L}^t (t-\tau)^{n-\nu-1} f(\tau) d\tau \right] \quad (12)$$

By applying the short memory principle, one gets:

$$D^\nu f(t) \approx \frac{1}{\Gamma(n-\nu)} \frac{d^n}{dt^n} \int_{t-L}^{t} (t-\tau)^{n-\nu-1} f(\tau) d\tau. \quad (13)$$

For a signal $f(t)$ with an upper bound, M, and a desired error, $\epsilon$, the required length to guarantee an accurate approximation by using the short memory principle is given by ((Podlubny, 1999, Ch. 7)):

$$L \ge \left[ \frac{M}{\epsilon |\Gamma(1-\nu)|} \right]^{\frac{1}{\nu}}. \quad (14)$$

Even though this expression will guarantee accuracy of the approximation for any given function $f(t)$, it does not take into account the shape of the given function or its frequency domain properties. Therefore, the inferior limit for this memory may be very large.

## 3.2. The cosine function

As periodic functions can be created from sine and cosine signals by means of their Fourier decomposition, let us consider the following signal:

$$f(t) = \cos(\omega t). \quad (15)$$

Its fractional derivative (according to Riemann-Liouville's definition 2.1 is:

$$D^\nu f(t) = \frac{1}{\Gamma(n-\nu)} \frac{d^n}{dt^n} \left[ \int_0^{t-L} \frac{\cos(\omega\tau)}{(t-\tau)^{n+\nu+1}} d\tau + \int_{t-L}^t \frac{\cos(\omega\tau)}{(t-\tau)^{n+\nu+1}} d\tau \right] \quad (16)$$

From (12) and according to Riemann-Liouville's definition 2.1, it is possible to decompose a fractional derivative with $\nu > 1$ like the fractional derivative of an integer one. Therefore, the analysis will be limited to the range $0 \le \nu < 1$ and $n = 1$. Under such assumptions, the fractional derivative becomes:

$$D^\nu f(t) = \frac{1}{\Gamma(1-\nu)} \frac{d}{dt} \left[ \int_0^{t-L} \frac{\cos(\omega\tau)}{(t-\tau)^\nu} d\tau + \int_{t-L}^t \frac{\cos(\omega\tau)}{(t-\tau)^\nu} d\tau \right] \quad (17)$$
By applying the short memory principle (see section 3.1):

$$ \int_{t-L}^{t} D^{\nu} f(t) \approx \frac{1}{\Gamma(1-\nu)} \frac{d}{dt} \int_{t-L}^{t} \frac{\cos(\omega\tau)}{(t-\tau)^{\nu}} d\tau. \quad (18) $$

Let us define an error signal $\epsilon(t)$ as being the difference between the exact derivative (17) and its truncated derivative (18):

$$ \begin{aligned} \epsilon(t) &= \int_{0}^{L} D^{\nu} f(t) - \int_{0}^{L} D^{\nu} f(t) \\ &= \frac{1}{\Gamma(1-\nu)} \frac{d}{dt} \int_{0}^{t-L} \frac{\cos(\omega\tau)}{(t-\tau)^{\nu}} d\tau. \end{aligned} \quad (19) $$

As can be seen, the time variable is included in the limits of the integral. Therefore, to apply the derivative, one uses the following property.

From the primitive $F(t)$

$$ F(t) = \int_{a(t)=0}^{b(t)=t-L} f(t, \tau) d\tau \quad (20) $$

its time derivative is:

$$ \frac{d}{dt} F(t) = f(t, t - L) + \int_{0}^{t-L} \frac{\partial f}{\partial t} d\tau. \quad (21) $$

The error signal (19) becomes:

$$ \epsilon(t) = \frac{1}{\Gamma(1-\nu)} \left[ \begin{aligned} & L^{-\nu} \cos(\omega(t-L)) \\ & -\nu \int_{0}^{t-L} \frac{\cos(\omega\tau)}{(t-\tau)^{\nu+1}} d\tau \end{aligned} \right]. \quad (22) $$

As the cosine is less than 1, an upper bound of the error (22) can be found as:

$$ \epsilon(t) \le \frac{1}{\Gamma(1-\nu)} \left[ L^{-\nu} - \nu \int_{0}^{t-L} \frac{\cos(\omega\tau)}{(t-\tau)^{\nu+1}} d\tau \right]. \quad (23) $$

This expression suggests that the error signal is conditioned by the cosine signal spectrum. However, it should be noted that an analytical expression for the integral term in relation (23) does not exist. Therefore, numerical simulations will be carried out to pursue the truncation error analysis.

### 3.3. Simulations and results

#### 3.3.1. Simulations with the fractional derivative by Grünwald-Letnikov definition

Approximation error for a truncation window $L$. The Grünwald-Letnikov's definition 2.2 being the most commonly used for numerical simulation, its definition is used to provide a truncation window of length $L$:

$$ L p^{\nu} f(t) \approx \frac{1}{T_s^{\nu}} \sum_{j=0}^{\lfloor \frac{L}{T_s} \rfloor} (-1)^j c_j f(t - jT_s). \quad (24) $$

The maximal error is then defined as follows:

$$ \epsilon_{\max} = \max |p^{\nu} f(t) - L p^{\nu} f(t)|. \quad (25) $$

A simple cosine function will be used for simulations. The sampling period is taken as $T_s = 0.01$s for the whole section.

Approximation error for different truncation windows $L$ with $\nu = 0.5$ and $T_s = 0.01$s. Figure 1 shows the maximum error for different frequencies $\omega$ within the band [1-100]rad/s. The maximal error has been computed and plotted for different window lengths $L$ by keeping the differentiation order $\nu = 0.5$ constant. As can be observed, an increasing window length leads to lower errors. However, the signal frequency also plays an important role in the approximation error. As frequency increases, there is a significant error reduction for all tested window length. It could be considered that a higher frequency $\omega$ implies a higher number of periods of the cosine function available for fractional differentiation approximations. Therefore, an error reduction according to the input signal frequency $\omega$ may be considered as a natural consequence of an increased knowledge of the function past.

<figure><img src="image_4.png" /><figcaption>Figure 1: Maximal approximation error vs frequency $\omega$ for different truncation windows $L$ with $\nu = 0.5$ and $T_s = 0.01$s</figcaption></figure>

Approximation error for a same number of periods with $\nu = 0.5$ and $T_s = 0.01$s. As a consequence, such a result suggests that the approximation error stays relatively the same if the same number of periods are taken to approximate its fractional derivative. Therefore, the truncation window $L$ is now taken as:

$$ L = \frac{X}{\omega}, \quad (26) $$

where $X$ is an arbitrary factor and $\omega$ is the input signal frequency. This choice of the truncation window length leads to an approximation based upon a fixed number of the signal periods rather than a time length. The simulation results are presented in figure 2.

Higher factors $X$ do lead to lower approximation errors of the truncated fractional derivative, but the maximal error is not preserved by keeping a constant number of periods. Moreover, as the input signal frequency $\omega$ increases, keeping a same number of periods for a truncated fractional derivative actually leads to an important increase of the maximal error. However, it should be noted that a high enough number of previous periods may severely limit the frequency effect, as is the case for $X = 50$.
<figure><img src="image_5.png" /><figcaption>Figure 2: Maximal approximation error vs frequency ω for different factors X with $v = 0.5$ and $T_s = 0.01s$</figcaption></figure>

An important remark should be done. In figure 1, as the input signal frequency ω increases, the same time length L is kept. This means that more periods of the cosine signal are included beneath the memory length as the signal frequency increases. This increasing in knowledge of the function past naturally leads to an approximation error reduction. On the other hand, figure 2 always includes a fixed number of periods. As the signal frequency increases, the time length L decreases in this case. Keeping an equal number of periods is not enough to avoid a significant increase in the approximation error.

Approximation error for different truncation windows L with $v = 0.5$ and $T_s = 0.001s$. Another interesting analysis is to focus on the higher frequencies and reproduce these simulations with a smaller sampling time $T_s = 0.001s$. The same window lengths are kept and the results are shown in figures 3.

<figure><img src="image_6.png" /><figcaption>Figure 3: Maximal approximation error vs the signal frequency ω for different truncation windows L and $T_s = 0.001s$</figcaption></figure>

As can be seen, amplitudes of maximum error are reduced for both cases even though the lengths are the same as in the

<figure><img src="image_7.png" /><figcaption>Figure 4: Maximal approximation error vs frequency ω for different factors X and $T_s = 0.001s$</figcaption></figure>

preceding example. This implies that there is a non-negligible part of the error due to sampling time. On the other hand, the behaviour as frequency increases stays the same.

Approximation error for different differentiation orders $v$ and $T_s = 0.01s$. Let us now analyze the influence of the differentiation order $v$. By taking a constant $ω = 1rad/s$, the maximal error is estimated for truncated fractional derivatives for differentiation orders in [0.1, 2.0) and different X factors. The simulation results are plotted in figure 5.

<figure><img src="image_8.png" /><figcaption>Figure 5: Maximal error vs differentiation order $v$ for different factors X</figcaption></figure>

As expected, higher X factors lead to lower approximation errors. On the other hand, it can also be seen that the errors are minimal as the differentiation order $v$ tends to an integer ($v = 0, 1$ or $2$). This aspect stems from the local property of classic (integer) derivatives as only one (respectively two) instant of the past is necessary for the computation of a derivative of order 1 (respectively order 2). The maximal approximation errors however, are non-trivial. For lower X factors, the maximal approximation errors get closer and are obtained for fractional orders
in [0.4 – 0.6]. The error curves become skewed to the left as the truncation window length L increases, leading to maximal errors within the [0.2 – 0.4] range.

### 3.3.2. Simulations with the fractional derivative by using the Oustaloup method

Recall of the Oustaloup method. Oustaloup’s method (see Oustaloup (1995)) is a well known approximation of the fractional derivative synthesized in the frequency domain. First of all, a fractional order derivative is truncated in the frequency domain. A low ($\omega_b$) and a high ($\omega_h$) frequency corners are chosen and the approximation is valid within this frequency range. The frequency response of a truncated fractional derivative is given by the following expression:

$$D(j\omega) = \left( C_0 \frac{1 + \frac{\omega}{\omega_b}}{1 + \frac{\omega}{\omega_h}} \right)^{\gamma} \quad (27)$$

By applying recursivity properties, this frequency-truncated derivative may be approximated by a series of N integer poles and zeros:

$$D(j\omega) = \lim_{N \to \infty} D_N(j\omega). \quad (28)$$

The integer approximation yields:

$$D_N(j\omega) = \left(\frac{\omega_u}{\omega_h}\right)^\gamma \prod_{k=1}^{N} \frac{1 + \frac{j\omega}{\omega_k}}{1 + \frac{j\omega}{\omega_k}}. \quad (29)$$

Two recursive factors $\alpha$ and $\eta$ allow to estimate the frequencies for both the poles and zeros:

$$\begin{align}
\left(\frac{\omega_u}{\omega_h}\right)^\gamma &= \left(\frac{\omega_b}{\omega_h}\right)^\gamma = \frac{1}{\alpha^{N+\frac{1}{2}}} \nonumber \\
\omega_0' &= \alpha^{-\frac{1}{2}} \omega_u \quad \omega_0 = \alpha^{\frac{1}{2}} \omega_u \nonumber \\
\frac{\omega_{k+1}'}{\omega_k'} &= \frac{\omega_{k+1}}{\omega_k} = \alpha \eta > 1 \nonumber \\
\frac{\omega_{k+1}}{\omega_k'} &= \alpha > 0 \nonumber \\
\frac{\omega_{k+1}}{\omega_k'} &= \eta > 0 \nonumber \\
N &= \frac{\log(\alpha)}{\log(\alpha\eta)} \tag{30}
\end{align}$$

Approximation error for a truncation window L. A time truncation of Oustaloup derivative is now compared with the results of the previous section. This was done by calculating the convolution product between Oustaloup derivative and the input signal for a given length window L:

$$p_{\text{oust}}^{\gamma} = D_N(p)u(t) \quad (31)$$

For $t^* - L < t < t^*$ and $t^*$ being defined as the present time. The maximal approximation error is then defined as:

$$\epsilon_{\max} = \max |p_{\text{oust}}^{\gamma} f(t) - L p^{\gamma} f(t)|. \quad (32)$$

Approximation error for different truncation windows L with $v = 0.5$ and $T_s = 0.01s$. Figures (6) shows the maximum error for different frequencies $\omega$ within the band [1–100]rad/s. In this case, the corner frequencies were chosen as $\omega_b = 0.001$ rad/s and $\omega_h = 1000$ rad/s. The number of cells $N = 20$. The maximal error has been computed and plotted for different window lengths L by keeping the differentiation order $v = 0.5$ constant. As can be noticed, for a fixed memory length errors on the low-frequency limit are similar to those obtained by Grünwald-Letnikov’s method presented in the previous subsection. On the range of [1–10] rad/s there is a reduction in the approximation error. However, the approximation error increases in a significant manner as the signal frequency tends to their higher limit of 100 rad/s. Another interesting remark is that the memory length doesn’t seem to have an important influence on the approximation error on the higher frequency limit, as all curves almost coincide.

<figure><img src="image_9.png" /><figcaption>Figure 6: Maximal error vs frequency $\omega$ for different truncation windows L by using Oustaloup method with $v = 0.5$, $T_s = 0.01s$, $\omega_b = 0.001$ rad/s and $\omega_h = 1000$ rad/s</figcaption></figure>

Approximation error for a same number of periods with $v = 0.5$ and $T_s = 0.01s$. For a fixed number of signal periods, figure 7 shows that the influence of the X factor is negligible for this case on the whole frequency range. For both figures, maximum error at high frequency is considerably higher than for the Grünwald-Letnikov’s method.

As can be noticed, for a fixed memory length errors on the low-frequency limit are similar to those obtained by Grünwald-Letnikov’s method presented in the previous subsection. On the range of [1–10] rad/s there is a reduction in the approximation error. However, the approximation error increases in a significant manner as the signal frequency tends to their higher limit of 100 rad/s. Another interesting remark is that the memory length doesn’t seem to have an important influence on the approximation error on the higher frequency limit, as all curves almost coincide. Figure 7 shows that the influence of the X factor is negligible for this case on the whole frequency range. For both figures, maximum error at high frequency is considerably higher than for the Grünwald-Letnikov’s method.
<figure><img src="image_10.png" /><figcaption>Figure 7: Maximal error vs frequency ω for different factors X with Oustaloup method</figcaption></figure>

A possible reason for the approximation error increase in high frequency could be a high frequency $ω_h$ that is fixed too low for the Oustaloup method. The same simulations are now performed by changing the corner frequencies to $ω_b = 0.0001 \text{ rad/s}$ and $ω_h = 10000 \text{ rad/s}$ and by increasing the number of cells to $N = 30$. Results are shown in figures 8 and 9.

<figure><img src="image_11.png" /><figcaption>Figure 8: Maximal error vs frequency ω for different truncation windows L with Oustaloup's method with an increased frequency range; $ω_b = 0.0001 \text{ rad/s}$ and $ω_h = 10000 \text{ rad/s}$</figcaption></figure>

Even though there is effectively some reduction in the approximation error when enlarging the bandwidth, it seems relatively limited if one considers that a whole additional decade was taken into account for higher frequencies. The global curve tendency stays the same as in the previous simulations.

<figure><img src="image_12.png" /><figcaption>Figure 9: Maximal error vs frequency ω for different factors X with Oustaloup's method with an increased frequency range, $ω_b = 0.0001 \text{ rad/s}$ and $ω_h = 10000 \text{ rad/s}$</figcaption></figure>

Approximation error for different truncation windows L with $v = 0.5$ and $T_s = 0.001$ s. Another factor that may have an impact on the approximation error concerns the sampling time. Figures 10 and 11 show the maximum errors for a sampling time of $T_s = 0.001$ s.

<figure><img src="image_13.png" /><figcaption>Figure 10: Maximal error vs frequency ω for different truncation windows L and $T_s = 0.001$ s with Oustaloup's method</figcaption></figure>

Just as for Grünwald method, there is a significant decrease in the approximation error amplitude. Just as in the previous cases, the length of truncation window seems to have little impact on the approximation error and both curves exhibit an error increase towards the higher frequency limit. Both figures show how a smaller sampling time helps to get a higher frequency range with an acceptable approximation error, as in both cases the error continues to decrease for most of the frequency range.

Given the results of this section, from now on, all simulations will be carried out by considering the Grünwald derivative in order to lower the approximation errors when considering the truncated window fractional derivatives.
<figure><img src="image_14.png" /><figcaption>Figure 11: Maximal error vs frequency ω for different factors X and T_s = 0.001 s with Oustaloup's method</figcaption></figure>

## 4. Fractional order system relaxation

### 4.1. Output estimation in system identification: the initial condition problem

System identification algorithms require estimating the system output $\hat{y}$. From the LMRPEM-2 method (details are provided in Duhe et al. (2021) for the coefficient estimation and in Victor et al. (2022b) for both coefficient and differentiation order estimation and briefly recalled in section 5.1), the estimated output is obtained from the following convolution:

$$ \hat{y}(t) = G(p)u(t). \qquad (33) $$

Such an online algorithm has proven to provide estimates without bias and with very low variance towards high noise level. Fractional order operators have the long-memory property, meaning that the past of the function is needed in order to have its good characterization. Therefore, the whole convolution in (33) is required, starting from $t = 0$ to the present time $t$, at each iteration. As a consequence, when the time gets large enough, the system identification algorithm would progressively become slower as the computational cost will increase due to data growth, and the real-time estimation might get impossible to exploit in spite of its recursive formulation. This imposes a need to obtain an accurate output estimation from a suitable time truncation window so that the fractional derivatives may be accurate enough.

For a truncation window $L$, the estimated output calculations are based upon the assumption of zero initial conditions. At a given time $t$, the output is estimated by assuming that at $t - L$, the system output is null, which is not necessarily true. Moreover, the non-local property of fractional derivatives implies a long-lasting influence of the system past.

In order to analyze such a challenge, a simple yet powerful property of fractional order operators may be applied: linearity. Figure 12 shows a simple block diagram for a system with transfer function $G(s)$ and two different inputs $u_1$ and $u_2$.

<figure><img src="image_15.png" /><figcaption>Figure 12: Fractional order system with two inputs</figcaption></figure>

Assuming that $u_1$ and $u_2$ respectively are the recent and distant past of a global input function $u$ in terms of current time $t^*$, then:

$$ u(t) = \begin{cases} u_1(t) & \text{for } 0 < t < t^* - L \\ u_2(t) & \text{for } t^* - L < t < t^*. \end{cases} \qquad (34) $$

The system output can therefore be divided into the recent and distant past by means of the linearity property, namely:

$$ y(t) = G(p)u_1(t) + G(p)u_2(t) = y_{old}(t) + y_{new}(t). \qquad (35) $$

If $y_{old}$ is negligible at the current time $t^*$, it may be possible to successfully estimate a fractional order model output. For a stable system and a sufficient truncation window $L$, one may then consider that the error related to an inaccurate initialization will be small. Therefore, there may be an interest in the analysis of the neglected term in order to better understand the choice of the window $L$.

### 4.2. Transient response: an elementary review

Let be an elementary fractional order model defined by:

$$ G(s) = \frac{1}{\lambda + s^{\nu}}. \qquad (36) $$

For null initial conditions, the step response in the Laplace domain writes:

$$ C(s) = \frac{1}{s(\lambda + s^{\nu})} = \frac{1}{\lambda} \left[ \frac{1}{s} + \frac{s^{\nu-1}}{s^{\nu} + \lambda} \right], \qquad (37) $$

or, in time domain, gives:

$$ c(t) = \frac{1}{\lambda} [1 - E_{\nu}(-\lambda t^{\nu})] u(t) \qquad (38) $$

where $E_ν$ is the well-known Mittag-Leffler function defined by:

$$ E_{\nu}(z) = \sum_{k=0}^{\infty} \frac{z^{\nu}}{\Gamma(k\nu + 1)}, \quad z \in \mathbb{C}. \qquad (39) $$

This function provides the basis for an elementary transient response of fractional order systems.

### 4.3. Fractional settling time and approximation of the Mittag-Leffler function

The settling time may be a useful reference to analyze the time response of a fractional order system. For integer order systems, charging and relaxation times are related.
System (36) is considered as settled if its response is within a tolerance $\delta$ from the steady state value:

$$1 - E_v(-\lambda t^\nu) = 1 - \delta. \quad (40)$$

Theoretically, the settling time $t_s$ satisfies:

$$E_v(-\lambda t_s^\nu) = \delta. \quad (41)$$

In order to have an approximated solution to the former equation, Mittag-Leffler's function needs to be approximated. Short and long time approximations are available in the literature (Mainardi et al. (2006)) for the case $\lambda = 1$. We propose to take the well-known Mainardi approximation for a general $\lambda$ parameter.

By developing the Mittag-Leffler function from its definition, one gets:

$$E_v(-\lambda t^\nu) = 1 + \frac{(\lambda t^\nu)}{\Gamma(\nu + 1)} + \frac{(\lambda t^\nu)^2}{\Gamma(2\nu + 1)} + \dots \quad (42)$$

Recall the series expansion of the exponential function:

$$e^x \approx 1 + x + \frac{x^2}{2!} + \dots \quad (43)$$

If one retains only the first term for both expansions, Mittag-Leffler function could be approximated by an exponential when $t \to 0$:

$$E_v(-\lambda t^\nu) \approx e^{-\frac{\lambda t^\nu}{\Gamma(1 + \nu)}}. \quad (44)$$

Therefore, the settling time may be approximated by:

$$t_{set}^{short} \approx - \frac{[\ln(\delta)\Gamma(1 + \nu)]^{\frac{1}{\nu}}}{\lambda}. \quad (45)$$

This approximation is called the short time approximation and is limited to cases with relatively fast settling times. However, fractional calculus is commonly used in order to take into account slow processes such as fluid mechanics (as the porous dyke relaxation (Oustaloup (2014))) or in bio-engineering scenarios (Magin (2004); Victor et al. (2020)). As a consequence, an approximation based upon infinitely long times is sought. Mainardi provides the following:

$$E_v(-z^\nu) \approx \frac{z^{-\nu}}{\Gamma(1 - \nu)}. \quad (46)$$

By taking $z = at$ and $a = \lambda^{\frac{1}{\nu}}$, the Mittag-Leffler's function for $t \to \infty$ yields:

$$E_v(-\lambda t^\nu) \approx \frac{t^{-\nu}}{\lambda \Gamma(1 - \nu)}. \quad (47)$$

This leads to the following settling time approximation:

$$t_{set}^{long} \approx [\delta \lambda \Gamma(1 - \nu)]^{-\frac{1}{\nu}}. \quad (48)$$

### 4.4. Charge and discharge of a fractional order system of the first kind

Suppose the following fractional order system of the first kind:

$$G(s) = \frac{1}{s^{0.5} + 1}. \quad (49)$$

Figure 13 is intended to illustrate how simulations are performed. By taking zero initial conditions, a unit step of variable time length $L$ is applied as the input. The sampling time is $T_s = 1s$ for this section. From time $L_u$, the input is removed and the system is allowed to relax. The time to relax the imposed initial condition is then estimated within a threshold level $\delta = 5\%$.

<figure><img src="image_16.png" /><figcaption>Figure 13: Fractional order system relaxation time $L_u$ after step application for system with $\nu = 0.5(49)$</figcaption></figure>

Relaxation time of a first kind fractional order model with $\nu = 0.5$. Figure 14 illustrates the variations of the relaxation time with respect to the charging time, which clearly shows that fractional operators have non-local behavior: the relaxation time is conditioned by the charging time.

<figure><img src="image_17.png" /><figcaption>Figure 14: fractional order system relaxation vs charging time</figcaption></figure>

For this simple scenario, $\nu = 0.5$ and $\lambda = 1$, which gives the following short and long time relaxation times (from equations (45) and (48)): $t_{set}^{short} = 7.05s$ and $t_{set}^{long} = 127.32s$.

Relaxation time of a first kind fractional order model with $\nu = 0.4$. A second charge-discharge test is done in order to appreciate the impact of the fractional order $\nu$ on the system dynamics. Figure 15 illustrates the charge-discharge curve when the order of system (49) is set to $\nu = 0.4$.
<figure><img src="image_18.png" /><figcaption>Figure 15: fractional order system relaxation time $L_u$ vs charging time with $v = 0.4$</figcaption></figure>

For this case, the long time approximation gives $t_{set}^{long} = 660$s. However, this requires an extremely long charging period. For figure 15, the length of the input step has been taken up to 25000$s and the relaxation time is still far from the long time approximation.

Relaxation time of a first kind fractional order model with $v = 0.7$. A third test is shown in figure 16 where $v = 0.7$. As this system exhibits a significantly faster dynamic behavior, sampling time was modified to $T_s = 0.1$s. For this case, $t_{set}^{long} = 15.1$s, which is inferior to the steady state value shown in the figure (around 17.5$s). This result is not surprising, as the relaxation time starts approaching short times. Even though the approximation is not exact, it stays relatively close in magnitude. Further studies could be carried in order to determine the validity range of the long time approximation.

<figure><img src="image_19.png" /><figcaption>Figure 16: fractional order system relaxation time $L_u$ vs charging time with $v = 0.7$</figcaption></figure>

Effectively, for both systems, an insufficient charging time (because of a shorter step input) leads to a shorter relaxation time. However, as the step length increases its relaxation time gets progressively closer to that of the long time approximation. The worst case scenario would be a fully charged system

with a step applied at $t = -\infty$. It should be noted that such a scenario would be unrealistic and would require an infinite amount of energy. This may allow to safely take the long time approximation as a safe limit for memory selection. However, it requires a rough knowledge of the system behavior before its identification. This is already the case of state variable filters for recursive least squares identification technique (Padilla (2017)).

# 5. System Identification application

## 5.1. Long-Memory recursive prediction error method

In this section, a recursive algorithm for fractional order system identification will be used in order to estimate its parameters. The long-memory recursive prediction error method (LMRPEM) (Victor et al. (2022a)) has proven to provide consistent estimates with low variance and without bias when considering the whole data set from $t = 0$ to the current instant time $kT_s$. A brief recall of this method is presented below.

Classic prediction error method minimizes a quadratic error criterion:

$$J = \frac{1}{2} \sum_{k=1}^{t} \epsilon(k)^2 \quad (50)$$

where the error function is originally defined as:

$$\epsilon(t) = y^*(t) - \hat{y}(t), \quad (51)$$

where $y^*(t)$ stands for the noisy output and the estimated output $\hat{y}(t)$ is computed as:

$$\hat{y}(t) = G(p, \hat{\theta})u(t). \quad (52)$$

The long memory property of fractional order derivatives have led to the choice of an extended error. The error signal $\tilde{\epsilon}$

$$\tilde{\epsilon}(k) = [\epsilon(0), \epsilon(T_s), \epsilon(2T_s), \dots, \epsilon(kT_s)]^T \quad (53)$$

gathers the errors at all instants from $t = 0$ to the current time $t = kT_s$ and is increased by one data-point per iteration. Additionally, the gradient $\tilde{\psi}_{\rho}(k, \hat{\theta})$ is a matrix defined by:

$$\tilde{\psi}_{\rho}(kT_s) = - \begin{bmatrix} \frac{\partial \epsilon(0)}{\partial b_0} & \frac{\partial \epsilon(T_s)}{\partial b_0} & \frac{\partial \epsilon(2T_s)}{\partial b_0} & \dots & \frac{\partial \epsilon(kT_s)}{\partial b_0} \\ \frac{\partial \epsilon(0)}{\partial b_1} & \frac{\partial \epsilon(T_s)}{\partial b_1} & \frac{\partial \epsilon(2T_s)}{\partial b_1} & \dots & \frac{\partial \epsilon(kT_s)}{\partial b_1} \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ \frac{\partial \epsilon(0)}{\partial a_{m_A}} & \frac{\partial \epsilon(T_s)}{\partial a_{m_A}} & \frac{\partial \epsilon(2T_s)}{\partial a_{m_A}} & \dots & \frac{\partial \epsilon(kT_s)}{\partial a_{m_A}} \end{bmatrix} \quad (54)$$

where $k$ indicates the present iteration.

By introducing the extended measured output $\tilde{Y}^*(k)$ and the extended estimated output $\tilde{Y}(k)$ as:

$$\tilde{Y}^*(kT_s) = [y^*(0), y^*(T_s), y^*(2T_s), \dots, y^*(kT_s)]^T \quad (55)$$

$$\tilde{Y}(kT_s) = [\hat{y}(0), \hat{y}(T_s), \hat{y}(2T_s), \dots, \hat{y}(kT_s)]^T, \quad (56)$$

the long-memory recursive prediction-error method (LMRPEM) is proposed for fractional order systems (see Victor et al. (2022b)):

$$\begin{cases} \tilde{\epsilon}(k) = \tilde{Y}^*(kT_s) - \tilde{Y}^*(kT_s) \\ R(k) = R(k-1) + \gamma_{\rho} [\tilde{\psi}_{\rho}(kT_s) \tilde{\psi}_{\rho}^T(kT_s) - R(k-1)] \quad (57) \\ \hat{\theta}(k) = \hat{\theta}(k-1) + \gamma_{\rho} R^{-1}(k) \tilde{\psi}_{\rho}(kT_s) \tilde{\epsilon}(kT_s) \end{cases}$$
<figure><img src="image_20.png" /><figcaption>Figure 17: Input-output data with SNR = 15dB</figcaption></figure>

where $γ_ρ$ is a refining gain analogous to the step in a gradient descent.

It is necessary to analyze the property of this algorithm proposed in Victor et al. (2022a) to see if the truncated approximation of the derivatives can be used for online estimation of the parameters.

## 5.2. Simulation results

Suppose a first kind fractional order system:

$$G(s) = \frac{K}{1 + \tau s^{\nu}} \qquad (58)$$

with the following true values: $K = 1$, $τ = 1$ and $ν = 0.5$. The input is a pseudo-random binary sequence oscillating between -5 and 5 and containing 2000 points with a sampling time of $T_s = 1.00s$. All simulations are initialized with $\hat{K}_0 = 2$ and $\hat{\tau}_0 = 4$, and $ν = 0.5$ is assumed known. The forgetting factor $γ_ρ = 0.01$ and the signal to noise ratio is $SNR = 15 dB$. Figure 17 shows the input-output data for the simulations.

### Case 1: without truncation

Parameter estimation as well as calculation time are presented in figures 18 and 19, respectively. As expected, an accurate estimation of the parameters is obtained. However, as expected, the calculation time dramatically increases with each iteration. The sampling time is set to $T_s = 1s$ and after 2000s, the calculation time has almost reached one third of the sampling time. Consequently, using the whole data set from $t = 0$ to the current time $t = kT_s$ for more longer scenario will lead to an unimplementable algorithm.

<figure><img src="image_21.png" /><figcaption>Figure 18: LMRPEM system identification without truncation</figcaption></figure>

<figure><img src="image_22.png" /><figcaption>Figure 19: Computation time vs iteration number without truncation</figcaption></figure>

### Case 2: truncation window equal to the relaxation time

Knowing the true system dynamics, if the truncation time window to approximate the fractional derivative is set to the system relaxation time, the parameter estimation as well as the calculation time are presented in figures 20 and 21, respectively. Truncation length is chosen to be $L = 128s$ based on long-time approximation. For this case, the windowing truncation leads to slight oscillations in the parameter estimation. However, it can be seen that they are still correctly estimated. Calculation time has been dramatically decreased and is found to have a mean value of $t_{calc} = 14.7ms$ per iteration (see figure 21). This is far beyond the sampling time limit, which enables real-time implementation. Parameter fluctuations may be due to rather short memory length $L$, as the LMRPEM method uses an extended error and the whole relaxation error is non-negligible.
<figure><img src="image_23.png" /><figcaption>Figure 20: LMRPEM system identification with a truncation window $L = 128$ s</figcaption></figure>

<figure><img src="image_24.png" /><figcaption>Figure 21: Computation time vs iteration number with a truncation window $L = 128$ s</figcaption></figure>

### Case 3: truncation window longer than the relaxation time
Given the true system dynamics, if the truncation time window to approximate the fractional derivative is now set longer than the system relaxation time, the parameter estimation as well as calculation time are presented in figures 22 and 23, respectively. The truncation window length is chosen to be $L = 256$ s (twice the relaxation time). In this case, oscillations are indeed significantly reduced. Even though the performance is not identical to that of the non-truncated case, it remains very close. Moreover, the calculation cost is not significantly increased with respect to the previous case. The mean calculation time is around $t_{calc} = 17.5$ ms, which remains far lower than the sampling time (see figure 23). The calculation cost curve shows a slight increase in the calculation time during the first iterations and then the truncation window imposes a limit to this increasing of the calculation time.

## 6. Conclusions and perspectives

The non-locality property or long-memory property of fractional derivatives constitutes both a main advantage and drawback in mathematical modelling. It allows to capture long-memory behavior, but may impose implementation limits in

<figure><img src="image_25.png" /><figcaption>Figure 22: LMRPEM system identification with a truncation window $L = 256$ s</figcaption></figure>

<figure><img src="image_26.png" /><figcaption>Figure 23: Computation time vs iteration number with a truncation window $L = 256$ s</figcaption></figure>

real-time scenarios. Podlubny's Short Memory Principle provides a mathematical proof of a truncation length $L$ to be used in order to get an approximated fractional derivative without the whole past of the function. However, even though it effectively provides a worst-case scenario, Podlubny's limit does not take into account the signal spectrum. A simple cosine function was analyzed by starting with Riemann-Liouville's definition, leading to an unsolvable integral for the analytic error expression. As a consequence, numerical simulations using Grünwald-Letnikov definition have been carried out. It is then seen that a fixed memory length $L$ will lead to smaller errors as the signal frequency $\omega$ increases. Even though this may suggest that the error reduction is due to an increased number of periods of the signal taken into account, this was shown to be not accurate enough. A fixed number of periods of a periodic function past increases the truncation error as the signal frequency increases. On the other hand, the influence of the fractional order $\nu$ was also studied. As the differentiation order approaches the integer cases, the truncation errors tend to zero. This is a logical result as integer order derivatives are local operators which do not require an extensive knowledge of the past.

The second part of this paper considered the initialization
problem when dealing with real-time (online) system identifi-
cation. Using only a truncated part of the system input *u* in or-
der to estimate its output *ŷ* implies an initialization error at the
chosen instant to start the convolution product. If the system
is stable, error due to a wrong initialization and neglecting the
distant past will progressively disappear. Transient response for
a simple fractional order system as well as approximations for
Mittag-Leffler's function were analyzed. The influence of the
step length before releasing a fractional order system was also
studied through a simulation example. Mittag-Leffler's long-
time approximation proved to be a safe limit to take into ac-
count for the choice of the truncation window *L*.

Finally, the concept of system relaxation was applied to sys-
tem identification. The LMRPEM algorithm was used with
different truncation windows in order to evaluate its real-time
feasability. Effectively, the most accurate estimation is obtained
for the non-truncated case. The truncation with a memory *L*
similar to the system relaxation time proved to correctly esti-
mate parameters, but fluctuations are present. Increasing the
memory length reduces such fluctuations and the calculation
time is still satisfying. This implies that a real-life case implies
a rough approximation of the system relaxation time: some
level of previous empirical knowledge of the plant may be re-
quired to correctly perform the system identification.

There are two main perspectives for this study. The first one
is to further analyze the link between parameter fluctuation and
calculation time. Depending on the scenario, optimizing the
truncation length *L* could lead to satisfying system identifica-
tion without exceeding the sampling time. On the other hand, a
more global fractional order system identification could be an-
alyzed. If the fractional differentiation orders need to be identi-
fied, it may have a severe influence on the choice of truncation
window. The Influence of the truncation window over differ-
entiation order estimation could be further analyzed. Adaptive
methods could also be considered. An improved simulation
scenario may also be proposed by optimizing the experience
and applying a more optimal type of input, as in (cite).

References

Battaglia, J., Cois, O., Puigsegur, L., Oustaloup, A., 2001. Solving an inverse heat conduction problem using a non-integer identified model. International Journal of Heat and Mass Transfer 44, 2671–2680. doi:10.1016/S0017-9310(00)00310-0.

De Wit, C., 1986. Recursive estimation of the continuous-time process param-
eters, in: 1986 25th IEEE Conference on Decision and Control, pp. 2016-
2020. doi:10.1109/CDC.1986.267390.

Duhe, J., Victor, S., Melchior, P., Abdelmoumen, Y., Roubertie, F., 2021. Recursive system identification for coefficient estimation of continuous-time fractional order systems, in: 19th IFAC Symposium on System Identification (SYSID 2021), Padova, Italy. pp. 114–119. doi:10.1016/j.ifacol.2021.08.344.

Duhé, J., Victor, S., Melchior, P., Abdelmoumen, Y., Roubertie, F., 2022.
Modeling thermal systems with fractional models: human bronchus
application. Nonlinear Dynamics 108, 579–595. doi:10.1007/s11071-022-07239-3.

Garnier, H., Wang, L., 2008. Identification of continuous-time models from
sampled data. Springer-Verlag.

Garrappa, R., Kaslik, E., Popolizio, M., 2019. Evaluation of fractional integrals and derivatives of elementary functions: Overview and tutorial. Mathematics 7. doi:10.3390/math7050407.

Grünwald, A., 1867. Über begrenzte Derivationen und deren Anwendung.
Zeitschrift für Mathematik und Physik, 441–480.

Ionescu, C., Copot, D., De Keyser, R., 2014. A fractional order impedance model to capture the structural changes in lungs. IFAC Proceedings Volumes 47, 5363–5368. doi:10.3182/20140824-6-ZA-1003.01124.

Ionescu, C., De Keyser, R., Sabatier, J., Oustaloup, A., Lebron, F., 2011. Low frequency constant-phase behavior in the respiratory impedance. Biomedical Signal Processing and Control 6, 197–208. doi:10.1016/j.bspc.2010.10.005. special Issue: The Advance of Signal Processing for Bioelectronics.

Kilbas, A.A., Srivastava, H.M., Trujillo, J.J., 2006. Theory and Applications of Fractional Differential Equations, Volume 204 (North-Holland Mathematics Studies). Elsevier Science Inc., New York, NY, USA.

Le Lay, L., 1998. Identification temporelle et fréquentielle par modèle non entier. Ph.D. thesis. Université Bordeaux 1. France.

Lenka, B., 2019. Fractional comparison method and asymptotic stability results for multivariable fractional order systems. Communications in Nonlinear Science and Numerical Simulation 69, 398–415. doi:10.1016/j.cnsns.2018.09.016.

Letnikov, A., 1868. Theory of differentiation of arbitrary order (russian). Matematiceskij Sbornik (Moscow) 3, 1–68.

Ljung, L., 1999. System identification – Theory for the user. 2 ed., Prentice-Hall, Upper Saddle River, N.J., USA.

Magin, R., 2004. Fractional calculus in bioengineering. Critical reviews in biomedical engineering 32, 1–104. doi:10.1615/critrevbiomedeng. v32.i1.10.

Magin, R., Ovadia, M., 2006. Modeling the cardiac tissue electrode interface using fractional calculus. 2nd IFAC Workshop on Fractional Differentiation and its Applications 39, 302–307. doi:10.3182/20060719-3-PT-4902. 00056.

Mainardi, F., Mura, A., Pagnini, G., Gorenflo, R., 2006. Fractional relax-
ation and time-fractional diffusion of distributed order, in: 2nd IFAC
Workshop on Fractional Differentiation and its Applications, pp. 1–21.
URL: https://www.sciencedirect.com/science/article/pii/
S147466701536465X, doi:10.3182/20060719-3-PT-4902.00002.

Malti, R., Mayoufi, A., Victor, S., 2022. Experiment design for elementary fractional models. Communications in Nonlinear Science and Numerical Simulation doi:10.1016/j.cnsns.2022.106337.

Malti, R., Moreau, X., Khemane, F., Oustaloup, A., 2011. Stability and res-
onance conditions of elementary fractional transfer functions. Automatica
47, 2462–2467. doi:10.1016/j.automatica.2011.08.029.

Malti, R., Sabatier, J., Akçay, H., 2009. Thermal modeling and identifica-
tion of an aluminium rod using fractional calculus, in: 15th IFAC Sympo-
sium on System Identification (SYSID’2009), St Malo, France. pp. 958–
963. doi:10.3182/20090706-3-FR-2004.00159.

Malti, R., Victor, S., Oustaloup, A., 2008a. Advances in system identification using fractional models. Journal of Computational and Nonlinear Dynamics 3, 021401.1-021401.7. doi:10.1115/1.2833910.

Malti, R., Victor, S., Oustaloup, A., Garnier, H., 2008b. An optimal instrumen-
tal variable method for continuous-time fractional model identification, in:
The 17th IFAC World Congress (IFAC’08), Seoul, Korea. pp. 14379–14384.
doi:10.3182/20080706-5-KR-1001.02436.

Matignon, D., d'Andréa-Novel, B., 1996. Some results on controllability and observability of finite-dimensional fractional differential systems, in: IMACS, IEEE-SMC, Lille, France. pp. 952–956.

Melchior, P., Victor, S., Pellet, M., Petit, J., Cabelguen, J., Oustaloup, A., 2019.
Handbook of Fractional Calculus with Applications. Applications in Engi-
neering, Life and Social Sciences, Part A. De Gruyter. volume 7. chapter
Skeletal muscle modeling by fractional multi-models: analysis of length ef-
fect. doi:10.1515/9783110571905-006.

Moze, M., Sabatier, J., 2005. LMI tools for stability analysis of frac-
tional systems, in: 20th ASME International Design Engineering Technical
Conferences and Computers and Information in Engineering Conference,
IDETC/CIE’05, Long Beach, CA. pp. 1–9.

Oustaloup, A., 1995. La dérivation non-entière : théorie, synthèse et applica-
tions. Hermès, Paris.

Oustaloup, A., 2014. Diversity and Non-integer Differentiation for System Dynamics. Wiley-ISTE.

Oustaloup, A., Le Lay, L., Mathieu, B., 1996. Identification of non integer order systems in the time domain, in: IEEE-CESA'96, SMC IMACS Multi-conference. SMC IMACS Multiconference - Computational Engineering in
Systems Application - Symposium on Control, Optimisation and Supervision, Lille, FRANCE. pp. 843-847.

Oustaloup, A., Moreau, X., Nouillant, M., 1993. From the second generation crone control to the crone suspension, in: Proceedings of IEEE Systems Man and Cybernetics Conference - SMC, Le Touquet, France. pp. 143-148.
doi:10.1109/ICSMC.1993.384863.

Padilla, A., 2017. Recursive identification of continuous-time systems with time-varying parameters. Theses. Université de Lorraine.

Podlubny, I., 1999. Fractional Differential Equations: An Introduction to Fractional Derivatives, Fractional Differential Equations, to Methods of Their Solution and Some of Their Applications. Academic Press, San Diego.

Rivero, M., Rogosin, S., Tenreiro Machado, J., Trujillo, J., 2013. Stability of fractional order systems. Mathematical Problems in Engineering 2013, 356215. doi:10.1155/2013/356215.

Samko, S., Kilbas, A., Marichev, O., 1993. Fractional integrals and derivatives: theory and applications. Gordon and Breach Science.

Shah, P., Agashe, S., 2016. Review of fractional pid controller. Mechatronics 38, 29-41. doi:10.1016/j.mechatronics.2016.06.005.

Sommacal, L., Melchior, P., Cabelguen, J., Oustaloup, A., Ijspeert, A., 2006. Fractional multimodels of the gastrocnemius frog muscle. 2nd IFAC Workshop on Fractional Differentiation and its Applications.

Trigeassou, J.C., Maamri, N., 2011. Initial conditions and initialization of linear fractional differential equations. Signal Processing 91, 427-436. doi:10.1016/j.sigpro.2010.03.010.

Trigeassou, J.C., Maamri, N., 2020. Analysis, Modeling and Stability of Fractional Order Differential Systems 2: The Infinite State Approach. Wiley-ISTE.

Victor, S., Duhe, J., Melchior, P., Abdelmoumen, Y., Roubertie, F., 2022a. Long memory recursive prediction error method for identification of continuous-time fractional models. Nonlinear Dynamics 110, 635-648. doi:10.1007/s11071-022-07628-8.

Victor, S., Malti, R., Garnier, H., Oustaloup, A., 2013. Parameter and differentiation order estimation in fractional models. Automatica 49, 926-935. doi:10.1016/j.automatica.2013.01.026.

Victor, S., Mayoufi, A., Malti, R., Chetoui, M., Aoun, M., 2022b. System identification of MISO fractional systems: parameter and differentiation order estimation. Automatica 141. doi:10.1016/j.automatica.2022.110268.

Victor, S., Melchior, P., Malti, R., Oustaloup, A., 2016. Robust motion planning for a heat rod process. Journal of Nonlinear Dynamics 86, 1271-1283. doi:10.1007/s11071-016-2963-2.

Victor, S., Melchior, P., Pellet, M., Oustaloup, A., 2020. Lung thermal transfer system identification with fractional models. IEEE Transactions on Control Systems Technology 28, 172-182. doi:10.1109/TCST.2018.2877606.

Victor, S., Ruiz, K., Melchior, P., Chaumette, S., 2022c. Dynamical repulsive fractional potential fields in 3d environment. Fractional Calculus & Applied Analysis, 321-345. doi:10.1007/s13540-022-00015-5.
# Highlights

## Fractional derivative truncation approximation for real-time applications

Jean-François Duhé, Stéphane Victor, Pierre Melchior, Youssef Abdelmounen, François Roubertie

* Short Memory Principle

* Fractional Relaxation

* System identification