Published: 31 March 2016

# Effect of fractional-order PID controller with acceleration feedback on a linear single degree-of-freedom oscillator

Jiangchuan Niu1
Yongjun Shen2
Shaopu Yang3
Xianghong Li4
1, 2, 3School of Mechanical Engineering, Shijiazhuang Tiedao University, Shijiazhuang, China
4Department of Mathematics and Physics, Shijiazhuang Tiedao University, Shijiazhuang, China
Corresponding Author:
Yongjun Shen
Views 69

#### Abstract

A linear single degree-of-freedom (SDOF) oscillator with fractional-order PID controller of acceleration feedback is investigated by the averaging method, and the approximately analytical solution is obtained. Moreover, the numerical solution of the system is obtained by the step-down order method and the power series method progressively. The effects of the parameters in fractional-order PID controller on the dynamical properties are characterized by some new equivalent parameters. The proportional component of fractional-order PID controller is characterized in the form of equivalent mass. The integral component of fractional-order PID controller is denoted in the form of the equivalent linear damping and equivalent mass. The differential component of fractional-order PID controller is denoted in the form of the equivalent linear negative damping and equivalent mass. Those equivalent parameters could distinctly illustrate the effects of the parameters in fractional PID controller on the dynamical response. A comparison between the analytical solution with the numerical results is made, and their satisfactory agreement verifies the correctness of the approximately analytical results. The effects of the parameters in fractional-order PID controller on control performance are further analyzed by some performance parameters of the time response. Finally, the robustness of the fractional-order PID controller based on acceleration feedback is demonstrated through the control of a SDOF quarter vehicle suspension model.

## 1. Introduction

Fractional-order calculus includes the fractional-order integration and fractional-order differentiation, which is a natural generalization of traditional integer-order counterpart. The fractional-order calculus theory and numerical methods have made great progress in recent years [1-7]. Since many engineering problems can be described by fractional-order differential or integral equations, the application of fractional-order calculus in engineering is becoming more and more extensive. One of the important aspects is to add fractional-order feedback in the controlled system to improve the control performance, where the advantages of good robustness and strong anti-noise ability with fractional-order feedback are considered. For example, Chen et al. [8] studied the influence of fractional-order damping on the dynamic characteristics of van der Pol system. Song et al. [9] addressed the robust consensus problem of fractional-order linear multi-agent systems with positive real uncertainty. Wu et al. [10] studied the application of fractional-order calculus in the control of semi-active air suspension. Zeng et al. [11] applied fractional-order control to a fish processing machine and evaluated the performance through physical experiment.

PID (Proportional Integral Derivative) controllers are the most popular controllers used in industry because of their simplicity, performance robustness, and the availability of many effective and simple tuning methods based on minimum plant model knowledge. Podlubny [12] extended the traditional integer-order PID controller to the fractional-order counterpart. Compared with the traditional integer-order PID controller, the fractional-order PID controller has great advantages in stability and dynamic performance. At present, the research on fractional-order PID controller was mainly focused on the design method, stability analysis and parameters optimization of the controller [13-19].

Due to the application of acceleration feedback control is very extensive [20-24], a linear single degree-of-freedom (SDOF) oscillator with fractional-order PID controller of acceleration feedback is studied in this paper. The paper is organized as follow. Based on the averaging method, Section 2 presents the approximately analytical solution of the SDOF oscillator with fractional-order PID controller of acceleration feedback, where the equivalent linear damping, equivalent linear negative damping and three equivalent mass parameters are defined. These five equivalent parameters could completely reveal the effects of the coefficients and the orders in fractional-order PID controller on convergence speed and frequency characteristic of the approximate solution. In Section 3 the numerical solution of the system is obtained by the step-down order method and the power series method progressively. And the comparison between the analytical and numerical solution is presented. In Section 4, the effects of the coefficients and the orders in fractional-order PID controller on system control performance are further analyzed by some performance parameters of the time response. In Section 5, the parameter tuning process of fractional-order PID is demonstrated through the control of a SDOF quarter vehicle suspension model. In the vehicle suspension model, the fractional-order integral controller is arranged in the negative feedback loop, and the fractional-order PD controller is arranged in the positive feedback loop. The main conclusions are drawn at last.

## 2. Approximately analytical solution

The mathematical model for fractional-order PID controller in time domain is:

1
$u\left(t\right)={K}_{p}e\left(t\right)+{K}_{i}{D}_{t}^{-\lambda }e\left(t\right)+{K}_{d}{D}_{t}^{\delta }e\left(t\right),$

where $u\left(t\right)$ and $e\left(t\right)$ are the output and input for controller respectively. ${K}_{p}$, ${K}_{i}$ and ${K}_{d}$ are the proportional constant, integration constant and differentiation constant respectively. $\lambda$ and $\delta$ are the orders of fractional-order integral and differential component respectively. Fractional-order calculus is a generalization of integer-order integration and differentiation to non-integer order with fundamental operator ${D}_{t}^{r}$. And the $r$-order continuous operator is defined as:

