#### Abstract

Fractional Calculus (FC) goes back to the beginning of the theory of differential calculus. Nevertheless, the application of FC just emerged in the last two decades, due to the progress in the area of chaos that revealed subtle relationships with the FC concepts. In the field of dynamical systems theory some work has been carried out but the proposed models and algorithms are still in a preliminary stage of establishment. Having these ideas in mind, the paper discusses FC in the study of system dynamics and control. In this perspective, this paper investigates the use of FC in the fields of controller tuning, legged robots, redundant robots, heat diffusion, and digital circuit synthesis.

#### 1. Introduction

The generalization of the concept of derivative to noninteger values of goes back to the beginning of the theory of differential calculus. In fact, Leibniz, in his correspondence with Bernoulli, L’Hôpital and Wallis (1695), had several notes about the calculation of . Nevertheless, the development of the theory of Fractional Calculus (FC) is due to the contributions of many mathematicians such as Euler, Liouville, Riemann, and Letnikov [1–3].

The FC deals with derivatives and integrals to an arbitrary order (real or, even, complex order). The mathematical definition of a derivative/integral of fractional order has been the subject of several different approaches [1–3]. For example, the Laplace definition of a fractional derivative of a signal is
where , . The Grünwald-Letnikov definition is given by (*ℜ*):
where is the Gamma function and is the time increment. However, (1.2) shows that fractional-order operators are “global’’ operators having a memory of all past events, making them adequate for modeling memory effects in most materials and systems.

The Riemann-Liouville definition of the fractional-order derivative is (): where is the Gamma function of .

Based on the proposed definitions it is possible to calculate the fractional-order integrals/derivatives of several functions (Table 1). Nevertheless, the problem of devising and implementing fractional-order algorithms is not trivial and will be the matter of the following sections.

In recent years FC has been a fruitful field of research in science and engineering [1–6]. In fact, many scientific areas are currently paying attention to the FC concepts and we can refer its adoption in viscoelasticity and damping, diffusion and wave propagation, electromagnetism, chaos and fractals, heat transfer, biology, electronics, signal processing, robotics, system identification, traffic systems, genetic algorithms, percolation, modeling and identification, telecommunications, chemistry, irreversibility, physics, control systems as well as economy, and finance [7–18].

Bearing these ideas in mind, Sections 2–6 present several applications of FC in science and engineering. In Section 2, it is presented the application of FC concepts to the tuning of PID controllers and, in Section 3, the application of a fractional-order PD controller in the control of the leg joints of a hexapod robot. Then in Section 4, it is presented the fractional dynamics in the trajectory control of redundant manipulators. Next, in Section 5, it is introduced the fractional characteristics of heat diffusion along a media and, in Section 6 it is shown the application of FC to circuit synthesis using evolutionary algorithms. Finally, the main conclusions are presented in Section 7.

#### 2. Tuning of PID Controllers Using Fractional Calculus Concepts

The PID controllers are the most commonly used control algorithms in industry. Among the various existent schemes for tuning PID controllers, the Ziegler-Nichols (Z-N) method is the most popular and is still extensively used for the determination of the PID parameters. It is well known that the compensated systems, with controllers tuned by this method, have generally a step response with a high percent overshoot. Moreover, the Z-N heuristics are only suitable for plants with monotonic step response.

In this section, we study a methodology for tuning PID controllers such that the response of the compensated system has an almost constant overshoot defined by a prescribed value. The proposed method is based on the minimization of the integral of square error (ISE) between the step responses of a unit feedback control system, whose open-loop transfer function is given by a fractional-order integrator and that of the PID compensated system [7].

Figure 1 illustrates the fractional-order control system that will be used as reference model for the tuning of PID controllers. The open-loop transfer function is defined as : where is the gain crossover frequency, that is, . The parameter is the slope of the magnitude curve, on a log-log scale, and may assume integer as well as noninteger values. In this study we consider , such that the output response may have a fractional oscillation (similar to an underdamped second-order system). This transfer function is also known as the Bode’s ideal loop transfer function since Bode studies on the design of feedback amplifiers in the 1940s [19].

