November 25, 2015

International Journal of Control

FO˙tool˙eval

To appear in the International Journal of Control
Vol. 00, No. 00, Month 20XX, 1–25

SURVEY PAPER

arXiv:1511.07521v1 [cs.SY] 24 Nov 2015

A review and evaluation of numerical tools for fractional calculus and
fractional order controls
Zhuo Lia Lu Liuab Sina Dehghana YangQuan Chena∗ and Dingyü Xueb
a

School of Engineering, University of California, Merced. 5200 N. Lake Road, Merced, CA
95343, USA b College of Information Science and Technology, Northeastern University,
Shenyang 110819, China
(June 2015)

∗ Corresponding author. Mechatronics, Embedded Systems and Automation Laboratory (MESA Lab) http://mechatronics.
ucmerced.edu. Email:ychen53@ucmerced.edu

November 25, 2015

International Journal of Control

FO˙tool˙eval

Abstract
In recent years, as fractional calculus becomes more and more broadly used in
research across different academic disciplines, there are increasing demands for
the numerical tools for the computation of fractional integration/differentiation,
and the simulation of fractional order systems. Time to time, being asked about
which tool is suitable for a specific application, the authors decide to carry out this
survey to present recapitulative information of the available tools in the literature,
in hope of benefiting researchers with different academic backgrounds. With this
motivation, the present article collects the scattered tools into a dashboard view,
briefly introduces their usage and algorithms, evaluates the accuracy, compares the
performance, and provides informative comments for selection.
Keywords:

Fractional calculus, fractional order controls, numerical tools.

2

November 25, 2015

1.

International Journal of Control

FO˙tool˙eval

Introduction

The fractional calculus (FC) got birth 300 years ago, and the research on fractional calculus experienced its boom in the past decades [1, 2, 3]. Besides the fundamental mathematical study, more
and more researchers from different academic disciplines begin to utilize it in a varieties of subjectassociated research, such as in biology and biomedical [4, 5], sociology [6, 7], economics [8, 9], and
control engineering [10, 11, 12], etc. Along with the rapid development of theoretical study, the
numerical methods and practical implementation also made considerable progress [13, 14, 15].
Sharp tools are prerequisite to a successful job. In this paper, an extensive collection of Matlab based tools are exhibited for the numerical computation of fractional order (FO) integration/differentiation, as well as some toolboxes for engineering applications, with an emphasis on
fractional order controls. A comprehensive table, table 1, is created to list the recapitulative information of these tools in a dashboard view. Brief description and basic evaluation of these numerical
algorithms are presented, in terms of usage, accuracy, unique features, advantages and drawbacks.
Through such efforts, it is hoped that an informative guidance is provided to the readers when they
face to the problem of selecting a numerical tool for a specific application. Thanks to the authors
of these tools. It is these pioneers who bring great convenience for the practical use of FC and FO
control. While a text descriptive survey on some of the tools under discussion can be found in book
[16], and 28 alternatives for the time-domain implementation of FO derivatives are documented in
[17], the present paper addresses more quantitative comparison and practical usage.
The rest of the paper are organized as follows: section 2 reviews 20 selected numerical tools
through brief description; section 3 evaluates and compares the quantitative performance of the
tools in three categories; section 4 gives comments for tool selection based on empirical use.
2.

Review and description

This article mainly covers the tools for fundamental fractional calculus, such as the numerical
computation of fractional integration/differentiation of a function or a signal, the Laplace transform of fractional differential equations [18], etc. Since automatic control is one of the engineering
disciplines that got the earliest exposure to fractional calculus [19, 20, 21, 22], the tools for the
application of fractional order controls are given more focus, associated with the authors’ expertise.
2.1

@fotf

@fotf (fractional order transfer function) is a control toolbox for fractional order systems developed
by Xue et al. Most of the functions inside are extended from the Matlab built-in functions. In
[23], the code and usage of the @fotf toolbox are described in very detail. It uses the overload
programming technique to enable the related methods of the Matlab built-in functions to deal
with FO models. The transfer function objects generated from it can be interactive with those
generated from the Matlab transfer function class. Yet, the overloading of associated functions such
as impulse(), step(), etc, lost the plotting functionality. As a work around, users can simply define
a time vector as the second input to these functions. fotf toolbox supports time delay in the TF,
e.g. fotf(a,na,b,nb,delay). It does not directly support transfer function matrix, hence, MIMO
systems cannot be simulated directly. However, since it provides Simulink block encapsulation of
the involved function fotf(), multiple input/output relationship can be established by manually
adding loop interactions in Simulink block diagrams. Therefore, the remark “could” is put in the
“MIMO” column in table 1, (where the ‘Delay’ column denotes if the script/toolbox is able to
handle time delay in the FO model; and the ‘MIMO’ column denotes if the script/toolbox is able
to handle MIMO FO models.).
A small drawback with @fotf is that the sampling time has relatively big impact on the accuracy,
3

November 25, 2015

International Journal of Control

FO˙tool˙eval

which has been remarked in the validation comments in [23]. Encouragingly, an update is upcoming
according to the author.
2.2

Ninteger

Ninteger, non-integer control toolbox for Matlab, is a toolbox intended to help with developing
fractional order controllers and assessing their performance, [24]. It uses integer order transfer
functions to approximate the fractional order integrator/differentiator, C(s) = ksν , ν ∈ R. It
offers three frequency domain approximation methods,
(1) The CRONE methods, which uses a recursive distribution,
C(s) = k′

N
Y
1 + s/ωzn
;
1 + s/ωpn

n=1

(2) The Carlson’s method that solves C α (s) using Newton’s iterative method,
Cn (s) = Cn−1 (s)

α (s) + (α + 1)g(s)
(α − 1)Cn−1
α (s) + (α − 1)g(s) ;
(α + 1)Cn−1