2
${D}_{t}^{r}=\left\{\begin{array}{ll}\frac{{d}^{r}}{d{t}^{r}},& \mathrm{R}\mathrm{e}\left(r\right)>0,\\ 1,& \mathrm{R}\mathrm{e}\left(r\right)=0,\\ {\int }_{0}^{t}\left(d\tau {\right)}^{-r},& \mathrm{R}\mathrm{e}\left(r\right)<0,\end{array}\right\$

where $r$ could be a real number, and $t$ is the upper bound of the integral operation.

Supposing that a fractional-order PID controller is arranged in the negative feedback loop, the considered SDOF oscillator with fractional-order PID controller of acceleration feedback is shown as:

3
$m\stackrel{¨}{x}\left(t\right)+kx\left(t\right)+c\stackrel{˙}{x}\left(t\right)+{K}_{p}\stackrel{¨}{x}\left(t\right)+{K}_{i}{D}_{t}^{-\lambda }\left[\stackrel{¨}{x}\left(t\right)\right]+{K}_{d}{D}_{t}^{\delta }\left[\stackrel{¨}{x}\left(t\right)\right]=0,$

where $m$ is the mass, $k$ is the linear stiffness, and $c$ is the linear damping. The initial values are generally considered as $x\left(0\right)={a}_{0}$ and $\stackrel{˙}{x}\left(0\right)=0$. The fractional orders are restricted as $\text{0}\le \lambda \le \text{1}$ and $\text{0}\le \delta \le \text{1}$. The necessity, advantage, and realization methods of fractional-order PID controller could be found in [12] and other references. There are many definitions about fractional-order derivative, and the Caputo’s definition is adopted here. Then the fractional-order components of Eq. (3) are defined as:

4
${D}_{t}^{\delta }\left[\stackrel{¨}{x}\left(t\right)\right]={D}_{t}^{2+\delta }\left[x\left(t\right)\right]=\frac{1}{\mathrm{\Gamma }\left(1-\delta \right)}{\int }_{0}^{t}\frac{x\mathrm{\text{'}}\mathrm{\text{'}}\mathrm{\text{'}}\left(\tau \right)}{\left(t-\tau {\right)}^{\delta }}d\tau ,$
5
${D}_{t}^{-\lambda }\left[\stackrel{¨}{x}\left(t\right)\right]={D}_{t}^{2-\lambda }\left[x\left(t\right)\right]=\frac{1}{\mathrm{\Gamma }\left(\lambda \right)}{\int }_{0}^{t}\frac{{x}^{\mathrm{\text{'}}\mathrm{\text{'}}}\left(\tau \right)}{\left(t-\tau {\right)}^{1-\lambda }}d\tau ,$

where $\mathrm{\Gamma }\left(z\right)$ is Gamma function satisfying $\mathrm{\Gamma }\left(z+1\right)=z\mathrm{\Gamma }\left(z\right)$.

Eq. (3) can be transformed into:

6
$\stackrel{¨}{x}\left(t\right)+{\omega }_{0}^{2}x\left(t\right)+2\epsilon \mu \stackrel{˙}{x}\left(t\right)+\epsilon {k}_{i}{D}_{t}^{2-\lambda }\left[x\left(t\right)\right]+\epsilon {k}_{d}{D}_{t}^{2+\delta }\left[x\left(t\right)\right]=0,$

where:

${\omega }_{0}=\sqrt{\frac{k}{m+{K}_{p}}},\mathrm{}\mathrm{}\mathrm{}2\epsilon \mu =\frac{c}{m+{K}_{p}},\epsilon {k}_{i}=\frac{{K}_{i}}{m+{K}_{p}},\epsilon {k}_{d}=\frac{{K}_{d}}{m+{K}_{p}}.$

Then Eq. (6) could be solved by averaging method.

The solution for Eq. (6) could be supposed as:

7
$x=a\left(t\right)\mathrm{c}\mathrm{o}\mathrm{s}\phi ,\stackrel{˙}{x}=-{\omega }_{0}a\left(t\right)\mathrm{s}\mathrm{i}\mathrm{n}\phi ,$

where $\phi ={\omega }_{0}t+\theta$. Based on the averaging method, one could obtain the standard equations:

8a
$\stackrel{˙}{a}\left(t\right)=-\frac{\epsilon }{{\omega }_{0}}\left\{{P}_{1}\left[a\left(t\right),\theta \right]+{P}_{2}\left[a\left(t\right),\theta \right]+{P}_{3}\left[a\left(t\right),\theta \right]\right\}\mathrm{s}\mathrm{i}\mathrm{n}\phi ,$
8b
$a\left(t\right)\stackrel{˙}{\theta }=-\frac{\epsilon }{{\omega }_{0}}\left\{{P}_{1}\left[a\left(t\right),\theta \right]+{P}_{2}\left[a\left(t\right),\theta \right]+{P}_{3}\left[a\left(t\right),\theta \right]\right\}\mathrm{c}\mathrm{o}\mathrm{s}\phi ,$

where:

${P}_{1}\left[a\left(t\right),\theta \right]=2\mu {\omega }_{0}a\left(t\right)\mathrm{s}\mathrm{i}\mathrm{n}\phi ,$
${P}_{2}\left[a\left(t\right),\theta \right]=-{k}_{i}{D}_{t}^{2-\lambda }\left[a\left(t\right)\mathrm{c}\mathrm{o}\mathrm{s}\phi \right],$
${P}_{3}\left[a\left(t\right),\theta \right]=-{k}_{d}{D}_{t}^{2+\delta }\left[a\left(t\right)\mathrm{c}\mathrm{o}\mathrm{s}\phi \right].$

The amplitude $a$ and the phase $\theta$ are the same small order variables as $\epsilon$, and they change slower than $\phi$. Accordingly, one could obtain the approximate forms of the amplitude $a$ and the phase $\theta$ in certain time interval. From the averaging method [25], one could select the time terminal $T$ as $T=2\pi$ when ${P}_{i}\left[a\left(t\right),\theta \right]$$\text{(}i=$1, 2, 3) is periodic function, or $T=\mathrm{\infty }$ when ${P}_{i}\left[a\left(t\right),\theta \right]$$\text{(}i=$1, 2, 3) is aperiodic one. Expanding the integrals, one could obtain another form of Eq. (8a):

9
$\stackrel{˙}{a}\left(t\right)={\stackrel{˙}{a}}_{1}\left(t\right)+{\stackrel{˙}{a}}_{2}\left(t\right)+{\stackrel{˙}{a}}_{3}\left(t\right),$

where:

10a
${\stackrel{˙}{a}}_{1}\left(t\right)=-\frac{\epsilon }{2\pi {\omega }_{0}}{\int }_{0}^{2\pi }{P}_{1}\left[a\left(t\right),\theta \right]\mathrm{s}\mathrm{i}\mathrm{n}\phi d\phi =-\epsilon \mu a\left(t\right),$
10b
${\stackrel{˙}{a}}_{2}\left(t\right)=-\underset{T\to \mathrm{\infty }}{\mathrm{l}\mathrm{i}\mathrm{m}}\frac{\epsilon }{T{\omega }_{0}}{\int }_{0}^{T}{P}_{2}\left[a\left(t\right),\theta \right]\mathrm{s}\mathrm{i}\mathrm{n}\phi d\phi ,$
10c
${\stackrel{˙}{a}}_{3}\left(t\right)=-\underset{T\to \mathrm{\infty }}{\mathrm{l}\mathrm{i}\mathrm{m}}\frac{\epsilon }{T{\omega }_{0}}{\int }_{0}^{T}{P}_{3}\left[a\left(t\right),\theta \right]\mathrm{s}\mathrm{i}\mathrm{n}\phi d\phi .$

In order to calculate Eq. (10b) and Eq. (10c), some important results about the definite integrals in references [26-27] are presented as:

11a
${B}_{1}=\underset{T\to \mathrm{\infty }}{\mathrm{l}\mathrm{i}\mathrm{m}}{\int }_{0}^{T}\frac{\mathrm{s}\mathrm{i}\mathrm{n}\left(\omega t\right)}{{t}^{p}}dt={\omega }^{p-1}\mathrm{\Gamma }\left(1-p\right)\mathrm{c}\mathrm{o}\mathrm{s}\left(\frac{p\pi }{2}\right),$
11b
${B}_{2}=\underset{T\to \mathrm{\infty }}{\mathrm{l}\mathrm{i}\mathrm{m}}{\int }_{0}^{T}\frac{\mathrm{c}\mathrm{o}\mathrm{s}\left(\omega t\right)}{{t}^{p}}dt={\omega }^{p-1}\mathrm{\Gamma }\left(1-p\right)\mathrm{s}\mathrm{i}\mathrm{n}\left(\frac{p\pi }{2}\right)$

Then one could obtain:

12a
${\stackrel{˙}{a}}_{2}\left(t\right)=\underset{T\to \mathrm{\infty }}{\mathrm{l}\mathrm{i}\mathrm{m}}\frac{\epsilon }{T{\omega }_{0}}{\int }_{0}^{T}{k}_{i}{D}_{t}^{2-\lambda }\left[a\left(t\right)\mathrm{c}\mathrm{o}\mathrm{s}\phi \right]\mathrm{s}\mathrm{i}\mathrm{n}\phi d\phi$
$\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}=-\frac{\epsilon {k}_{i}{\omega }_{0}a\left(t\right)}{\mathrm{\Gamma }\left(\lambda \right)}\underset{T\to \mathrm{\infty }}{\mathrm{l}\mathrm{i}\mathrm{m}}\frac{1}{T}{\int }_{0}^{T}\left[{\int }_{0}^{t}\frac{\mathrm{c}\mathrm{o}\mathrm{s}\left({\omega }_{0}\tau +\theta \right)}{\left(t-\tau {\right)}^{1-\lambda }}d\tau \right]\mathrm{s}\mathrm{i}\mathrm{n}\left({\omega }_{0}t+\theta \right)dt$
$\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}=-\frac{\epsilon {k}_{i}a\left(t\right)}{2}{\omega }_{0}^{1-\lambda }\mathrm{c}\mathrm{o}\mathrm{s}\left(\frac{1-\lambda }{2}\pi \right),$
12b
${\stackrel{˙}{a}}_{3}\left(t\right)=\underset{T\to \mathrm{\infty }}{\mathrm{l}\mathrm{i}\mathrm{m}}\frac{\epsilon }{T{\omega }_{0}}{\int }_{0}^{T}\left[{k}_{d}{D}_{t}^{2+\delta }\left[a\left(t\right)\mathrm{c}\mathrm{o}\mathrm{s}\phi \right]\right]\mathrm{s}\mathrm{i}\mathrm{n}\phi d\phi$
$\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}=\frac{\epsilon {k}_{d}{\omega }_{0}^{2}a\left(t\right)}{\mathrm{\Gamma }\left(1-\delta \right)}\underset{T\to \mathrm{\infty }}{\mathrm{l}\mathrm{i}\mathrm{m}}\frac{1}{T}{\int }_{0}^{T}\left[{\int }_{0}^{t}\frac{\mathrm{s}\mathrm{i}\mathrm{n}\left({\omega }_{0}\tau +\theta \right)}{\left(t-\tau {\right)}^{\delta }}d\tau \right]\mathrm{s}\mathrm{i}\mathrm{n}\left({\omega }_{0}t+\theta \right)dt$
$\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}=\frac{\epsilon {k}_{d}a\left(t\right)}{2}{\omega }_{0}^{\delta +1}\mathrm{s}\mathrm{i}\mathrm{n}\left(\frac{\delta }{2}\pi \right).$