The Bode diagrams of amplitude and phase of are illustrated in Figure 2. The amplitude curve is a straight line of constant slope dB/dec, and the phase curve is a horizontal line positioned at rad. The Nyquist curve is simply the straight line through the origin, arg rad.

This choice of gives a closed-loop system with the desirable property of being insensitive to gain changes. If the gain changes, the crossover frequency will change, but the phase margin of the system remains rad, independent of the value of the gain. This can be seen from the curves of amplitude and phase of Figure 2.

The closed-loop transfer function of fractional-order control system of Figure 1 is given by

The unit step response of is given by the expression:

For the tuning of PID controllers, we address the fractional-order transfer function (2.2) as the reference system [8]. With the order and the crossover frequency we can establish the overshoot and the speed of the output response, respectively. For that purpose we consider the closed-loop system shown in Figure 3, where and are the PID controller and the plant transfer functions, respectively.

The transfer function of the PID controller is where is the error signal and is the controller’s output. The parameters , , and are the proportional gain, the integral time constant, and the derivative time constant of the controller, respectively.

The design of the PID controller will consist on the determination of the optimum PID set gains that minimize , the integral of the square error (ISE), defined as where is the step response of the closed-loop system with the PID controller (Figure 3) and is the desired step response of the fractional-order transfer function (2.2) given by (2.3).

To illustrate the effectiveness of proposed methodology we consider the third-order plant transfer function: with nominal gain .

Figure 4 shows the step responses and the Bode diagrams of phase of the closed-loop system with the PID for the transfer function for gain variations around the nominal gain () corresponding to , that is, for a variation up to 40% of its nominal value. The system was tuned for (), rad/s. We verify that we get the same desired iso-damping property corresponding to the prescribed () values.

**(a)**

**(b)**

In fact, we observe that the step responses have an almost constant overshoot independent of the variation of the plant gain around the gain crossover frequency . Therefore, the proposed methodology is capable of producing closed-loop systems robust to gain variations and step responses exhibiting an iso-damping property. The proposed method was tested on several systems revealing good results. It was also compared with other tuning methods showing comparable or superior results [8].

#### 3. Fractional P Control of a Hexapod Robot

Walking machines allow locomotion in terrain inaccessible to other type of vehicles, since they do not need a continuous support surface, but at the cost of higher requirements for leg coordination and control. For these robots, joint level control is usually implemented through a PID-like scheme with position feedback. Recently, the application of the theory of FC to robotics revealed promising aspects for future developments [9]. With these facts in mind, this section compares different Fractional P robot controller tuning, applied to the joint control of a walking system (Figure 5) with legs, equally distributed along both sides of the robot body, having each three rotational joints (i.e., ) [10].

During this study leg joint can be either mechanical actuated or motor actuated (Figure 5). For the mechanical actuated case, we suppose that there is a rotational pre-tensioned spring-dashpot system connecting leg links and . This mechanical impedance maintains the angle between the two links while imposing a joint torque [10].

Figure 5 presents the dynamic model for the hexapod body and foot-ground interaction. It is considered robot body compliance because walking animals have a spine that allows supporting the locomotion with improved stability. The robot body is divided in identical segments (each with mass ) and a linear spring-damper system (with parameters defined so that the body behaviour is similar to the one expected to occur on an animal) is adopted to implement the intrabody compliance [10]. The contact of the th robot feet with the ground is modelled through a nonlinear system [11], being the values for the parameters based on the studies of soil mechanics [11].

The general control architecture of the hexapod robot is presented in Figure 6 [12]. In this study we evaluate the effect of different P, *ℜ*, controller implementations for , while is a proportional controller with gain . For the P algorithm, implemented through a discrete-time 4th-order Padé approximation (, *∈**ℜ*, ), we have
where and are the proportional and derivative gains, respectively, and is the fractional order, for joint . Therefore, the classical PD^{1} algorithm occurs when the fractional order .

It is analysed the system performance of the different tuning, during a periodic wave gait at a constant forward velocity , for two cases: two leg joints are motor actuated and the ankle joint is mechanical actuated and the three leg joints are fully motor actuated [10].

The analysis is based on the formulation of two indices measuring the mean absolute density of energy per traveled distance () and the hip trajectory errors () during walking, according to