(3) The Matsuda’s methods, that approximates C with a gain known at several frequencies.
(s − ωk−1 )/dk (ωk )]+∞
k−1 ,
ω − ωk
.
d0 (ω) = |C(jω)|, dk+1 (ω) =
dk (ω) − dk (ωk )
C(s) = [d0 (ω0 );

It also provides Simulink block encapsulation of the involved functions, such as ‘nid’ and ‘nipid’
blocks. Moreover, it offers a user-friendly GUI for fractional order PID controller design.
There is a problem with ninteger toolbox in Matlab version 2013a or later. Without additional
editing, it has conflicts with some built-in functions due to the overload editing of the Matlab
built-in function isinteger(). For example, calling the mean() function will prompt an error.
2.3

ooCroneToolbox

The CRONE Toolbox, developed since the nineties by the CRONE team, is a Matlab and Simulink
toolbox dedicated to applications of non integer derivatives in engineering and science [25]. It
evolved from the original script version to the current object-oriented version. A good feature of
the Crone toolbox is that some of the methods are implemented for MIMO fractional transfer functions. For example, executing sysMIMO=[sys,sys;sys2,sys2] generates a two-input-two-output
TF matrix. Many simulation results in the literature are obtained using the CRONE toolbox such
as the design of centralized CRONE controller with the combination of the MIMO-QFT approach
in [26]. Several other toolboxes are inspired by CRONE, e.g. ninteger and FOMCON. A drawback
of the CRONE toolbox is that time delay cannot be incorporated into the generated FO TF. Manually multiplying the delay to the frac tf object does not work either because the exp() operation
is not overloaded by frac tf class. CRONE is a toolbox much more powerful than merely simulating fractional order systems. In spite of this basic functionality, it is also capable of fractional
order system identification and robust control analysis and design.

4

November 25, 2015

International Journal of Control

2.4

FO˙tool˙eval

FOMCON

The FOMCON (Fractional-Order Modeling and Control) toolbox is developed by Tepljakov et.
al, [27]. Its kernel utilizes the algorithms in FOTF, Ninteger and Crone. It encapsulates some
of the major functionalities of those three toolboxes, and builds a GUI shell on top, aiming at
extending classical control schemes for FO controller designs. The relation of FOMCON with the
three toolboxes is shown in figure 1. Some notable changes/patches to the original FOTF are:
• newfotf() uses the string parser to enable users to input TF as a string;
• tf2ss() is overloaded and foss() is added, which makes the conversion between an FO TF
object and an FO state space object possible. The CRONE toolbox is also able to do the
task, yet the script is encrypted in Matlab P code format.
2.5

M-L functions

M-L functions, as the name implies, are Matlab functions developed for numerically computing the
Mittag-Leffler function (definitions can be found in [4] etc). There are several versions of code by
different authors available in the literature. Five of them are listed in table 1, where
(1) mlf(α, β, x, p) is for the calculation of the 2-parameter M-L function in the form of Eα,β (x)
with the precision of p for each element in x;
(2) ml func([α, β, γ, q], z, n, ε0 ) is capable of computing the M-L function with either 1, 2, 3, or
4 parameters, and the script is available in the books [28] or [10]. It uses the fast truncation
algorithm to improve the efficiency, and embeds the mlf() in the file such that when the fast
truncation algorithm is not convergent, solution is guaranteed by trading off some efficiency;
(3) ml fun(α, β, x, n, ε0 ) (α > 0, β > 0) is also for 2-parameter M-L function with error tolerance
of ε0 , which is implemented using C-MEX .dll (dynamic-link library) technique and can be
used in Simulink through s-functions;
(4) gml fun(α, β, γ, x, ε0 ) calculates the generalized M-L function with 3 parameters in the form
γ
of Eα,β
(x), [29];
(5) ml(x, α, β, γ) can calculate the M-L function with either 1, 2, ro 3 parameters.
Alternatively, the generalized hypergeometric function [pfq]=genHyper (a,b,z, lnpfq, ix,
nsigfig) in [30], or [y,tt,nterms]=pfq (a,b,z,d) in [31] can also achieve the numerical computation of the generalized M-L functions under certain conditions. For more details, refer to [32].
2.6

NILT

The inversion of Laplace transform is fundamentally important in the applications of Laplace
transform method. It can be carried out with one of the following three approaches: 1). analytical
solution using definition and basic properties; 2). Laplace transform tables; and 3). numerical
computation. While analytical solutions are usually too hard to be obtained, and tables do not
cover arbitrary cases, the numerical computation becomes an inevitable way. Among the numerous
algorithms for numerical inversion of Laplace transform (NILT), NILT in [33, 34] and the “improved
NILT” in [35, 36, 37] have relatively bigger literature exposure. Lubomir’s NILT method applies
the fast Fourier Transformation (FFT) and the ε-algorithm to speed up the convergence of infinite
complex Fourier series. A very detailed description and performance evaluation of these methods
is available in [38]. Hence, repetitive comparison among different NILTs are not presented here.
Focus is mainly put on the comparison between NILT and other numerical methods.
A good feature of the two NILT code is that both support the direct input of time delay in
the form of exp(-Ls). Yet, INVLAP() gives some glitch at the end of the delay, for example,
[x,y]=INVLAP(’1/(s* (s^0.5+1)) * exp(-s)’, 0.01,10,1000). There is a tricky part need
5

November 25, 2015

International Journal of Control

FO˙tool˙eval

to be noted in evaluating the computational error of NILT. If the same initial, terminating and
sampling time (t0 , tf and Ts ) for other tools are used in the script, the NILT actually computes
one point less than the other tools which use regularly spaced time vector. That is because: let
0
represent the amount of points computed by NILT, then, the time interval is actually
M = tfT−t
s
t
−t
f
0
′
Ts = M −1 due to the script t=linspace(0,tm,M). Whereas the conventional assignment of time
vector (t=t0:Ts:tf) generates M+1 points. In order to compute the same amount of points aligned
to the time stamps used for baseline analytical solution, the time vector for analytical computation
needs to be adjusted so as to adapt to that used by NILT. This means to let analytical computation
use the time vector generated by NILT, which can be achieved by either 1). t=0:M*Ts/(M-1):M*Ts,
or 2). t=linspace(0,tm,M). This cannot be done the other way around, i.e. replaced by t =
0:Ts:M*Ts-Ts nor t=linspace(0,tm-Ts,M). Otherwise, cumulated computation error will cause
inaccuracy of the final simulation result. Alternatively, if tf is not a concern, user can assign one
point less to M in the NILT script while keeping Ts unchanged. Thus, NILT generates the same time
stamps except a tf shortened by one sampling period. The difference in dealing with time vectors
can be easily visualized if longer sampling time is assigned. An example of the resulting computation
error is demonstrated in figure 2. Similar time stamp assignment issue exists in INVLAP(). In
addition, the initial time stamp is not allowed to be 0 due to the constraint in the INVLAP() script.
2.7

dfod

DFOD (Digital Fractional Order Differentiator/integrator) is a set of Maltab functions writhen by
Petráš, for the approximation of fractional order differentiators and integrators. There are three
versions of dfod:
(1) dfod1() is the infinite impulse response (IIR) type based on continued fraction expansion
(CEF), shown in equation (1), of weighted operator with the mixed scheme of the trapezoidal
(Tustin) rule and the backward difference (Euler) rule, [39];
Z{D α x(t)} = CF E{(

1 Pp (z −1 )
1 − z −1 α
) }X(z) ≈ ( )α
X(z).
T
T Qq (z −1 )

(1)

(2) dfod2() is the finite impulse response (FIR) type based on power series expansion (PSE),
shown in equation (2), of the backward difference (Euler) rule, [40];
D ∓α (z) =

T ∓α
T ∓α
 
≈
Qq (z −1 )
±α −j
(−1)j
z
j
j=0

1
= ∞
P
(1 − z −1 )±α

(2)

(3) dfod3() is a new IIR type based on power series expansion of the trapezoidal (Tustin) rule,
[41].
Euler : sα ≈



1 − z−1
T

α

Tustin : sα ≈

,



2 1 − z−1
T 1 + z−1

α

.

(3)

There are other FO algorithms based on IIR, such as newfod() by Chen, [42].
Regarding discretization, besides the aforementioned methods used in the various tools, other
methods exist such as the Prony’s technique, direct discretization, the binomial expansion of the
backward difference, etc, [43].