Accordingly, one can get:

13
$\stackrel{˙}{a}\left(t\right)=-\epsilon \mu a\left(t\right)-\frac{\epsilon {k}_{i}a\left(t\right)}{2}{\omega }_{0}^{1-\lambda }\mathrm{c}\mathrm{o}\mathrm{s}\left(\frac{1-\lambda }{2}\pi \right)+\frac{\epsilon {k}_{d}a\left(t\right)}{2}{\omega }_{0}^{\delta +1}\mathrm{s}\mathrm{i}\mathrm{n}\left(\frac{\delta }{2}\pi \right).$

Similarly, one could get the simplified result about generalized phase as:

14
$a\left(t\right)\stackrel{˙}{\theta }=-\frac{\epsilon {k}_{i}a\left(t\right)}{2}{\omega }_{0}^{1-\lambda }\mathrm{s}\mathrm{i}\mathrm{n}\left(\frac{1-\lambda }{2}\pi \right)-\frac{\epsilon {k}_{d}a\left(t\right)}{2}{\omega }_{0}^{\delta +1}\mathrm{c}\mathrm{o}\mathrm{s}\left(\frac{\delta }{2}\pi \right).$

Furthermore, the amplitude and the phase are solved as:

15a
$a\left(t\right)={a}_{0}\mathrm{exp}\left\{\left[-\epsilon \mu -\frac{\epsilon {k}_{i}}{2}{\omega }_{0}^{1-\lambda }\mathrm{sin}\left(\frac{\lambda }{2}\pi \right)+\frac{\epsilon {k}_{d}}{2}{\omega }_{0}^{\delta +1}\mathrm{sin}\left(\frac{\delta }{2}\pi \right)\right]t\right\},$
15b
$\theta =\left[-\frac{\epsilon {k}_{i}}{2}{\omega }_{0}^{1-\lambda }\mathrm{c}\mathrm{o}\mathrm{s}\left(\frac{\lambda }{2}\pi \right)-\frac{\epsilon {k}_{d}}{2}{\omega }_{0}^{\delta +1}\mathrm{c}\mathrm{o}\mathrm{s}\left(\frac{\delta }{2}\pi \right)\right]t.$