To tune the different controller implementations we adopt a systematic method, testing and evaluating several possible combinations of parameters, for all controller implementations. Therefore, we adopt the parameters that establish a compromise in what concerns the simultaneous minimisation of and . Moreover, it is assumed high-performance joint actuators, with a maximum actuator torque of Nm, and the desired angle between the foot and the ground (assumed horizontal) is established as . We tune the P joint controllers for different values of the fractional order while making .

We start by considering that leg joints 1 and 2 are motor actuated and joint 3 is mechanical actuated. For this case we tune the P joint controllers for different values of the fractional order , with step , namely, . Afterwards, we consider that joint 3 is also motor actuated, and we repeat the controller tuning procedure versus .

For the first situation under study, we verify that the value of (Figure 7), being the gains of the P controller , , , and the parameters of the mechanical spring-dashpot system for the ankle actuation , , presents the best compromise situation between the simultaneous minimisation of and .

Regarding the case when all joints are motor actuated, Figure 8 presents the best controller tuning for different values of . The experiments reveal the superior performance of the P controller for , with , , , , and , .

For the results are very poor and for , the hexapod locomotion is unstable. Furthermore, we conclude that the best case corresponds to all leg joints being motor actuated.

In conclusion, the experiments reveal the superior performance of the FO controller for and a robot with all motor actuated joints, as can be concluded analysing the curves for the joint actuation torques (Figure 9) and for the hip trajectory tracking errors and (Figure 10).

**(a)**

**(b)**

**(c)**

**(a)**

**(b)**

Since the objective of the walking robots is to walk in natural terrains, in the sequel it is examined how the different controller tunings behave under different ground properties, considering that all joints are motor actuated. For this case, and considering the previously tuning controller parameters, the values of are varied simultaneously through a multiplying factor that is varied in the range . This variation for the ground model parameters allows the simulation of the ground behaviour for growing stiffness, from peat to gravel [11].

The performance measure versus the multiplying factor of the ground parameters is presented on Figure 11. Analysing the system performance from the viewpoint of the index , it is possible to conclude that the best P implementation occurs for the fractional order . Moreover, it is clear that the performances of the different controller implementations are almost constant on all range of the ground parameters, with the exception of the fractional order . For this case, presents a significant variation with . Therefore, we conclude that the controller responses are quite similar, meaning that these algorithms are robust to variations of the ground characteristics [12].

#### 4. Fractional Dynamics in the Trajectory Control of Redundant Manipulators

A redundant manipulator is a robotic arm possessing more degrees of freedom (*dof*) than those required to establish an arbitrary position and orientation of the end effector. Redundant manipulators offer several potential advantages over non-redundant arms. In a workspace with obstacles, the extra degrees of freedom can be used to move around or between obstacles and thereby to manipulate in situations that otherwise would be inaccessible [20–23].

When a manipulator is redundant, it is anticipated that the inverse kinematics admits an infinite number of solutions. This implies that, for a given location of the manipulator’s gripper, it is possible to induce a self-motion of the structure without changing the location of the end effecter. Therefore, the arm can be reconfigured to find better postures for an assigned set of task requirements.

Several kinematic techniques for redundant manipulators control the gripper through the rates at which the joints are driven, using the pseudoinverse of the Jacobian [22, 24]. Nevertheless, these algorithms lead to a kind of chaotic motion with unpredictable arm configurations.

Having these ideas in mind, Section 4.1 introduces the fundamental issues for the kinematics of redundant manipulators. Based on these concepts, Section 4.2 presents the trajectory control of a three *dof* robot. The results reveal a chaotic behavior that is further analyzed in Section 4.3.

##### 4.1. Kinematics of Redundant Manipulators

A kinematically redundant manipulator has more *dof* than those required to define an arbitrary position and orientation of the gripper. In Figure 12 is depicted a planar manipulator with *∈**ℵ* rotational () joints that is redundant for . When a manipulator is redundant it is anticipated that the inverse kinematics admits an infinite number of solutions. This implies that, for a given location of the manipulator’s gripper, it is possible to induce a self-motion of the structure without changing the location of the gripper. Therefore, redundant manipulators can be reconfigured to find better postures for an assigned set of task requirements but, on the other hand, have a more complex structure requiring adequate control algorithms.