6

November 25, 2015

International Journal of Control

2.8

FO˙tool˙eval

IRID

The impulse response invariant discretization (IRID) is a family of functions designed by Chen,
Li, Sheng et al. [44, 45], for different approximation purposes based on the algorithm as its name
implies. It includes the following members:
(1) irid fod() is designed to compute a discrete-time finite dimensional (z) transfer function to
approximate a continuous irrational transfer function sα where ‘s’ is the Laplace transform
variable and −1 < α < 1. It has been tested that the algorithm still works for α > 1 and
α < −1, by removing the input checking statement.
Rb
(2) irid doi() is for the approximation of distributed order integrator a s1α dα, where ‘a’ and
‘b’ are arbitrary real numbers in the range of (0.5 ,1), and a < b.
(3) irid dolp() is for the approximation of a continuous-time fractional order low-pass filter in
the form of 1/(τ s + 1)α
(4) irid fsof() is for the approximation of fractional second order filter in the form of 1/(s2 +
as + b)α where 0 < α < 1.
(5) BICO irid() is for thep
approximation of BICO (Bode’s Ideal Cut-Off) transfer function in
the form of 1/(s/w0 + (s/w0 )2 + 1)α , where α > 0.
2.9

ora foc

ora foc() is for the approximation of fractional order differentiators,
Recursive-Approximation method described in [47].
2.10

1
sα ,

[46], using the Oustaloup-

fderiv

fderiv() calculates the fractional derivative of order α for the given function r(t) using the
Grünwald-Letnikov (G-L) definition, [48]. The input of the given function is represented by a
vector of signal values. There is an improved implementation of this function, fgl deriv(), by
Jonathan, which uses vectorization for faster computation with Matlab, [49].
2.11

glfdiff

glfdiff(y,t,α) (G-L finite diff) is a Matlab function written by Xue et al. [50] for calculating
the αth derivative of a given function, whose inputs y, t are the signal and time vectors. It is based
on the forward finite difference approximation of the G-L definition,
α
a Dt f (t)

(t−1)/h
1 X (α)
ωj f (t − jh),
≈ α
h

(4)

j=0

where the binomial coefficients are recursively calculated, [50]:
(α)
ω0

2.12

= 1,

(α)
ωj

=



α+1
1−
j



(α)

ω(j−1) ,

j = 1, 2, . . .

(5)

Fractional differentiation and integration

Many of the above functions approximate the fractional order integral or derivative operator.
This Matlab function calculates the αth order derivative or integral of a function, defined in a
7

November 25, 2015

International Journal of Control

FO˙tool˙eval

given range through Fourier series expansion. The necessary integrations are performed with the
Gauss-Legendre quadrature rule, [51]. Three examples are offered in this package, namely FO
differ/integral of identity, cubic polynomial and tabular functions, respectively. The main call
function is fourier diffint().
2.13

FIT

FIT is the Fractional Integration Toolbox developed by Santamaria Laboratory at the University
of Texas at San Antonio, [52]. It is for the numerical computation of fractional integration and
differentiation of the Riemann-Liouville (R-L) type, and is designed for large data size, which allows
parallel computing of multiple fractional integration/differentiation on GPUs (graphical processing
units). The extrapolation and interpolation algorithms used by this toolbox are implemented in
C++ and are integrated with Matlab via MEX mechanism. Detailed explanation can be found in
[53].
2.14

DFOC

DFOC, written by Petráš, is a digital version of the Fractional-Order PID Controller of the form:
C(s) = K + Ti

1
+ Td sd .
sm

(6)

It provides a transfer function of the FO PID controller for given parameters, [54].
2.15

FOPID

The fractional order PID (FOPID) controller toolbox, presented by Lachhab et al., is for the design
of robust fractional order P I α D β controllers, [55]. The tuning rules for the parameters follow those
promoted in [56] and [57]. Thus, the fractional order PID tuning is converted to a 5-parameter
optimization problem. This toolbox utilize the “non-smooth” H∞ synthesis in [58] to perform the
minimization. For now, there is not a publicly available source for download.
2.16

Sysquake FO PID

In [59], Pisoni et al. presented an interactive tool for fractional order PID controllers developed on
the Sysquake software environment, which is a similar effort with that for integer order PIDs done
by Åström et al. in [60]. Sysquake is a numerical computing environment based on a programming
language mostly-compatible with Matlab. However, the interactive tool for FO PID runs in the
Sysquake environment instead of Matlab. Hence, it is not reviewed in detail here.
2.17

FOCP

In [61], Tricaud and Chen et al. formulated the Fractional Optimal Control Problems (FOCP) into
the integer order format by using a rational approximation of the fractional derivative obtained
from the singular value decomposition (SVD) of the Hankel matrix of the impulse response. Then,
RIOTS 95 [62, 63] is used to perform the optimization. The scheme is potentially able to solve any
type of FOCPs and is implemented in Matlab for public accessibility, [64]. It supports MIMO FO
optimal control, but does not handle time delay due to the limitation of RIOTS.

8

November 25, 2015

2.18

International Journal of Control

FO˙tool˙eval

FSST

FSST is a simulation toolkit in Matlab/Simulink for the fractional order discrete state-space system
education. The toolkit consists of a set of C-MEX s-functions which are encapsulated in Simulink
blocks. Several typical fractional order system simulation examples are provided as shown in figure
3, such as the fractional order state-space model and the fractional Kalman filter (FKF), [65].
The version 1.7 is available for free download at [66]. Two of the superior strengths of FSST are:
1.) it can directly simulate MIMO systems since it is a Simulink block kit handling state space
representations; 2.) it is able to incorporate the initial conditions into the dynamic equations to be
simulated, which is a unique feature among all the aforementioned tools. The drawback of FSST
is that the step size has large impact on the simulation results, even larger than the impact by
“cilcular” buffer size. A sample illustration is plotted in figure 4.
2.19

Fractional variable orders

All the above tools/toolboxes (except irid doi()) deal with constant fractional orders. Yet, there
exists a type of differentiations that have fractional variable orders (FVO). The definitions in the
G-L format are given as follows, [67]:
Theorem 1 (The 1st type FVO):
α(t)
f (t)
0 Dt

= lim

h→0

1
hα(t)

n
X
r=0



α(t)
(−1)
f (t − rh).
r
r

(7)

The 2nd and 3rd types can be found in the same reference.
Regarding the fractional variable order differentiation, there are dedicated tools. Podlubny et al.
offers a matrix approach that unifies the numerical differentiation of integer order and the n-fold
integration, using the so-called triangular strip matrices, [68]. It is available for download at [69]
and can be applied on the solution to FODEs and FPDEs.
Sierociuk et al. provides a C-MEX s-function based Simulink toolkit, “fvoderiv”, for this purpose,
[70]. It supports real-time-workshop.
The toolbox “vod” created by Valério et al. calculates variable fractional or complex order derivatives. R-L, Caputo and G-L definitions are provided; the three types of definitions in [67] are all
considered. Fuzzy supervised implementations in Simulink are also provided, [71].
2.20

FO root locus