Substituting the original system parameters into Eq. (15), one can get:

16a
$a\left(t\right)={a}_{0}\mathrm{exp}\left[\left\{-\frac{1}{2\left[m+M\left(P\right)\right]}\left[c+C\left(I\right)+C\left(D\right)\right]\right\}t\right],$
16b
$\theta =-\frac{M\left(I\right)+M\left(D\right)}{2\left[m+M\left(P\right)\right]}{\omega }_{0}t,$

where the five new equivalent parameters are defined as:

17a
$M\left(P\right)={K}_{p},$
17b
$C\left(I\right)={K}_{i}{\omega }_{0}^{1-\lambda }\mathrm{s}\mathrm{i}\mathrm{n}\left(\frac{\lambda }{2}\pi \right),$
17c
$C\left(D\right)=-{K}_{d}{\omega }_{0}^{\delta +1}\mathrm{s}\mathrm{i}\mathrm{n}\left(\frac{\delta }{2}\pi \right),$
17d
$M\left(I\right)={K}_{i}{\omega }_{0}^{-\lambda }\mathrm{c}\mathrm{o}\mathrm{s}\left(\frac{\lambda }{2}\pi \right),$
17e
$M\left(D\right)={K}_{d}{\omega }_{0}^{\delta }\mathrm{c}\mathrm{o}\mathrm{s}\left(\frac{\delta }{2}\pi \right),$

where ${\omega }_{0}=\sqrt{k/\left(m+{K}_{p}\right)}=\sqrt{k/\left[m+M\left(P\right)\right]}$.

Accordingly, the approximately analytical solution for the SDOF oscillator with fractional-order PID controller of acceleration feedback could be established as:

18
$x\left(t\right)={a}_{0}\mathrm{exp}\left[\left\{-\frac{1}{2\left[m+M\left(P\right)\right]}\left[c+C\left(I\right)+C\left(D\right)\right]\right\}t\right]×\mathrm{c}\mathrm{o}\mathrm{s}\left\{{\omega }_{0}-\frac{M\left(I\right)+M\left(D\right)}{2\left[m+M\left(P\right)\right]}{\omega }_{0}\right\}t.$

From Eq. (18) it could be found that $C\left(I\right)$ increases the damping of Eq. (3) so that $C\left(I\right)$ is defined as the equivalent damping. $C\left(D\right)$ decreases the damping of Eq. (3) and it could be defined as the equivalent negative damping. Moreover, $M\left(P\right)$, $M\left(I\right)$ and $M\left(D\right)$ will increase the mass of Eq. (3), and decrease the natural frequency of Eq. (3), so that $M\left(P\right)$, $M\left(I\right)$ and $M\left(D\right)$ are all defined as the equivalent masses. The effects of the fractional-order PID controller on the dynamical properties of SDOF oscillator could be summarized as follow.

1) The proportional component of fractional-order PID controller could increase the mass and reduce the natural frequency of the system.

2) The integral component of fractional-order PID controller would increase the damping and limit the amplitude of the system by ${K}_{i}{\omega }_{0}^{1-\lambda }\mathrm{s}\mathrm{i}\mathrm{n}\left(\lambda \pi /2\right)$, and decrease the natural frequency of the system by $-{K}_{i}{\omega }_{0}^{1-\lambda }\mathrm{c}\mathrm{o}\mathrm{s}\left(\lambda \pi /2\right)/2\left(m+{K}_{p}\right)$.

3) The differential component of fractional-order PID controller would decrease the damping and increase the amplitude of the system by $-{K}_{d}{\omega }_{0}^{\delta +1}\mathrm{s}\mathrm{i}\mathrm{n}\left(\delta \pi /2\right)$, and decrease the natural frequency of the system by $-{K}_{d}{\omega }_{0}^{\delta +1}\mathrm{c}\mathrm{o}\mathrm{s}\left(\delta \pi /2\right)/2\left(m+{K}_{p}\right)$.

From Eq. (17), it could also be found that the function of fractional-order PID controller is not only related to the feedback of the system, but also on the coefficients ${K}_{p}$, ${K}_{i}$, ${K}_{d}$ and the orders $\lambda$ and $\delta$. That is listed here.

1) The proportional component of fractional-order PID controller regulates the mass and the natural frequency of the system by ${K}_{p}$.

2) The integral component of fractional-order PID controller affects the dynamic characteristics of the system by ${K}_{i}$ and $\lambda$. When $\lambda =\text{0}$, the integral component degrades as linear mass ${K}_{i}$, and decreases the natural frequency of the system by $-{K}_{i}{\omega }_{0}/2\left(m+{K}_{p}\right)$. When $\lambda$ changes from 0 to 1, the integral component will gradually degrade into linear damping. When $\lambda =\text{1}$, the integral component degrades into linear damping ${K}_{i}$.

3) The differential component of fractional-order PID controller affects the dynamic characteristics of the system by ${K}_{d}$ and $\delta$. When $\delta =\text{0}$, the differential component degrades into linear mass ${K}_{d}$, and decreases the natural frequency of the system by $-{K}_{d}{\omega }_{0}/2\left(m+{K}_{p}\right)$. When $\delta$ changes from 0 to 1, the differential component gradually changes into linear negative damping. When $\delta =\text{1}$, the differential component degrades into linear negative damping, and decreases the damping of the system by $-{K}_{d}{\omega }_{0}^{2}$.

Supposing that the fractional-order PID controller is arranged in the positive feedback loop, the effect of fractional-order PID controller on the dynamic characteristics will be reversed, which will yield the negation of the equivalent parameters. This conclusion is easy to verify.

## 3. Comparisons between the approximate analytical solution and numerical solution