We consider a manipulator with degrees of freedom whose joint variables are denoted by . We assume that a class of tasks, we are interested in can be described by variables, and that the relation between and is given by where is a function representing the direct kinematics.

Differentiating (4.1) with respect to time yields where , , and . Hence, it is possible to calculate a path in terms of a prescribed trajectory in the operational space. We assume that the following condition is satisfied:

Failing to satisfy this condition usually means that the selection of manipulation variables is redundant and the number of these variables can be reduced. When condition (4.3) is verified, we say that the degree of redundancy of the manipulator is . If, for some we have then the manipulator is in a singular state. This state is not desirable because, in this region of the trajectory, the manipulating ability is very limited.

Many approaches for solving redundancy [25, 26] are based on the inversion of (4.2). A solution in terms of the joint velocities is sought as
where is one of the generalized inverses of the [26–28]. It can be easily shown that a more general solution to (4.2) is given by
where is the identity matrix and is a arbitrary joint velocity vector and is the pseudoinverse of the . The solution (4.6) is composed of two terms. The first term is relative to minimum norm joint velocities. The second term, the *homogeneous solution*, attempts to satisfy the additional constraints specified by . Moreover, the matrix allows the projection of in the null space of . A direct consequence is that it is possible to generate internal motions that reconfigure the manipulator structure without changing the gripper position and orientation [27–30]. Another aspect revealed by the solution of (4.6) is that repetitive trajectories in the operational space do not lead to periodic trajectories in the joint space. This is an obstacle for the solution of many tasks because the resultant robot configurations have similarities with those of a chaotic system.

##### 4.2. Robot Trajectory Control

The direct kinematics and the Jacobian of a 3-link planar manipulator with rotational joints ( robot) has a simple recursive nature according with the expressions: where is the length of link , , , and .

During all the experiments it is considered seconds, , and *. *

In the closed-loop pseudoinverse’s method the joint positions can be computed through the time integration of the velocities according with the block diagram of the inverse kinematics algorithm depicted in Figure 13, where represents the vector of reference coordinates of the robot gripper in the operational space.

Based on (4.7) we analyze the kinematic performances of the -robot when repeating a circular motion in the operational space with frequency rad s^{-1}, centre at distance and radius .

Figure 14 shows the joint positions for the inverse kinematic algorithm (4.5) for and . We observe that the following hold.

**(a)**

**(b)**

**(c)**

**(d)**

^{-1}.

##### 4.3. Analysis of the Robot Trajectories

In the previous subsection we verified that the pseudoinverse-based algorithm leads to unpredictable arm configurations. In order to gain further insight into the pseudoinverse nature several distinct experiments are devised in the sequel during a time window of 300 cycles. Therefore, in a first set of experiments we calculate the Fourier transform of the -robot joints velocities for a circular repetitive motion with frequency rad s^{-1}, radius , and radial distances .

Figure 15 shows versus the frequency ratio and the distance where represents the Fourier operator. Is verified an interesting phenomenon induced by the gripper repetitive motion because a large part of the energy is distributed along several subharmonics. These fractional-order harmonics (*foh*) depend on and making a complex pattern with similarities with those revealed by chaotic systems. Furthermore, we observe the existence of several distinct regions depending on .

**(a)**

**(b)**

**(c)**

**(d)**

For example, selecting in Figure 15 several distinct cases, namely for , we have the different signal Fourier spectra clearly visible in Figure 16.

**(a)**

**(b)**

**(c)**

In the authors’ best knowledge the *foh* are aspects of fractional dynamics [14, 15, 31], but a final and assertive conclusion about a physical interpretation is a matter still to be explored.

For joints velocities 1 and 3 the results are similar to the verified ones for joint velocity 2.

#### 5. Heat Diffusion

The heat diffusion is governed by a linear one-dimensional partial differential equation (PDE) of the form: where is the diffusivity, is the time, is the temperature, and is the space coordinate. However, (5.1) involves the solution of a PDE of parabolic type for which the standard theory guarantees the existence of a unique solution [16].