Three Matlab based scripts for plotting root locus of fractional order TFs are available. Two early
works are frlocus() in [72], and the code attached in the paper [73] by Machado et al.. The
other is forlocus() developed by the author which is listed in the last row in table 1 and can be
downloaded from [74]. Besides, the newest version of @fotf toolbox also features the root locus plot
of FO systems. Figure 5 shows a demonstrating plot of the root locus of the following fractional
order transfer function,
G(s) =

1.2s1.3 + 1
,
0.8s2.6 + 0.6s1.3 + 1

(8)

where figure 5a shows the plot on the Laplace s-plane and figure 5b shows the plot on the w = s1.3
plane. A closer view of the 2nd quadrant in figure 5b tells that the root locus in this example has
two branches on the first layer of the Riemann sheet [75, 76]. One starts from the pole marked in
9

November 25, 2015

International Journal of Control

FO˙tool˙eval

green, and the other is from the next Riemann sheet. As the system gain varies, they aggregate at
(−1.25 + 1.1i) and then bifurcate. One approaches to the open loop zero marked in red and the
other goes to infinity.
2.21

Other tools

Text description of a few tools listed above can also be found in [77]. There are other fractional
calculus related tools or Matlab scripts available for specific applications, such as the fractional
Fourier transform (FrFT) [78, 79], closed-form solutions to linear fractional order differential equations, fode sol() [10], the M-L random number generator mlrnd() [80], digital fractional order
Savitzky-Golay differentiator [81], and the functions for simulating fractional-order chaotic systems
[82], etc. Considering the scope of research, they are not enumerated here and only fundamental
FC and FO control related tools are reviewed.

3.
3.1

Evaluation and comparison
Comparison I

To evaluate the collected tools, several groups of benchmark problems and inputs are designed. For
the FO control toolboxes, the following problems are used,
(1) Baseline model: first order transfer function,
gb (s) =

1
,
s+1

whose time domain analytical solution of its step response is: y(t) = 1 − e−t ;
(2) Impulse response of half order integrator:
1
ghint (s) = √ ,
s
whose time domain analytical solution is:
(3) TF with a half order pole:

√1 ;
πt

ghp (s) = √

1
,
s+1

whose time domain analytical solution is:
√
√
1
1
√ E 1 , 1 (− t), or equivalently, √ − et erfc[ t];
t 2 2
πt
(4) The commensurate order TF:
gcom (s) =

6s1.2 + s0.8 + 2s0.4 + 3
;
5s1.6 + s0.8 + 2

10

November 25, 2015

International Journal of Control

FO˙tool˙eval

(5) Step response of the irrational order TF:
gir (s) =

s

√

2s
5

√

3

+1

+ 3s

√

2

+1

.

The accuracy is quantified by the conventional integral absolute error (IAE) criteria, S =
RT
0 |e(t)|dt. All comparison have been kept as fair as possible. The numerical values of the time domain analytical solution using Matlab built-in functions are assumed to be accurate and is adopted
as the baseline. The computational errors when Ts = 0.05 are summarized in table 2, where the
row indices represent the methods numbered in table 1, and the column indices represent the test
problems respectively. Besides, ‘M’ denotes the Matlab built-in TF and ‘-’ means the underlying
method is not applicable for the test problem. Two sample plots of the step responses of problems
1 and 5 are shown in figures 6 and 7. For problems 4 and 5, since analytical solution is hard to
obtain, all methods are compared to the values computed by fotf.
For the impulse response of the half-order integrator, the first point is ignored for error calculating
because it is infinity. Two graphic views of the comparison are shown in figures 8a and 8b, with
Ts=0.01 and Ts=0.1 respectively.
As stated in [45], irid fod() uses finite dimensional (z) TF as the approximation method.
Hence, the order of the (z) TF has impact on the approximation accuracy. The error listed in
table 2 is based on the 10th order approximation. An illustrative plot is shown in figure 11. The
sampling time also has impact on its accuracy. The anti-intuitive fact is that relatively greater Ts
gives higher accuracy. A heat map of the error on the field of Ts=0.01:0.001:0.1 and order=3:30
is plotted in figure 12. At some particular high orders, “rank deficient” would occur during the
call of prony() inside irid fod(). Users can choose appropriate orders according to their specific
accuracy requirement.
The analytical expression of M-L function is a summation of infinite terms. Hence, it is not
surprising to see the numerical computation induced error in the results.
3.2

Comparison II

fderiv(), glfdiff(), fourier diffint() and FIT are integration /differentiation tools for functions. For this group of tools, the following two problems are designed to compare the performance.
(1) Half order derivative of the function y(t) = 3t on the interval of [0, 5], whose analytical
solution is,
0.5
0 Dt y

=

3Γ(2) √
t.
Γ(1.5)

(2) 0.75 order integration of the function y(t) =
−0.75
y
0 Dt

=

√

(9)

t, whose analytical solution is,

Γ(1.5) 1.25
t .
Γ(2.25)

(10)

The time steps are all set to 0.01 sec. It can be seen that fourier diffint() performs not
as well as other methods although a big number of Fourier and Gaussian coefficients have been
assigned (default values are 260 and 520 for identity polynomial). Its performance on a 3rd order
polynomial is better. The results are plotted in figures 9 and 10. Quantitative comparison including
computational error and averaged elapsed time (for 20 runs each) are listed in table 3, for the above
two problems respectively.

11

November 25, 2015

3.3

International Journal of Control

FO˙tool˙eval

Comparison III

Although the simulation of fractional order pseudo state space (SS) models can be achieved indirectly, some toolboxes do provide the direct simulation capability, such as the CRONE toolbox and
FSST. Since the function frac ss in CRONE toolbox only adopts the input of commensurate order
systems, for comparison purposes, the following commensurate order pseudo state space model is
selected,

(0.7) 

  
x1
0
1
x1
0
=
+
u
x2
−0.1
−0.2
x
1
2



 x1
y = 0.1 0.3
x2


(11)

To involve more tools into comparison, the FO integrator blocks in the FOTF and Ninteger toolboxes are used to represent the above fractional differential equations in Simulink, as shown in
figure 13. The comparison of the unit step responses computed by the four toolboxes are plotted
in figure 14, from which it can be seen that the result obtained using FSST (1 sec for step size)
has bigger difference from the others. However, since analytical solution is not easy to obtain, it is
insufficient to claim which method gives highest accuracy. Hence, quantitative comparison is not
provided. As an alternative, users can transform the above FO SS model to an FO transfer function
model, assuming zero initial conditions,
G(s) = C(sα I − A)−1 B =

3s0.7 + 1
.
10s1.4 + 2s0.7 + 1

(12)

Thus, the NILT scripts can be used to compute the numerical solution, which has relatively higher
reliability according to the authors observation.

4.

Comments for selection

A tricky part for the simulation of fractional order systems is that even if the system is broken down
to the bottom layer, i.e. the analytical solution, it usually still involves the computation of M-L
functions, which still needs to rely on the numerical tools or scripts. From the comparison, it can be
seen that in the category of integrating/differentiating a function, glfdiff and FIT outperform other
tools in terms of accuracy; in the category of control system simulation, NILT always provides
higher accuracy. However, other toolboxes has advantages, for example, ninteger and CRONE
toolbox provide integrator blocks in Simulink, which makes the simulation of nonlinear systems
possible.