An illustrative example is studied herein as defined by the system parameters: ${a}_{0}=$–1, $m=\mathrm{}$5, $k=\mathrm{}$10, $c=\text{0.3,}$${K}_{p}=\text{0.5,}$${K}_{i}=\text{1,}$${K}_{d}=\text{0.5,}$$\lambda =\text{0.6}$ and $\delta =\text{0.4.}$ Based on the approximately analytical solution as Eq. (18), one could analytically obtain the time history of the displacement for Eq. (3) denoted by the solid line shown in Fig. 1.

Fig. 1Comparison between the approximate analytical solution and numerical solution

In order to verify the precision of the approximately analytical solution, we also present the numerical results. The relationship to the explicit numerical approximation of the power series method in reference [1-2] is:

19
${D}^{p}\left[y\left({t}_{n}\right)\right]\approx {h}^{-p}\sum _{j=0}^{n}{C}_{j}^{p}y\left({t}_{n-j}\right),$

where ${t}_{n}=n×h$ is the sample points, $h$ is the sample step, and ${C}_{j}^{p}$ is the binomial coefficient with the iterative relationship as:

20
${C}_{0}^{p}=1,\mathrm{}\mathrm{}\mathrm{}{C}_{j}^{p}=\left(1-\frac{1+p}{j}\right){C}_{j-1}^{p}.$

However, we can’t directly use the power series method for Eq. (3). Because the fractional-order differentiation and integration are linear operations:

21
${D}_{t}^{p}\left[\gamma f\left(t\right)+\mu g\left(t\right)\right]=\gamma {D}_{t}^{p}f\left(t\right)+\mu {D}_{t}^{p}g\left(t\right).$

Eq. (3) can be transformed into:

22
$\left(m+{K}_{p}\right){D}_{t}^{2-\delta }\left[x\left(t\right)\right]+c{D}_{t}^{1-\delta }\left[x\left(t\right)\right]+k{D}_{t}^{-\delta }\left[x\left(t\right)\right]+{K}_{i}{D}_{t}^{2-\lambda -\delta }\left[x\left(t\right)\right]+{K}_{d}{D}_{t}^{2}\left[x\left(t\right)\right]=0,$

by step-down order.

According to Eqs. (19), (20) and (22), the numerical scheme for Eq. (3) can be shown as:

23a
$x\left({t}_{n}\right)=y\left({t}_{n-1}\right)h-\sum _{j=1}^{n}{C}_{j}^{1}x\left({t}_{n-j}\right),$
23b
$y\left({t}_{n}\right)=-\left[\left(m+{K}_{p}\right){z}_{1}\left({t}_{n-1}\right)+c{z}_{2}\left({t}_{n-1}\right)+k{z}_{3}\left({t}_{n-1}\right)+{K}_{i}{z}_{4}\left({t}_{n-1}\right)\right]\frac{h}{{K}_{d}}$
$\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}-\sum _{j=1}^{n}{C}_{j}^{1}y\left({t}_{n-j}\right),$
23c
${z}_{1}\left({t}_{n}\right)=-\left[\left(m+{K}_{p}\right){z}_{1}\left({t}_{n-1}\right)+c{z}_{2}\left({t}_{n-1}\right)+k{z}_{3}\left({t}_{n-1}\right)+{K}_{i}{z}_{4}\left({t}_{n-1}\right)\right]\frac{{h}^{\delta }}{{K}_{d}}$
$\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}-\sum _{j=1}^{n}{C}_{j}^{\delta }{z}_{1}\left({t}_{n-j}\right),$
23d
${z}_{2}\left({t}_{n}\right)=\left[y\left({t}_{n}\right)\right]{h}^{\delta }-\sum _{j=1}^{n}{C}_{j}^{\delta }{z}_{2}\left({t}_{n-j}\right),$
23e
${z}_{3}\left({t}_{n}\right)=\left[x\left({t}_{n}\right)\right]{h}^{\delta }-\sum _{j=1}^{n}{C}_{j}^{\delta }{z}_{3}\left({t}_{n-j}\right).$

For $\lambda +\delta \le \text{1}$, one could obtain:

23f
${z}_{4}\left({t}_{n}\right)=-\left[\left(m+{K}_{p}\right){z}_{1}\left({t}_{n}\right)+c{z}_{2}\left({t}_{n}\right)+k{z}_{3}\left({t}_{n}\right)+{K}_{i}{z}_{4}\left({t}_{n-1}\right)\right]\frac{{h}^{\lambda +\delta }}{{K}_{d}}$
$\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}-\sum _{j=1}^{n}{C}_{j}^{\lambda +\delta }{z}_{4}\left({t}_{n-j}\right).$

If $1<\lambda +\delta \le \text{2}$, one should use:

23g
${z}_{4}\left({t}_{n}\right)=\left[y\left({t}_{n}\right)\right]{h}^{\lambda +\delta -1}-\sum _{j=1}^{n}{C}_{j}^{\lambda +\delta -1}{z}_{4}\left({t}_{n-j}\right).$

Here we select $h=$0.001. And the initial values for the fractional-order terms could be roughly selected as [27]:

24
${D}_{t}^{p}\left[x\left(0\right)\right]={a}_{0}\mathrm{c}\mathrm{o}\mathrm{s}\left(\frac{p}{2}\pi \right).$

As a comparison, the time history by numerical solution is also shown in Fig. 1 by the circles. From the observation of Fig. 1, it could be concluded that the approximately analytical solution agrees very well with the numerical results.

## 4. Analysis of the control performance

In this section, we will study the effects on system control performance of the coefficients ${K}_{p}$, ${K}_{i}$ and ${K}_{d}\text{,}$ and the orders $\lambda$ and $\delta$ in fractional-order PID controller. The effects of fractional-order PID controller parameters on the rise time and peak value of system will be emphatically analyzed.

### 4.1. Rise time

The rise time ${t}_{r}$ which is the time for the system response curve to reach the steady-state amplitude. Based on Eq. (18), one could get:

25
${a}_{0}\mathrm{exp}\left[\left\{-\frac{1}{2\left[m+M\left(P\right)\right]}\left[c+C\left(I\right)+C\left(D\right)\right]\right\}t\right]×\mathrm{c}\mathrm{o}\mathrm{s}\left\{{\omega }_{0}t-\frac{M\left(I\right)+M\left(D\right)}{2\left[m+M\left(P\right)\right]}{\omega }_{0}t\right\}=0.$