For the case of a planar perfectly isolated surface we usually apply a constant temperature at and analyzes the heat diffusion along the horizontal coordinate . Under these conditions, the heat diffusion phenomenon is described by a noninteger-order model: where is the space coordinate, is the boundary condition, and is the system transfer function.

In our study, the simulation of the heat diffusion is performed by adopting the Crank-Nicholson implicit numerical integration based on the discrete approximation to differentiation as [16, 17] where , , and are the increments and the integration indices for space and time, respectively.

##### 5.1. Control Strategies

The generalized PID controller has a transfer function of the form where and are the orders of the fractional integrator and differentiator, respectively. The constants , , and are correspondingly the proportional gain, the integral time constant, and the derivative time constant.

Clearly, taking we get the classical controllers, respectively.

The controller is more flexible and gives the possibility of adjusting more carefully the closed-loop system characteristics.

In the following two subsections, we analyze the system of Figure 17 by adopting the classical integer-order PID and a fractional PI, respectively.

##### 5.2. PID Tuning Using the Ziegler-Nichols Rule

In this subsection, we analyze the closed-loop system with a conventional PID controller given by the transfer function (5.4) with . Usually, the PID parameters () are tuned by using the so-called Ziegler-Nichols open loop (ZNOL) method [17]. The ZNOL heuristics are based on the approximate first-order plus dead-time model:

For the heat system, the resulting parameters are leading to the PID constants .

A step input is applied at m and the closed-loop response is analyzed for m, without actuator saturation (Figure 18). We verify that the system with a PID controller, tuned through the ZNOL heuristics, does not produce satisfactory results giving a significant overshoot and a large settling time , namely, , where represents the peak time and the rise time. We consider two indices that measure the response error, namely, the integral square error (ISE) and the integral time square error (ITSE) criteria defined as

We can use other performance criteria such as the integral absolute error (IAE) or the integral time absolute error (ITAE); however, in the present case, the ISE and the ITSE criteria have produced the best results and are adopted in the study.

In this case, the ZNOL PID tuning leads to the values . The poor results indicate again that the method of tuning may not be the most adequate for the control of the heat system.

In fact, the inherent fractional dynamics of the system lead us to consider other configurations. In this perspective, we propose the use of fractional controllers tuned by the minimization of the indices ISE and ITSE.

##### 5.3. PID^{β} Tuning Using Optimization Indices

In this subsection, we analyze the closed-loop system under the action of the PI controller given by the transfer function (5.4) with and . The fractional derivative term in (5.4) is implemented through a fourth-order Padé discrete rational transfer function. It used a sampling period of second.

The PI controller is tuned by the minimization of an integral performance index. For that purpose, we adopt the ISE and ITSE criteria.

A step reference input is applied at m and the output is analyzed for m, without actuator saturation. The heat system is simulated for 3000 seconds. Figure 19 illustrates the variation of the fractional PID parameters () as function of the order’s derivative , for the ISE and the ITSE criteria. The dots represent the values corresponding to the classical PID (ZNOL-tuning) addressed in the previous section.

**(a)**

**(b)**

**(c)**

The curves reveal that for the parameters () are slightly different, for the two ISE and ITSE criteria, while for they lead to almost similar values. This fact indicates a large influence of a weak-order derivative on system’s dynamics.

To further illustrate the performance of the fractional-order controllers a saturation nonlinearity is included in the closed-loop system of Figure 17 and inserted in series with the output of the controller . The saturation element is defined as

The controller performance is evaluated for and which corresponds to a system without saturation. We use the same fractional-PID parameters obtained without considering the saturation nonlinearity.

Figures 20 and 21 show the step responses of the closed-loop system and the corresponding controller output, for the PI tuned in the ISE and ITSE perspectives for and , respectively. The controller parameters correspond to the minimization of those indices leading to the values ISE: and ITSE: .

**(a)**

**(b)**

**(a)**

**(b)**

The step responses reveal a large diminishing of the overshoot and the rise time when compared with the integer PID, showing a good transient response and a zero steady-state error. The PI leads to better results than the classical PID controller tuned through the ZNOL rule. These results demonstrate the effectiveness of the fractional algorithms when used for the control of fractional-order systems. The step response and the controller output are also improved when the saturation level is diminished.