5.

Conclusion

In this paper, a comprehensive review of the Matlab based numerical tools for fractional calculus
and fractional order controls is presented. Quantitative evaluation of the selected tools is conducted.
The summarized description and numerical comparison are designated to serve as a reference and
guidance for readers when selecting tools for specific applications.

12

November 25, 2015

International Journal of Control

FO˙tool˙eval

References
[1] J. T. Machado, V. Kiryakova, and F. Mainardi, “Recent history of fractional calculus,” Commun Nonlinear Sci Numer Simulat, vol. 16, pp. 1140–1153, 2011.
[2] K. S. Miller and B. Ross, An Introduction to the Fractional Calculus and Fractional Differential
Equations, 1st ed. New York, NY, USA: Wiley-Interscience, 1993.
[3] J. Sabatier, O. Agrawal, and J. T. Machado, Eds., Advances in Fractional Calculus: Theoretical
Developments and Applications in Physics and Engineering. Springer, 2007.
[4] R. L. Magin, Fractional calculus in bioengineering. Begell House Inc., 2006.
[5] B. J. West, Where madicine went wrong. World Scientific Publishing Company, 2006.
[6] B. J. West, M. Turalska, and P. Grigolini, Complex Networks: From Social Crises to Neuronal
Avalanches. Wiley-VCH Verlag GmbH & Co. KGaA, 2014, pp. 509–524.
[7] T. G. Lewis, Book of Extremes: Why the 21st Century isnt Like the 20th Century. Switzerland:
Springer, 2014.
[8] Z. Ding, C. W. Granger, and R. F. Engle, “A long memory property of stock market returns
and a new model,” Journal of Empirical Finance 1, pp. 83–106, 1993.
[9] B. G. Malkiel, A Random Walk Down Wall Street, 7th ed. W. W. Norton & Company, 1999.
[10] C. A. Monje, Y. Q. Chen, B. M. Vinagre, D. Xue, and V. Feliu, Fractional order Systems and
Controls: Fundamentals and Applications. Begell House Inc., 2006.
[11] C. Yin, Y. Q. Chen, and S. M. Zhong, “Fractional-order sliding mode based extremum seeking
control of a class of nonlinear systems,” Automatica, vol. 50, pp. 3173–3181, 2014.
[12] Z. Li and Y. Q. Chen, “Identification of linear fractional order systems using the relay feedback
approach,” in 2014 American Control Conference (ACC), Portland, OR, USA, 2014.
[13] R. Barbosa and J. A. T. Machado, “Implementation of discrete-time fractional-order controllers based on LS approximations,” Acta Polytechn Hung, vol. 3, no. 4, pp. 5–22, 2006.
[14] G. W. Bohannan, “Analog fractional order controller in temperature and motor control applications,” Journal of Vibration and Control, vol. 14, no. 9, pp. 1487–1498, 2008.
[15] C. Jiang, T. T. Hartley, J. Carletta, and R.J.Veillette, “A systematic approach for implementing fractional-order operators and systems,” IEEE Journal on Emerging and Selected Topics
in Circuits and Systems, vol. 3, no. 3, pp. 301–312, 2013.
[16] S. Das and I. Pan, Fractional order signal processing: Introductory Concepts and Applications.
Springer, 2012.
[17] D. Valério and J. S. da Costa, “Time-domain implementatio of fractional order controllers,”
Control Theory and Applications, IEE Proceedings -, vol. 152, no. 5, pp. 539–552, Sept 2005.
[18] I. Podlubny, Fractional differential equations. Academic Press, 1999.
[19] H. Bode, Network Analysis and Feedback Amplifier Design. New York, USA: D. Van Nostrand
Company, Inc., 1945.
[20] R. L. Bagley and R. A. Calico, “Fractional order state equations for the control of viscoelastically damped structure,” Journal of Guidance, vol. 14, no. 2, pp. 304–310, 1989.
[21] M. Axtell and M. E. Bise, “Fractional calculus applications in control systems,” in Nat.
Aerospace and Electronic Conf., New York, NY, USA, 1990, pp. 563–566.
[22] I. Podlubny, “Fractional-order systems and PIλ Dµ controllers,” IEEE Transactions on Automatic Control, vol. 44, no. 1, pp. 208–214, 1999.
[23] Y. Q. Chen, I. Petráš, and D. Xue, “Fractional order control - a tutorial,” in American Control
Conference, 2009. ACC ’09., St. Louis, MO, USA, June 2009, pp. 1397–1411.
[24] D. Valério and J. S. da Costa, “Ninteger: a non-integer control toolbox for Matlab,” in 1st
IFAC Workshop on Fractional Differentiation and Applications, Bordeaux, France, 2004, pp.
208–213.
[25] A. Oustaloup, P. Melchior, P. Lanusse, O. Cois, and F. Dancla, “The CRONE toolbox for
matlab,” Computer-Aided Control System Design, pp. 190–195, 2000.
[26] N. Yousfi, P. Melchior, C. Rekik, N. Derbel, and A. Oustaloup, “Design of centralized CRONE
13

November 25, 2015

[27]

[28]
[29]
[30]
[31]
[32]
[33]
[34]
[35]

[36]

[37]

[38]
[39]
[40]
[41]
[42]
[43]
[44]

[45]

[46]
[47]

International Journal of Control

FO˙tool˙eval