In order to get the simplified form, another two equivalent parameters are defined, i.e. the equivalent negative stiffness $K\left(I\right)=-M\left(I\right){\omega }_{0}^{2}$ and $K\left(D\right)=-M\left(D\right){\omega }_{0}^{2}$. The effects on system dynamics of the equivalent negative stiffness and equivalent mass are consistent, which can be proved by the following method. By using the Taylor series expansion one could obtain:

26a
${F}_{1}\left[M\left(I\right)+M\left(D\right)\right]=\sqrt{\frac{k}{m+M\left(P\right)+M\left(I\right)+M\left(D\right)}}\approx {\omega }_{0}-\frac{M\left(I\right)+M\left(D\right)}{2\left[m+M\left(P\right)\right]}{\omega }_{0},$
26b
${F}_{2}\left[K\left(I\right)+K\left(D\right)\right]=\sqrt{\frac{k+K\left(I\right)+K\left(D\right)}{m+M\left(P\right)}}\approx {\omega }_{0}+\frac{K\left(I\right)+K\left(D\right)}{2\left[m+M\left(P\right)\right]{\omega }_{0}},$

where ${F}_{1}\left(z\right)=\sqrt{k/\left[m+M\left(P\right)+z\right]}$ and ${F}_{2}\left(z\right)=\sqrt{\left(k+z\right)/\left[m+M\left(P\right)\right]}$ and $z$ is a variable.

Letting $M$ be the equivalent mass of system and $M=m+M\left(P\right)$, $K$ be the equivalent stiffness of system and $K=k+K\left(I\right)+K\left(D\right)$, ${\omega }_{n}$ be the equivalent natural frequency of system and ${\omega }_{n}=\sqrt{K/M}$, from Eqs. (25) and (26) one could obtain:

27
${t}_{r}=\frac{\pi }{2\left\{{\omega }_{0}-\frac{M\left(I\right)+M\left(D\right)}{2\left[m+M\left(P\right)\right]}{\omega }_{0}\right\}}\approx \frac{\pi }{2\sqrt{\frac{K}{M}}}=\frac{\pi }{2{\omega }_{n}}.$

Obviously, the rise time of system is related with the stiffness and the equivalent mass. The rise time will become longer with the increase of ${K}_{p}$. The fractional-order PID controller will be changed into the traditional integer-order counterpart when $\lambda =\text{1}$ and $\delta =\text{1}$. And the rise time will remain unchanged with the change of ${K}_{i}$ and ${K}_{d}$. Conversely, when $\lambda \ne \text{1}$ and $\delta \ne \text{1}$, the rise time will last longer with the increase of ${K}_{i}$ or ${K}_{d}$.

In order to analyze the effect of the fractional orders more intuitively, other two parameters are fixed as ${K}_{i}=$1 and ${K}_{d}=$1. Henceforth we shall adopt the same initial value, except the parameter is specified as a variable. If the fractional order $\lambda$ in integral component is changed from 0 to 1, one could obtain the equivalent natural frequency curve which is shown as Fig. 2. It could be found that the equivalent natural frequency tends to become larger with the increase of $\lambda$. Then the rise time will become shorter with the increase of $\lambda$.

Fig. 2Equivalent natural frequency curve with λ

Fig. 3Equivalent natural frequency curve with δ

When the fractional order $\delta$ in differential component is changed from 0 to 1, one could obtain the equivalent natural frequency curve shown as Fig. 3. It could be found that the equivalent natural frequency remains almost unchanged for small $\delta$ and then tends to become larger with the increase of $\delta$. Then the rise time will remain approximately invariant at first and then become shorter with the increase of $\delta$.

### 4.2. Peak value

The time of peak value ${t}_{p}$ is the time for the system response curve to reach the first peak amplitude. From Eq. (18), one could get:

28
$\frac{dx\left(t\right)}{dt}\left|{\mathrm{}}_{t={t}_{p}}\right=0.$

After expanding and simplifying Eq. (28), one could obtain:

29
$\mathrm{s}\mathrm{i}\mathrm{n}\left\{{\omega }_{0}{t}_{p}-\frac{M\left(I\right)+M\left(D\right)}{2\left[m+M\left(P\right)\right]}{\omega }_{0}{t}_{p}\right\}=0.$

Then, one could get:

30
${t}_{p}=\frac{\pi }{\left\{{\omega }_{0}-\frac{M\left(I\right)+M\left(D\right)}{2\left[m+M\left(P\right)\right]}{\omega }_{0}\right\}}\approx \frac{\pi }{\sqrt{\frac{K}{M}}}=\frac{\pi }{{\omega }_{n}}.$

Substituting Eq. (30) into Eq. (18), the peak value could be obtained as:

31
${M}_{p}=-{a}_{0}\mathrm{exp}\left\{-\frac{C}{2\sqrt{MK}}\pi \right\}=-{a}_{0}{e}^{\left\{-\xi \pi \right\}},$

where $C=c+C\left(I\right)+C\left(D\right)$ and $\xi =C/2\sqrt{MK}$. Here $C$ is defined as the equivalent linear damping of system, and $\xi$ is defined as the equivalent damping ratio. Apparently, the peak value will become smaller with the increase of $\xi$.

The peak value will become larger with the increase of ${K}_{p}$. When the integral component coefficient ${K}_{i}$ is changed from 0 to 8, one could obtain the equivalent damping ratio curve of the fractional-order PID controller, which is shown in Fig. 4 by solid line. At the same time, the equivalent damping ratio curve of the traditional integer-order PID counterpart is shown by dotted line. From Fig. 4, it could be found that the equivalent damping ratio becomes larger with the increase of ${K}_{i}$. The peak values of the two kinds of controllers will be almost the same at the beginning stage, and then the peak value of the fractional-order PID controller will be smaller than the traditional integer-order counterpart.