Figure 22 depicts the ISE and ITSE indices for , when and . We verify the existence of a minimum for and for the ISE and ITSE cases, respectively. Furthermore, the higher the the lower the value of the index.

**(a)**

**(b)**

Figures 23 and 24 show the variation of the settling time , the peak time , the rise time , and the percent overshoot , for the closed-loop response tuned through the minimization of the ISE and the ITSE indices, respectively.

**(a)**

**(b)**

**(c)**

**(d)**

**(a)**

**(b)**

**(c)**

**(d)**

In the ISE case , , and diminish rapidly for , while for the parameters increase smoothly. For the ITSE, we verify the same behavior for . On the other hand, increases smoothly for , while for it decreases very quickly, both for the ISE and the ITSE indices.

In conclusion, for we get the best controller tuning, superior to the performance revealed by the classical integer-order scheme.

#### 6. Circuit Synthesis Using Evolutionary Algorithms

In recent decades evolutionary computation (EC) techniques have been applied to the design of electronic circuits and systems, leading to a novel area of research called Evolutionary Electronics (EE) or Evolvable Hardware (EH). EE considers the concept for automatic design of electronic systems. Instead of using human conceived models, abstractions, and techniques, EE employs search algorithms to develop implementations not achievable with the traditional design schemes, such as the Karnaugh or the Quine-McCluskey Boolean methods.

Several papers proposed designing combinational logic circuits using evolutionary algorithms and, in particular, genetic algorithms (GAs) [32, 33] and hybrid schemes such as the memetic algorithms (MAs) [34].

Particle swarm optimization (PSO) constitutes an alternative evolutionary computation technique, and this paper studies its application to combinational logic circuit synthesis. Bearing these ideas in mind, the organization of this section is as follows. Section 6.1 presents a brief overview of the PSO. Section 6.2 describes the PSO-based circuit design, while Section 6.3 exhibits the simulation results.

##### 6.1. Particle Swarm Optimization

In literature about PSO the term ‘swarm intelligence’ appears rather often and, therefore, we begin by explaining why this is so.

Noncomputer scientists (ornithologists, biologists, and psychologists) did early research, which led into the theory of particle swarms. In these areas, the term “swarm intelligence’’ is well known and characterizes the case when a large number of individuals are able of accomplish complex tasks. Motivated by these facts, some basic simulations of swarms were abstracted into the mathematical field. The usage of swarms for solving simple tasks in nature became an intriguing idea in algorithmic and function optimization.

Eberhart and Kennedy were the first to introduce the PSO algorithm [35], which is an optimization method inspired in the collective intelligence of swarms of biological populations, and was discovered through simplified social model simulation of bird flocking, fishing schooling, and swarm theory.

In the PSO, instead of using genetic operators, as in the case of GAs, each particle (individual) adjusts its flying according with its own and its companions experiences. Each particle is treated as a point in a D-dimensional space and is manipulated as described in what follows in the original PSO algorithm:
where and are positive constants, and are two random functions in the range , represents the *i*th particle, is the best previous position (the position giving the best fitness value) of the particle, the symbol represents the index of the best particle among all particles in the population, and is the rate of the position change (velocity) for particle .

However, (6.1a) and (6.1b) represent the flying trajectory of a population of particles. Also, (6.1a) describes how the velocity is dynamically updated and (6.1b) the position update of the “flying” particles. Moreover, (6.1b) is divided in three parts, namely the momentum, the cognitive and the social parts. In the first part the velocity cannot be changed abruptly: it is adjusted based on the current velocity. The second part represents the learning from its own flying experience. The third part consists on the learning group flying experience [36].

The first new parameter added into the original PSO algorithm is the inertia weigh. The dynamic equation of PSO with inertia weigh is modified to be where constitutes the inertia weigh that introduces a balance between the global and the local search abilities. A large inertia weigh facilitates a global search while a small inertia weigh facilitates a local search.

Another parameter, called constriction coefficient , is introduced with the hope that it can insure a PSO to converge. A simplified method of incorporating it appears in (6.3), where is function of and as it is presented as follows: where , .