controller combined with MIMO-QFT approach applied to non-square multivariable systems,”
Int. J. Comput. Appl, vol. 45, pp. 6–14, 2012.
A. Tepljakov, E. Petlenkov, and J. Belikov, “FOMCON: Fractional-order modeling and control
toolbox for MATLAB,” in 18th International Conference “Mixed Design of Integrated Circuits
and Systems”, Gliwice, Poland, 2011.
D. Xue and Y. Q. Chen, System Simulation Techniques with Matlab and Simulink. John
Wiley & Sons, Ltd., 2014.
T. Prabhakar, “A singular integral equation with a generalised Mittag-Leffler function in the
kernel,” Yokohama Math. J, vol. 19, pp. 7–15, 1971.
B. Barrowes, “Generalized hypergeometric function,” [Online] (Matlab Central) http://www.
mathworks.com/matlabcentral/fileexchange/5616, 2005.
J. Huntley, “Generation of random variates,” [Online] (Matlab Central) http://www.
mathworks.com/matlabcentral/fileexchange/35008, 2012.
V. Chaurasia and S. Pandey, “On the fractional calculus of generalized Mittag-Leffler function,” SCIENTIA Series A: Mathematical Sciences, vol. 20, pp. 113–122, 2010.
F. de Hoog, J. Knight, and A.N.Stokes, “An improved method for numerical inversion of
Laplace transforms,” SIAM Journal on Scientific Computing, vol. 3, no. 3, pp. 357–366, 1982.
J. Liang, “Control of linear time-invariant distributed parameter systems: from interger order
to fractional order,” Ph.D. dissertation, Utah State University, 2005.
J. Valsa and L. Branc̆ík, “Fractional order state equations for the control of viscoelastically
damped structure,” International Journal of Numerical Modelling: Electronic Networks, Devices and Fields, vol. 11, no. 3, pp. 153–166, 1998.
L. Branc̆ı́k, “An improved numerical inversion of two-dimensional Laplace transforms with
application to transient analysis of transmission lines,” in Proceedings EDS’99 Brno, Czech
Republic, 1999.
L. Branc̆ík, “Utilization of quotient-difference algorithm in FFT-based numerical ILT method,”
in Proceedings of the 11th International Czech-Slovak Scientific Conference Radioelektronika,
Czech Republic, 2001.
H. Sheng, Y. Li, and Y. Q. Chen, “Application of numerical inverse laplace transform algorithms in fractional calculus,” Journal The Franklin Institute, vol. 348, pp. 315–330, 2011.
I. Petráš, “Digital fractional order differentiator/integrator - IIR type,” [Online] (Matlab Central) http://www.mathworks.com/matlabcentral/fileexchange/3672, 2003 (updated 2011).
——, “Digital fractional order differentiator/integrator - FIR type,” [Online] (Matlab Central)
http://www.mathworks.com/matlabcentral/fileexchange/3673, 2003 (updated 2011).
——, “Digital fractional order differentiator/integrator - new IIR type,” [Online] (Matlab
Central) http://www.mathworks.com/matlabcentral/fileexchange/31358, 2011.
Y. Q. Chen and B. M. Vinagre, “A new IIR-type digital fractional order differentiator,” Signal
processing, vol. 83, pp. 2359–2365, 2003.
R. Caponetto, G. Dongola, L. Fortuna, and I. Petráš, Fractional Order Systems: Modelling
and Control Applications. World Scientific, 2010.
Y. Li, H. Sheng, and Y. Q. Chen, “Impulse response invariant discretization of a generalized
commensurate fractional order filter,” in Proceedings of the 8th World Congress on Intelligent
Control and Automation, Jinan, China, July 2010.
Y. Q. Chen, “Impulse response invariant discretization of fractional order integrators/differentiators,” [Online] (Matlab Central) http://www.mathworks.com/matlabcentral/
fileexchange/21342, 2008.
D. Xue, Y. Q. Chen, and D. Atherton, Linear feedback control - analysis and design with
Matlab 6.5. Society for Industrial and Applied Mathematics, 2009.
A. Oustaloup, F. Levron, B. Mathieu, and F. M. Nanot, “Frequency-band complex noninteger differentiator: Characterization and synthesis,” IEEE Trans. On Circ. And Sys. - I:
Fundamental Theory And Applications, vol. 47, no. 1, pp. 25–39, 2000.
14

November 25, 2015

International Journal of Control

FO˙tool˙eval

[48] F. M. bayat, “Fractional differentiator,” [Online] (Matlab Central) http://www.mathworks.
com/matlabcentral/fileexchange/13858, 2007.
[49] Jonathan, “Fractional derivative,” [Online] (Matlab Central) http://www.mathworks.com/
matlabcentral/fileexchange/45982, 2014.
[50] D. Xue and Y. Q. Chen, Modeling, Analysis and Design of Control Systems in Matlab and
Simulink. World Scientific, 2014.
[51] G. Papazafeiropoulos, “Fractional differentiation and integration,” [Online] (Matlab Central)
http://www.mathworks.com/matlabcentral/fileexchange/45877, 2014.
[52] T. M. Marinov, N. Ramirez, and F. Santamaria, “Fractional integration toolbox (FIT),” [Online] (Santamaria Lab) http://www.cbi.utsa.edu/FIT, 2013.
[53] ——, “Fractional integration toolbox,” Fractional Calculus & Applied Analysis, vol. 16, no. 3,
pp. 670–681, 2013.
[54] I. Petráš, “Discrete fractional-order PID controller,” [Online] (Matlab Central) http://www.
mathworks.com/matlabcentral/fileexchange/33761, 2011.
[55] N. Lachhab, F. Svaricek, F. Wobbe, and H. Rabba, “Fractional order PID controller (FOPID)
- toolbox,” in 2013 European Control Conference (ECC), Zuich, Switzerland, 2013, pp. 3694–
3699.
[56] H. Li, Y. Luo, and Y. Q. Chen, “A fractional order proportional and derivative (FOPD) motion
controller: Tuning rule and experiments,” IEEE Transactions on Control Systems Technology,
vol. 18, no. 2, pp. 516–520, 2010.
[57] Y. Luo and Y. Q. Chen, Fractional order motion controls. John Wiley & Sons, Inc., 2013.
[58] P. Apkarian and D. Noll, “Nonsmooth H∞ synthesis,” IEEE Transactions on Automatic Control, vol. 51, no. 1, pp. 71–86, 2006.
[59] E. Pisoni, A. Visioli, and S. Dormido, “An interactive tool for fractional order PID controllers,”
in IECON ’09. 35th Annual Conference of IEEE, 2009.
[60] J. L. Guzman, K. J. Astrom, S. Dormido, T. Hagglund, M. Berenguel, and Y. Piguet, “Interactive modules for PID control,” EEE Control Systems Magazine I, vol. 28, no. 5, pp. 118–134,
2008.
[61] C. Tricaud and Y. Q. Chen, “Solution of fractional order optimal control problems using SVDbased rational approximations,” in The 2009 American Control Conference (ACC), St. Louis,
MO, USA, 2009, pp. 1430–1435.
[62] T. Zhao, Y. Q. Chen, and Z. Li, “Fractional order nonlinear model predictive control using
RIOTS 95,” in The International conference on fractional differentiation and its applications
(ICFDA), Catania, Italy, 2014.
[63] C. Tricaud and Y. Q. Chen, “Solving fractional order optimal control problems in RIOTS 95 a general-purpose optimal control problem solver,” in Proceedings of the 3rd IFAC Workshop
on Fractional Differentiation and its Applications, Ankara, Turkey, 2008.
[64] C. Tricaud, “Solution of fractional optimal control problems,” [Online] (Matlab Central)
http://www.mathworks.com/matlabcentral/fileexchange/22196, 2009.
[65] A. Dzieliński and D. Sierociuk, “Simulation and experimental tools for fractional order control
education,” in Proceedings IFAC World Congress, Seoul, Korea, 2008, pp. 11 654–11 659.
[66] D. Sierociuk, “Fractional states-space toolkit (FSST),” [Online] http://www.ee.pw.edu.pl/
∼dsieroci/fsst/fsst.htm, 2003.
[67] C. F. Lorenzo and T. T. Hartley, “Variable order and distributed order fractional operators,”
Nonlinear Dynamics, vol. 29, no. 1, pp. 57–98, 2002.
[68] I. Podlubny, “Matrix approach to discrete fractional calculus,” Factional Calculus and Applied
Analysis, vol. 29, no. 4, pp. 281–296, 2000.
[69] ——, “Matrix approach to discretization of ODEs and PDEs of arbitrary real order,” [Online]
(Matlab Central) http://www.mathworks.com/matlabcentral/fileexchange/22071, 2008.
[70] D. Sierociuk, “Fractional variable order derivative Simulink toolkit,” [Online] (Matlab Central)
http://www.mathworks.com/matlabcentral/fileexchange/38801, 2012.
15