If the differential component coefficient ${K}_{d}$ is changed from 0 to 5, one could obtain the equivalent damping ratio curve of the fractional-order PID controller, which is shown in Fig. 5 by solid line. Meanwhile, the equivalent damping ratio curve of the traditional integer-order PID counterpart is shown by dotted line. From Fig. 5, it could be found that the equivalent damping ratios become progressively smaller with the increase of ${K}_{d}\text{.}$ The peak values of the fractional-order PID controller will be smaller than the traditional integer-order counterpart at the beginning stage, and then it will be larger than the traditional integer-order counterpart. Obviously, when $\xi <\text{0}$, the system will be unstable.

If the fractional order $\lambda$ in integral component is changed from 0 to 1, one could obtain the equivalent damping ratio curve, which is shown in Fig. 6. From Fig. 6, it can be observed that the equivalent damping ratio will increase to the maximum and then decrease slightly. Therefore, the peak value will decrease to the minimum and then increase slightly.

When the fractional order $\delta$ in differential component is changed from 0 to 1, one could obtain the equivalent damping ratio curve, which is shown in Fig. 7. From Fig. 7, it could be found that the equivalent damping ratio become smaller progressively with the increase of $\delta$. Accordingly, the peak value will become larger with the increase of $\delta$.

Fig. 4Equivalent damping ratio curves with Ki

Fig. 5Equivalent damping ratio curves with Kd

Fig. 6Equivalent damping ratio curve with λ

Fig. 7Equivalent damping ratio curve with δ

## 5. Control example

When designing a closed-loop control system, one should determine the type of feedback and then select the components of controller and tune parameters according to the control objectives, so as to improve the control performance. The vehicle suspension could be simplified to a SDOF control system under some assumptions. Next, we will study the active control of a SDOF quarter vehicle suspension model by fractional-order PID controller with acceleration feedback.

The selected basic system parameters are $m=\text{240}$, $k=\text{16000}$, $c=\text{1000}$, ${a}_{0}=\text{–0.01}$, and one can get the damping ratio $\xi =\text{0.2552}$ and the natural frequency ${\omega }_{n}=\text{8.165}$. A well-designed vehicle suspension should provide isolation by reducing forces transmitted from vehicle axle to vehicle body [28]. In order to improve the ride comfort and guarantee a certain rapidity of the dynamic response, the control objectives are designed for that, the damping ratio $\xi =\text{0.4}$, and the rise time remain unchanged, when the sprung mass is changed to ${m}_{f}=\text{320}$.

Through the analysis of Section 4, one should arrange the fractional-order integral controller in the negative feedback loop, and arrange the fractional-order PD controller in the positive feedback loop.

According to the control objectives of system, one could get:

32a
$\frac{c+{K}_{i}{\omega }_{0}^{1-\lambda }\mathrm{s}\mathrm{i}\mathrm{n}\left(\frac{\lambda }{2}\pi \right)+{K}_{d}{\omega }_{0}^{\delta +1}\mathrm{s}\mathrm{i}\mathrm{n}\left(\frac{\delta }{2}\pi \right)}{2\sqrt{\left({m}_{f}-{K}_{p}\right)\left[k-{K}_{i}{\omega }_{0}^{2-\lambda }\mathrm{c}\mathrm{o}\mathrm{s}\left(\frac{\lambda }{2}\pi \right)+{K}_{d}{\omega }_{0}^{\delta +2}\mathrm{c}\mathrm{o}\mathrm{s}\left(\frac{\delta }{2}\pi \right)\right]}}=0.4,$
32b
$\sqrt{\frac{k-{K}_{i}{\omega }_{0}^{2-\lambda }\mathrm{c}\mathrm{o}\mathrm{s}\left(\frac{\lambda }{2}\pi \right)+{K}_{d}{\omega }_{0}^{\delta +2}\mathrm{c}\mathrm{o}\mathrm{s}\left(\frac{\delta }{2}\pi \right)}{{m}_{f}-{K}_{p}}}=8.165,$
32c
$\sqrt{\frac{k}{{m}_{f}-{K}_{p}}}={\omega }_{0}.$

Solving Eq. (32), one can obtain infinite sets of solutions. For convenience of calculations, one could let $\lambda =\text{0.5}$, $\delta =\text{0.5}$, and ${K}_{i}={K}_{d}{\omega }_{0}$, then obtain a set of solutions, which is ${K}_{p}=\text{80}$, ${K}_{i}=\text{140.4774}$ and ${K}_{d}=\text{17.2049}$. The time histories of the vibration displacements before and after control are shown in Fig. 8. Obviously, the rise time of the system response is guaranteed, and the peak value is reduced, when the system is under the fractional-order PID control. The results show that, the fractional-order PID controller can achieve the control demands.

Fig. 8Comparison of before and after control

Fig. 9Comparison between fractional-order and integer-order PID controllers

Letting $\lambda =\text{1}$, $\delta =\text{1}$, and ${K}_{i}={K}_{d}{\omega }_{0}^{2}$, one could obtain another set of solutions, which is ${K}_{p}=\text{80}$, ${K}_{i}=\text{283.8367}$, ${K}_{d}=\text{4.2576}$. The traditional integer-order PID controller can achieve the same control demands.

Next, the control performances are compared between the fractional-order PID and the traditional integer-order PID controller. When the stiffness of vehicle suspension is changed as $k=\text{10000}$, one could obtain the time histories of the vibration displacements with the PID controller of the above parameters, which are shown in Fig. 9. Apparently, the fractional-order PID controller has better robustness for the change of the stiffness.

## 6. Conclusions

A linear single degree-of-freedom oscillator with fractional-order PID controller of acceleration feedback is studied in this paper. The approximately analytical solution is obtained by the averaging method. The numerical solution is obtained by the step-down order method and the power series method. The accuracy of the method and the analytical results are demonstrated by comparison with the numerical solution. The effects of the parameters in fractional PID controller on the dynamical properties are characterized by equivalent linear damping, equivalent linear negative damping and three equivalent mass parameters. The effects of the parameters in fractional-order PID controller on system control performance are analyzed in detail by some performance parameters. Finally, the robustness of the fractional-order PID controller based on acceleration feedback is demonstrated through a control example. The research method and results may be helpful to the design of fractional-order PID controller of acceleration feedback, and provide a foundation for the control of free vibration with fractional-order PID controller.

#### References