There are two different PSO topologies, namely, the global version and the local version. In the global version of PSO, each particle flies through the search space with a velocity that is dynamically adjusted according to the particle’s personal best performance achieved so far and the best performance achieved so far by all particles. On the other hand, in the local version of PSO, each particle’s velocity is adjusted according to its personal best and the best performance achieved so far within its neighborhood. The neighborhood of each particle is generally defined as topologically nearest particles to the particle at each side.

PSO is an evolutionary algorithm simple in concept, easy to implement and computationally efficient. Figures 25, 26, and 27 present a generic EC algorithm, a hybrid algorithm, more precisely a MA and the original procedure for implementing the PSO algorithm, respectively.

The different versions of the PSO algorithms are the real-value PSO, which is the original version of PSO and is well suited for solving real-value problems; the binary version of PSO, which is designed to solve binary problems; and the discrete version of PSO, which is good for solving the event-based problems. To extend the real-value version of PSO to binary/discrete space, the most critical part is to understand the meaning of concepts such as trajectory and velocity in the binary/discrete space.

Kennedy and Eberhart [35] use velocity as a probability to determine whether (a bit) will be in one state or another (zero or one). The particle swarm formula of (6.1a) remains unchanged, except that now and are integers in and a logistic transformation is used to accomplish this modification. The resulting change in position is defined by the following rule: where the function is a sigmoid limiting transformation and is a random number selected from a uniform distribution in the range .

##### 6.2. PSO Based Circuit Design

We adopt a PSO algorithm to design combinational logic circuits. A truth table specifies the circuits and the goal is to implement a functional circuit with the least possible complexity. Four sets of logic gates have been defined, as shown in Table 2, being *Gset 2* the simplest one (i.e., a RISC-like set) and *Gset 6* the most complex gate set (i.e., a CISC-like set). Logic gate named WIRE means a logical no-operation.

In the PSO scheme the circuits are encoded as a rectangular matrix **A** (*row* * column *) of logic cells as represented in Figure 28.

Three genes represent each cell: <*input1><input2><gate type>*, where *input1* and *input2* are one of the circuit inputs, if they are in the first column, or one of the previous outputs, if they are in other columns. The gate type is one of the elements adopted in the gate set. The chromosome is formed with as many triplets as the matrix size demands (e.g., triplets = ). For example, the chromosome that represents a matrix is depicted in Figure 29.

The initial population of circuits (particles) has a random generation. The initial velocity of each particle is initialized with zero. The following velocities are calculated applying (6.2a) and the new positions result from using (6.2b). This way, each potential solution, called particle, flies through the problem space. For each gene is calculated the corresponding velocity. Therefore, the new positions are as many as the number of genes in the chromosome. If the new values of the input genes result out of range, then a re-insertion function is used. If the calculated gate gene is not allowed a new valid one is generated at random. These particles then have memory and each keeps information of its previous best position (*pbest*) and its corresponding fitness. The swarm has the *pbest* of all the particles and the particle with the greatest fitness is called the global best (*gbest*).

The basic concept of the PSO technique lies in accelerating each particle towards its *pbest* and *gbest* locations with a random weighted acceleration. However, in our case we also use a kind of mutation operator that introduces a new cell in 10% of the population. This mutation operator changes the characteristics of a given cell in the matrix. Therefore, the mutation modifies the gate type and the two inputs, meaning that a completely new cell can appear in the chromosome.

To run the PSO we have also to define the number of individuals to create the initial population of particles. This population is always the same size across the generations, until reaching the solution.

The calculation of the fitness function in (6.6) has two parts, and , where measures the functionality and measures the simplicity. In a first phase, we compare the output produced by the PSO-generated circuit with the required values , according with the truth table, on a bit-per-bit basis. By other words, is incremented by one for each correct bit of the output until reaches the maximum value that occurs when we have a functional circuit. Once the circuit is functional, in a second phase, the algorithm tries to generate circuits with the least number of gates. This means that the resulting circuit must have as much genes <*gate type>* <*wire>* as possible. Therefore, the index , that measures the simplicity (the number of null operations), is increased by *one* (*zero*) for each *wire* (*gate*) of the generated circuit, yielding
where and represent the number of inputs and outputs of the circuit.

The concept of dynamic fitness function results from an analogy between control systems and the GA case, where we master the population through the fitness function. The simplest control system is the proportional algorithm; nevertheless, there can be other control algorithms, such as the proportional and the differential scheme.