November 25, 2015

International Journal of Control

FO˙tool˙eval

[71] D. Valério, “Variable order derivatives,” [Online] (Matlab Central) http://www.mathworks.
com/matlabcentral/fileexchange/24444, 2009.
[72] F. M. bayat, “Root-locus plot of fractional order systems,” [Online] (Matlab Central) http://
www.mathworks.com/matlabcentral/fileexchange/20577, 2008.
[73] J. A. T. Machado, “Root locus of fractional linear systems,” Commun Nonlinear Sci Numer
Simulat, vol. 16, pp. 3855–3862, 2011.
[74] Z. Li, “Fractional order root locus,” [online] (Matlab Central) http://www.mathworks.com/
matlabcentral/fileexchange/50458, 2015.
[75] R. M. Corless and D. J. Jeffrey, “Graphing elementary Riemann surface,” SIGSAM Bulletin,
vol. 32, no. 1, pp. 11–17, 1998.
[76] H. M. Farkas and I. Kra, Riemann Surfaces, 2nd ed. New York: Springer-Verlag, 1980.
[77] I. Petráš, Engineering Education and Research Using MATLAB. InTech, 2011, ch. Fractional
Derivatives, Fractional Integrals, and Fractional Differential Equations in Matlab, pp. 239–264.
[78] H. M. Ozaktas, Z. Zalevsky, and M. A. Kutay, The Fractional Fourier Transform. John Wiley
& Sons, 2001.
[79] V. A. Narayanana and K. Prabhu, “The fractional fourier transform: theory, implementation
and error analysis,” Microprocessors and Microsystems, vol. 27, pp. 511–521, 2003.
[80] G. Germano, “Mittag-Leffler random number generator,” [Online] (Matlab Central) http://
www.mathworks.com/matlabcentral/fileexchange/19392, 2008 (updated 2014).
[81] D. Chen, Y. Q. Chen, and D. Xue, “Digital fractional order Savitzky-Golay differentiator,”
IEEE Transactions on Circuits and Systems - II: Express Briefs, vol. 58, no. 11, pp. 758–762,
2011.
[82] I. Petráš, Fractional-Order Nonlinear Systems: modeling, analysis and simulation. Springer
Science & Business Media, Jointly published by Higher Education Press, 2011.
[83] The CRONE Team, “The CRONE toolbox homepage,” Internet, http://www.ims-bordeaux.
fr/CRONE/toolbox, Feb/2014.
[84] I. Podlubny, “Mittag-leffler function,” [Online] (Matlab Central) http://www.mathworks.com/
matlabcentral/fileexchange/8738, 2005.
[85] S. Mukhopadhyay, “Mittag-leffler function, M-file, cmex DLL, and S-function,” [Online] (Matlab Central) http://www.mathworks.com/matlabcentral/fileexchange/20731, 2008.
[86] Y. Q. Chen, “Generalized Mittag-Leffler function,” [Online] (Matlab Central) http://www.
mathworks.com/matlabcentral/fileexchange/20849, 2008.
[87] R. Garrappa, “The Mittag-Leffler function,” [Online] (Matlab Central) http://www.
mathworks.com/matlabcentral/fileexchange/48154, 2014.
[88] Juraj, “Numerical inversion of Laplace transforms in Matlab,” [Online] (Matlab Central)
http://www.mathworks.com/matlabcentral/fileexchange/32824, 2011.
[89] Y. Q. Chen, “Oustaloup-recursive-approximation for fractional order differentiators,” [Online]
(Matlab Central) http://www.mathworks.com/matlabcentral/fileexchange/3802, 2003.
[90] A. Tepljakov, “Fractional-order calculus based identification and control of linear dynamic
systems,” Master’s thesis, Tallinn University of Technology, 2011.

16

November 25, 2015

Name
fotf
ninteger
Crone
FOMCON
mlf
ml func
ml fun
gml fun
ml
NILT
INVLAP
dfod1,2,3
irid fod ...
ora foc
fderiv
glfdiff

Typical usage
FO control toolbox
FC and FOC toolbox
FO control toolbox
FO modeling & control
2-param M-L func
1 ∼ 4 param M-L func
2-param M-L func
Generalized M-L func
1,2,3-param M-L func
Num Inverse of Laplace
Num Inverse of Laplace
Digital FO diff/int
Impulse Resp Invariant
Oustaloup-Rec-Approx
FO diff of r(t)
Finite Diff of G-L
FO diff of f(x)
FO integration toolbox
Discrete FO PID
FO PID

Sample syntax
s=fotf(’s’)

Author(s)
Dingyü Xue
D Valério
CRONE team
A Tepljakov
I. Podlubny
Dingyü Xue