• Petras I. Fractional-Order Nonlinear Systems. Higher Education Press, Beijing, 2011.
• Podlubny I. Fractional Differential Equations, Mathematics in Science and Engineering. Academic Press, New York, 1999.
• Shen Y. J., Yang S. P., Xing H. J., et al. Primary resonance of Duffing oscillator with fractional- order derivative. Communication in Nonlinear Science and Numerical Simulation, Vol. 17, Issue 7, 2012, p. 3092-3100.
• Rossikhin Y. A., Shitikova M. V. Application of fractional derivatives to the analysis of damped vibrations of viscoelastic single mass systems. Acta Mechanica, Vol. 120, Issues 1-4, 1997, p. 109-125.
• Yang S. P., Shen Y. J. Recent advances in dynamics and control of hysteretic nonlinear systems. Chaos Solitons and Fractals, Vol. 40, Issue 4, 2009, p. 1808-1822.
• Li C. P., Deng W. H. Remarks on fractional derivatives. Applied Mathematics and Computation, Vol. 187, Issue 1, 2007, p. 777-784.
• Cao J. X., Ding H. F., Li C. P. Implicit difference schemes for fractional diffusion equations. Communication on Applied Mathematics and Computation, Vol. 40, Issue 4, 2013, p. 61-74.
• Chen J. H., Chen W. C. Chaotic dynamics of the fractionally damped van der Pol equation. Chaos, Solitons and Fractals, Vol. 35, Issue 1, 2008, p. 188-198.
• Song C., Cao J. D., Liu Y. Z. Robust consensus of fractional-order multi-agent systems with positive real uncertainty via second-order neighbors information. Neurocomputing, Vol. 165, 2015, p. 293-299.
• Wu G. Q., Huang H. J., Ye G. H. Semi-active control of automotive air suspension based on fractional calculus. Transactions of the Chinese Society for Agricultural Machinery, Vol. 45, Issue 7, 2014, p. 19-24.
• Zeng Q. S., De Silva L. W. The application of fractional order control in an industrial fish processing machine. Control and Intelligent Systems, Vol. 40, Issue 3, 2012, p. 177-185.
• Podlubny I. Fractional-order systems and PID – controllers. IEEE Transactions on Automatic Control, Vol. 44, Issue 1, 1999, p. 208-214.
• Zeng G. Q., Chen J., Dai Y. X., et al. Design of fractional order PID controller for automatic regulator voltage system based on multi-objective extremal optimization. Neurocomputing, Vol. 160, 2015, p. 173-184.
• Zhong J. P., Li L. C. Tuning fractional-order (PID mu)-D-lambda controllers for a solid-core magnetic bearing system. IEEE Transactions on Control Systems Technology, Vol. 23, Issue 4, 2015, p. 1648-1656.
• Chen Y. Q., Petras I., Xue D. Y. Fractional order control – a tutorial. American Control Conference Hyatt Regency Riverfront, St. Louis, MO, USA, 2009.
• Saidi B., Amairi M., Najar S., et al. Bode shaping-based design methods of a fractional order PID controller for uncertain systems. Nonlinear Dynamics, Vol. 80, Issue 4, 2015, p. 1817-1838.
• Chen L. C., Zhao L., Li W., et al. Bifurcation control of bounded noise excited Duffing oscillator by a weakly fractional-order PID feedback controller. Nonlinear Dynamics, Vol. 83, Issue 1, 2016, p. 529-539.
• Chen L. C., Zhu W. Q. Stochastic jump and bifurcation of Duffing oscillator with fractional derivative damping under combined harmonic and white noise excitations. International Journal of Non-Linear Mechanics, Vol. 46, Issue 10, 2011, p. 1324-1329.
• Chen L. C., Li Z. S., Zhuang Q. J., et al. First-passage failure of single-degree-of-freedom nonlinear oscillators with fractional derivative. Journal of Vibration and Control, Vol. 19, Issue 14, 2013, p. 2154-2163.
• Wang Z. H., Zheng Y. G. The optimal form of the fractional-order difference feedbacks in enhancing the stability of a sDOF vibration system. Journal of Sound and Vibration, Vol. 326, Issues 3-5, 2009, p. 476-488.
• Zhang T. M., Li G. T., Liang L. Switch variable structure control of a flexible manipulator with acceleration feedback. Chinese Journal of Mechanical Engineering, Vol. 38, Issue 3, 2002, p. 39-42.
• Dumetz E., Dieulot J. Y., Barre P. J., et al. Control of an industrial robot using acceleration feedback. Journal of Intelligent and Robotic Systems: Theory and Applications, Vol. 46, Issue 2, 2006, p. 111-128.
• Wang H. P., Li J., Zhang K. Sup-resonant response of a nonautonomous maglev system with delayed acceleration feedback control. IEEE Transactions on Magnetics, Vol. 44, Issue 10, 2008, p. 2338-2350.
• An F., Chen W. D., Shao M. Q. Dynamic behavior of time-delayed acceleration feedback controller for active vibration control of flexible structures. Journal of Sound and Vibration, Vol. 333, Issue 20, 2014, p. 4789-4809.
• Sanders J. A., Verhulst F., Murdock J. Averaging Methods in Nonlinear Dynamical Systems. Springer, New York, 2007.
• Shen Y. J., Yang S. P., Xing H. J. Dynamical analysis of linear single degree-of-freedom oscillator with fractional-order derivative. Acta Phsysica Sinica, Vol. 61, Issue 11, 2012, p. 110505.
• Shen Y. J., Yang S. P., Sui C. Y. Analysis on limit cycle of fractional-order van der Pol oscillator. Chaos, Solitons and Fractals, Vol. 67, 2014, p. 94-102.
• Rajamani R. Vehicle Dynamics and Control. Springer, NewYork, 2012.

30 September 2015
Accepted
18 January 2016
Published
31 March 2016
SUBJECTS
Vibration generation and control
Keywords
fractional-order PID controller
averaging method
approximately analytical solution
Acknowledgements

This work was financially supported by the National Natural Science Foundation of China (No. 11372198), the Cultivation Plan for Innovation Team and Leading Talent in Colleges and Universities of Hebei Province (No. LJRC018), the Higher Level Talents Scientific Research of Hebei Province (No. GCC2014053), and the High Level Talent of Hebei Province (No. A201401001).