In this line of thought, (6.6) is a static fitness function and corresponds to using a simple proportional algorithm. Therefore, to implement a proportional-derivative evolution the fitness function needs a scheme of the type [18] where is the differential fractional-order and is the “gain’’ of the dynamical term.

##### 6.3. Experiments and Results

A reliable execution and analysis of an EC algorithm usually requires a large number of simulations to provide a reasonable assurance that the stochastic effects are properly considered. Therefore, in this study are developed simulations for each case under analysis.

The experiments consist on running the three algorithms to generate a typical combinational logic circuit, namely, a 2-to-1 multiplexer (), a 1-bit full adder (*FA1*), a 4-bit parity checker (*PC4*) and a 2-bit multiplier (*MUL2*), using the fitness scheme described in (6.6) and (6.7). The circuits are generated with the gate sets presented in Table 2 and , , , and .

Figure 30 depicts the standard deviation of the number of generations to achieve the solution versus the average number of generations to achieve the solution for the algorithms , the circuits , and the gate sets . In these figure, we can see that the *MUL2* circuit is the most complex one, while the *PC4* and the are the simplest circuits. It is also possible to conclude that *Gset 6* is the less efficient gate set for all algorithms and circuits.

**(a)**

**(b)**

**(c)**

Figure 30 reveals that the plots follow a power law: Table 3 presents the numerical values of the parameters for the three algorithms.

In terms of versus , the MA algorithm presents the best results for all circuits and gate sets. In what concerns the other two algorithms, the PSO is superior (inferior) to the GA for complex (simple) circuits.

Figure 31 depicts the average processing time to obtain the solution versus the average number of generations to achieve the solution for the algorithms , the circuits and the gate sets . When analysing these charts it is clear that the PSO algorithm demonstrates to be around ten times faster than the MA and the GA algorithms.

**(a)**

**(b)**

**(c)**

These plots follow also a power law:

Table 3 shows parameters () and we can see that the PSO algorithm has the best values.

Figures 32 and 33 depict the standard deviation of the number of generations to achieve the solution and the average processing time to obtain the solution , respectively, versus the average number of generations to achieve the solution for the PSO algorithm using , the circuits *{**M2-1, FA1, PC4, MUL2**}*, and the gate sets . We conclude that leads to better results in particular for the *MUL2* circuit and for the .

Figures 34 and 35 present a comparison between and .

In terms of versus it is possible to say that the MA algorithm presents the best results. Nevertheless, when analysing Figure 31, that shows versus for reaching the solutions, we verify that the PSO algorithm is very efficient, in particular, for the more complex circuits.

The PSO-based algorithm for the design of combinational circuits follows the same profile as the other two evolutionary techniques presented in this paper.

Adopting the study of the versus for the three evolutionary algorithms, the MA algorithm presents better results over the GA and the PSO algorithms. However, in what concerns the processing time to achieve the solutions, the PSO outcomes clearly the GA and the MA algorithms. Moreover, applying the the results obtained are improved further in all gate sets and in particular for the more complex circuits.

#### 7. Conclusions

Fractional Calculus (FC) goes back to the beginning of the theory of differential calculus. Nevertheless, the application of FC just emerged in the last two decades, due to the progress in the area of chaos that revealed subtle relationships with the FC concepts.

Recently FC has been a fruitful field of research in science and engineering and many scientific areas are currently paying wider attention to the FC concepts. In the field of dynamical systems theory, some work has been carried out but the proposed models and algorithms are still in a preliminary stage of establishment. This article presented several case studies on the implementation of FC-based models and control systems, being demonstrated the advantages of using the FC theory in different areas of science and engineering. In fact, this paper studied a variety of different physical systems, namely

(i)tuning of PID controllers using fractional calculus concepts;(ii)fractional P control of a hexapod robot;(iii)fractional dynamics in the trajectory control of redundant manipulators;(iv)heat diffusion;(v)circuit synthesis using evolutionary algorithms.It has been recognized the advantageous use of this mathematical tool in the modeling and control of these dynamical systems, and the results demonstrate the importance of Fractional Calculus and motivate for the development of new applications.