Source
[50]
nid(k,a,[w1 w2],5,’crone’)
[24]
frac tf(1,frac poly exp(1,0.5)
[83]
sys foss = tf2ss(g)
[27]
y=mlf(a,b,-t)
[84]
y=ml func([a,b,r],-t)
[10]
y=ml fun(a,b,x,n,e)
S. Mukhopadhyay
[85]
gml fun(a,b,r,x,eps0)
YQ Chen
[86]
e= ML(x,a,b,r)
R Garrappa
[87]
Script based
L Branc̆ík
[36]
[t,y]=INVLAP(’1/s’,1,10,100)
Code by Juraj
[88]
sysdfod=dfod3(n,T,r)
I Petráš
[39]
df=irid fod(-.5,.1,5)
YQ Chen
[45]
ora foc(0.5,2,0.1,100)
YQ Chen
[89]
y=fderiv(0.5,r,Ts)
F. M. bayat
[48]
y1=glfdiff(y,t,r)
Dingyü Xue
[50]
fourier diffint
fourier diffint(f,x,..)
G Papazafeiropoulos
[51]
FIT
fracIntegrationSIM(...)
Marinov et al.
[52]
DFOC
DFOC(K,Ti,Td,m,d,Ts,n)
I Petráš
[54]
FOPID
—
Lachhab et al.
[55]
FOCP
Fractional optimal control
Calling RIOTS
C Tricaud et al.
[64]
FSST
FO S-S Toolkit
Simulink blocks
D. Sierociuk
[66]
FVO
Fractional variable order ban(alpha,N,h)
Podlubny et al.
[68]
forlocus
RL plot of FO TFs
forlocus(num,den,l)
Zhuo Li et al.
[74]
‘Delay’ column denotes if the script/toolbox is able to handle time delay in the FO model.
‘MIMO’ column denotes if the script/toolbox is able to handle MIMO FO models.

Delay
X
Could
×
X

MIMO
Could
Could
X
X

N/A

N/A

X
X
N/A
N/A
N/A
N/A
N/A
N/A
N/A
N/A
—
×
X
N/A
N/A

×
×
N/A
N/A
N/A
N/A
N/A
N/A
N/A
×
×
X
X
N/A
N/A

FO˙tool˙eval

17

#
1
2
3
4
5a
5b
5c
5d
5e
6a
6b
7
8
9
10
11
12
13
14
15
16
17
18
19
The
The

International Journal of Control

Table 1.: Matlab based numerical tools for computation of fractional operations and fractional order controls.

International Journal of Control

FO˙tool˙eval

ninteger

Crone

FOMCON
fotf
Figure 1.: FOMCON’s relation to other numerical tools, [90].

1
0.9

Amplitude (C)

November 25, 2015

0.8

1.01

0.7

1

0.6

0.99
0.98

0.5

0.97
0.4
0.96
0.3

4

6

0.2

8

10

M
M −1 Ts

Using
Using Ts

0.1
0
0

1

2

3

4

5
Time (sec)

6

7

8

9

10

Figure 2.: Computation error of NILT caused by mis-assignment of Ts.

Figure 3.: The Simulink block set provided in FSST.

18

November 25, 2015

International Journal of Control

FO˙tool˙eval

1.4
1.2
1
0.8
0.6
0.4

Ts=2
Ts=1
Ts=0.1
Ts=0.01

0.2
0

0

5

10

15

20

25

30

35

40

45

50

Figure 4.: The impact of simulation step size on the FSST toolbox.

1.5

2.5
2

1

1.5
1

0.5
Im{w}

Im{s}

0.5
0

0
−0.5

−0.5

−1
−1.5

−1

−2
−1.5
−1.5

−1

−0.5

0
Re{s}

0.5

1

−2.5
−5

1.5

−4

−3

−2

−1
Re{w}

(a) On s-plane.

(b) On w-plane.

Figure 5.: The RL plot of equation (8) on different planes.

19

0

1

2

International Journal of Control

FO˙tool˙eval

1
Time domain anlytical
Matlab built−in tf()
fotf()
ninteger
Crone
mlf
ml_fun
NILT
INVLAP
irid

Amplitude

0.8

0.6
0.98
0.96
0.4

0.94
0.92
0.9
0.88

0.2

0.86
0.84
0

2
0

1

2.2

2

3

2.4
4

2.6

5
Time (sec.)

2.8
6

3

7

3.2

8

9

Figure 6.: Comparison of the step responses of problems 1.

1.4

1.2

1

1.15
Amplitude

November 25, 2015

0.8

fotf()
ninteger
Crone
NILT
INVLAP

1.1

0.6

0.4

1.05

0.2

1
6

0

0

5

8

10

12

14

10

16
15

18

20

22

24

20

Time (sec.)

Figure 7.: Comparison of the step responses of problem 5.

20

25

International Journal of Control

FO˙tool˙eval

20
18

Anlytical
fotf()
ninteger
irid_fod()
NILT

1.9
1.8

16

1.7
14
1.6
Amplitude

12

1.5
1.4

10

0.11 0.12 0.13 0.14 0.15 0.16 0.17
8
6

0.25

4

0.15

0.2

0.1

2

7
0

0

1

2

3

4

7.5

8

5
Time (sec.)

8.5

6

9
7

9.5
8

10
9

(a) Ts=0.01
10
Anlytical
fotf()
ninteger
irid_fod()
NILT

9
1.6
8

1.4
1.2

7

1
0.8

6
Amplitude

November 25, 2015

0.6
0.2

5

0.4

0.6

0.8

1

4
0.24
0.22
0.2
0.18
0.16
0.14

3
2
1

9
0

0

1

2

3

4

9.2

5
Time (sec.)

9.4
6

9.6
7

9.8
8

10
9

(b) Ts=0.1

Figure 8.: Comparison of the impulse responses of the half order integrator.

21

November 25, 2015

International Journal of Control

FO˙tool˙eval

15
y
D0.5y Analytical
0.5

D

y fderiv

0.5

Half derivative of y = 3t

D

y glfdiff

0.5

D

10

y fourier_diffint

0.5

D

y FIT

5

0
0

1

2

3

4

5

t [sec]

Figure 9.: Comparison of half derivative of function y(t) = 3t, using fderiv(), glfdiff(),
fourier diffint() and FIT respectively.

7
y
0.75

I

6
t
√
0 .7 5 order integ ra l o f y =

y Analytical

0.75

I

y fderiv

0.75

I

5

y glfdiff

0.75

I

y fourier_diffint

0.75

I

y FIT

4

3

2

1

0

0

1

2

3

4

5

t [sec]

Figure 10.: Comparison of 0.75th order integration of function y(t) =
fderiv(), glfdiff(), fourier diffint() and FIT respectively.

22

√

t, using analytical solution,

November 25, 2015

International Journal of Control

FO˙tool˙eval

Ts=0.01
Impulse resp

10
8

0.2

6

0.1

4

0

Analytical
th
5 order
th
10 order
th
20 order

9.5

2

9.6

9.7

9.8

9.9

10

0
0

1

2

3

4

5

6

7

8

9

10

Ts=0.1

Impulse resp

3
0.21
0.2
0.19
0.18
0.17

2

Analytical
th
5 order
th
10 order
th
20 order

1
9.5

9.6

9.7

9.8

9.9

0
0

1

2

3

4

5

6

7

8

9

10

Time [sec]

Figure 11.: The Ts and order impact on irid fod().

Order of the z TF

5
10
15
20
25
30
0.01

0.02

0.03

0.04

0.05

0.06

0.07

0.08

0.09

0.1

Ts [sec]

Figure 12.: The heat map of approximation error of irid fod() versus the order and sampling
time.

23

November 25, 2015

International Journal of Control

FO˙tool˙eval

Figure 13.: The Simulink block diagrams for simulating the FO pseudo state space model in equation
(11).

1.4
1.2
1
0.8
0.6
CRONE
FSST
FOTF
Ninteger

0.4
0.2
0

0

5

10

15

20

25

30

35

40

45

50

Figure 14.: Comparison of the simulation results of the FO SS model obtained with different
toolboxes.

24

November 25, 2015

International Journal of Control

FO˙tool˙eval

Table 2.: Quantitative evaluation results for the test problems 1∼5.
❍❍
Error
❍❍
Method
❍
M
1
2
3
5a
5c
6a
6b
8

1

2

3

4

5

0
1.4955
3.18×10−13
0.4956
0
8.62×10−12
0.0016
0.0059
0.5327

8.4176
2.5287
2.9627
0.0236
0.0012
0.0071

6.6813
0.3831
1.4254
4.69×10−4
1.08×10−10
0.0206
2.49×10−5
0.2189

“0”
9.8434
10.454
6.1528
3.4722
-

“0”
2.5519
3.1857
2.4477
1.5042
-

Table 3.: Quantitative comparison of function int/diff tools.
PP Methods
PP
Analytical
CriteriaPPP
Error 1
Elapsed T1
0.0001
Error 2
Elapsed T2
0.0001

fderiv()

glfdiff()

fourier diffint()

FIT

140.5000
1.4028
339.7973
1.4105

1.8232
0.0029
2.1208
0.0029

792.4660
0.0874
250.5433
0.0893

0.0000
0.0209
0.0743
0.0201

25

