# Differential Evolution-Based PID Control of Nonlinear Full-Car Electrohydraulic Suspensions

Authors: Jimoh O. Pedro; Muhammed Dangor; Olurotimi A. Dahunsi; M. Montaz Ali
Journal: Mathematical Problems in Engineering
 (2013)

---

## Abstract

This paper presents a differential-evolution- (DE-) optimized, independent multiloop proportional-integral-derivative (PID) controller design for full-car nonlinear, electrohydraulic suspension systems. The multiloop PID control stabilises the actuator via force feedback and also improves the system performance. Controller gains are computed using manual tuning and through DE optimization to minimise a performance index, which addresses suspension travel, road holding, vehicle handling, ride comfort, and power consumption constraints. Simulation results showed superior performance of the DE-optimized PID-controlled active vehicle suspension system (AVSS) over the manually tuned PID-controlled AVSS and the passive vehicle suspension system (PVSS).


---

## Body

## 1. Introduction

The evolution of modern instrumentation and control techniques has made semi-active vehicle suspension systems (SAVSS) [1] and AVSS designs more promising. However, the biggest challenge remains the attainment of trade-offs between the conflicting performance criteria like suspension travel, ride comfort, road holding, and vehicle handling.Although there is a vast amount of literature documented about AVSS design, most are affected by at least one of the following limitations: use of linear suspension models, ignoring actuator dynamics and performance evaluation based on the 2 degree-of-freedom (dof) quarter-car model [2, 3]. This work presents a 7 dof full-car nonlinear suspension system that accounts for system complexities associated with inherent nonlinearities, actuator dynamics, and coupling interactions [4].Linear optimal control schemes like linear quadratic regulator (LQR), linear quadratic gaussian (LQG),ℋ∞, mixed ℋ2/ℋ∞, and linear parameter varying (LPV) control methods are well-developed control schemes that have been employed in SAVSS and AVSS designs [1, 4–8]. Their stability and robustness properties are more readily established but always limited when employed for complex nonlinear control schemes; they also normally assume time-invariant situation. Minimising the chosen performance objectives when using these control methods requires that some optimal feedback gains be obtained. Need for measurement of all state variables is another additional challenge. Moreover, solving AVSS control problem by the use of linear control schemes is usually based on the assumption of a broad bandwidth actuator, whose response is fast enough and its parameters can be linearised (the Jacobian way) within some operating regions [9].Application of nonlinear control methods like feedback linearization (FBL) [4], backstepping, and sliding mode control [10] methods is also well documented among the AVSS control design literature. These methods are able to overcome most of the challenges encountered in the implementation of the linear control methods, but challenges like difficulty in establishing system stability and repeated differentiation of nonlinear functions in backstepping lead to additional system complications, and exact FBL is not always possible because full knowledge about the system is unavailable. Moreover, FBL and sliding mode are usually susceptible to chattering that degrades system performance. Availability of good system dynamic model is usually necessary and performance is also readily affected by measurement noise [4, 10–12].Intelligent controllers provide new alternatives with better prospects but they inherit the setbacks of the nonlinear controllers since they are implemented in most cases in combination with the nonlinear control techniques [12–15].PID control is well known for its simplicity and ease with which its gains can be adjusted. It is also known to have setbacks in its lack of robustness to parameter variations and physical implementation limitations due to high-loop gains [2]. Evidence in the literature shows that tuning of the majority of the industrial PID controllers was either done manually or based on its default setting [16].Previous works have demonstrated the prospect for superior AVSS performance when PID control is combined with appropriate computational intelligence or global optimization methods, especially for the selection of PID optimal gains. Ziegler Nichols tuning method has been successfully employed over the years in tuning PID controllers but it often requires further fine tuning which has been done manually or intuitively [9, 13, 15, 17].Evolutionary-based optimization techniques such as genetic algorithms (GA) and particle swarm optimization (PSO) have been successfully applied for controller tuning of linear PID control of quarter-car AVSS [18–20]. Similarly, optimizations of ℋ∞ and LQR-based control methods were performed with success [21, 22]. Optimization has also been combined with intelligent control methods to obtain superior performance in some cases [23–25]. With the aid of these evolutionary-based optimization algorithms, the AVSS was able to meet hard design specifications that could not be achieved by conventional tuning methods. The drawback of the proposed methods was the weakness of the cost functions as  they were only based on improving ride comfort or vertical body displacement. There was no attempt made to resolve the other conflicting design requirements of an AVSS.Differential evolution is a heuristic, population set-based direct search global optimization algorithm with numerous applications [13, 17, 26–28]. The global optimization algorithm is inherently flexible and relatively simpler in comparison with other techniques. It gives better search space exploration characteristics with similar or even better results than previously employed optimization routines [26, 29, 30]. Secondly, it applies predefined conditions where new generations are only accepted if they produce enriched and better solutions. Those characteristics slow the routine down but helps to prevent early convergence [31].Four independent multi-loop PID controllers have been employed in this work to simultaneously control the full-car AVSS and to ensure actuator stability. Multiobjective DE direct search global optimization algorithm was employed in tuning the multi-loop PID controllers’ gains. The objective function simultaneously addressed the conflicting trade-off challenge in the AVSS design.


## 2. System Overview and Modelling

AVSS responds dynamically to road disturbance inputs by inducing relative motion between the body and the wheel through the force generated by the servo-hydraulic actuator. Obtaining the appropriate control voltage for the actuator requires an optimal trade-off between the design objectives in the presence of road disturbance inputs. The success of this process yields a suspension system that is adaptive to the road disturbance and other operating conditions.Figure1 presents the generic AVSS feedback control loop. The system consists of a controller issuing the command input signal to the actuator to generate a manipulating signal. Better command input signal is sent because the controller output is optimized by the use of evolutionary algorithms.Figure 1
AVSS feedback control loop.


### 2.1. Physical Modelling

Figure2 shows the seven dof full-car suspension system used in this study. Ms, Iθ, and Iϕ are the sprung mass, pitch inertia, and roll inertia of the vehicle, respectively; lf, lr, and tf are the distance from the vehicle front axle to the centre of gravity, the distance from the vehicle rear axle to the centre of gravity, and the half width of the vehicle, respectively. z, θ, and ϕ are the vehicle body heave, pitch, and roll, respectively; the displacement and mass of each wheel are denoted as ztij and muij respectively with (i,j) representing the position of the tyre, where i refers to the lateral position front f or rear r, and j to the longitudinal position right r or left l, and the respective road disturbance profile at each wheel is given as wij. The spring and damping force contributions of the suspension system are denoted by Fkij and Fbij, respectively, and those of the tyre are denoted by Fktij and Fbtij, respectively. Faij is the actuator force supplied between the chassis and each wheel. The suspension travels are denoted by zij [4].Figure 2
Simplified full-car model of the active vehicle suspension system.


### 2.2. Mathematical Modelling

The governing equations are derived using Newton-Euler approach, with the resultant force acting in each suspension system given as(1)Fij=Fkij+Fbij-Faij.The heave, pitch, and roll dynamics of the vehicle are [4]
(2)Msz¨=Ffr+Ffl+Frr+Frl,Iθθ¨=-Ffrlf-Ffllf+Frrlr+Frllr,Iϕϕ¨=tf[Ffr-Ffl+Frr-Frl].The dynamics of each wheel are determined as follows:(3)muijz¨tij=-Fkij-Fbij+Faij+Fktij+Fbtij.The suspension travel at each wheel is determined as(4)zfl=z-lfsinθ-tfsinϕ,zrl=z+lrsinθ-tfsinϕ,zfr=z-lfsinθ+tfsinϕ,zrr=z+lfsinθ+tfsinϕ.The suspension spring and damping forces are computed as(5)Fkij=ksijl(ztij-zij)+ksijnl(ztij-zij)3,Fbij=bsijl(z˙tij-z˙ij)+bsijnl|z˙tij-z˙ij|sgn(z˙tij-z˙ij)=bsijsym|z˙tij-z˙ij|with linear and nonlinear suspension spring coefficientsksijl and ksijnl, and linear, nonlinear, and symmetric damping coefficients bsijl, bsijnl, and bsijsym, respectively. The dynamics of the tyres are modelled as
(6)Fktij=ktij(wij-ztij),Fbtij=btij(wij-ztij),
where ktij and btij are the spring and damping coefficients of the wheel. The pressure developed by the servo-hydraulic actuators is modelled as [4, 32]
(7)P˙Lij=γijΦijxvij-βijPLij+αijAijz˙ij,
where
(8)αij=4βeijVtij,βij=αijCtpij,γij=CdijSij1ρij,Φij=ϕ1ij×ϕ2ij,ϕ1ij=sgn(Psij-sgn(xvij)PLij),ϕ2ij=|Psij-sgn(xvij)PLij|,wherexvij are the spool-valves’ displacements, Psij are the supply pressures, Aij are the pistons’ cross-sectional areas, Vtij are the volumes of the cylinders, βeij are the bulk moduli of the working fluids, Φij are the hydraulic load flows, Ctpij are the leakage coefficients.C
d
i
j are the discharge coefficients between the supply lines and the hydraulic cylinders and Sij are spool-valves’ area gradients. The valves’ positions are modelled with first-order lag as follows:
(9)x˙vij=1τij(Kvijui-xvij),whereKvij are the servo-valves’ gains and τij are the time constants. The disturbance is a deterministic road bump. The vehicle approaches a sinusoidal bump at 40 km/h. The heights of the bump on the right and left sides are 4 and 6 cm, respectively, such that roll motion is induced. The road profile at each wheel is given as
(10)wfr={a1(1-cos2π(V/λ)t)20.45≤t≤0.90otherwise,wfl={a2(1-cos2π(V/λ)t)20.45≤t≤0.90otherwise,wrr=wrl={a1(1-cos2π(V/λ)t)20.45+λV≤t≤0.9+λV0otherwise,
where a1 and a2 are the amplitudes of the bump, V is the horizontal speed, and λ is the wavelength of the bump. Values of the system parameters are given in Table 1.Table 1
Parameters of the full-car model.
Parameters
Value
Sprung mass (Ms)
1060 kg
Unsprung masses (muij)
40, 40, 35, 35 kg
Pitch moment of inertia (Iθ)
2200kgm2
Roll moment of inertia (Iϕ)
460 kgm2
Distance from vehicle front axle to its centre of gravity (lf)
1 m
Distance from vehicle rear axle to its centre of gravity (lr)
1.5 m
Vehicle width (tf)
1.5 m
Linear suspension stiffness at each wheel (ksijl)
2.35
×
1
0
4 N/m
Nonlinear suspension stiffness at each wheel (ksijnl)
2.35
×
1
0
6 N/m
Tyre stiffness at each wheel (ktij)
1.9
×
1
0
5 N/m
Linear suspension damping at each wheel (bsijl)
700 Ns/m
Nonlinear suspension damping at each wheel (bsijnl)
400 Ns/m
Asymmetric suspension damping at each wheel (bsijsym)
400 Ns/m
Tyre damping at each wheel (btij)
80, 70 Ns/m
Actuator parameter (αij)
4.515
×
1
0
13
Actuator parameter (βij)
1
Actuator parameter (γij)
1.545
×
1
0
9
Piston area (Aij)
3.35
×
1
0
-
4 m2
Supply pressure (Psij)
10,342,500Pa
Time constant (τij)
3.33
×
1
0
-
2 s
Servo-valve gains (Kvij)
0.001m/V
Bump amplitudes (a1, a2)
4, 6cm
Bump wavelength (λ)
5m
Vehicle speed(V)
40km/h


## 2.1. Physical Modelling

Figure2 shows the seven dof full-car suspension system used in this study. Ms, Iθ, and Iϕ are the sprung mass, pitch inertia, and roll inertia of the vehicle, respectively; lf, lr, and tf are the distance from the vehicle front axle to the centre of gravity, the distance from the vehicle rear axle to the centre of gravity, and the half width of the vehicle, respectively. z, θ, and ϕ are the vehicle body heave, pitch, and roll, respectively; the displacement and mass of each wheel are denoted as ztij and muij respectively with (i,j) representing the position of the tyre, where i refers to the lateral position front f or rear r, and j to the longitudinal position right r or left l, and the respective road disturbance profile at each wheel is given as wij. The spring and damping force contributions of the suspension system are denoted by Fkij and Fbij, respectively, and those of the tyre are denoted by Fktij and Fbtij, respectively. Faij is the actuator force supplied between the chassis and each wheel. The suspension travels are denoted by zij [4].Figure 2
Simplified full-car model of the active vehicle suspension system.


## 2.2. Mathematical Modelling

The governing equations are derived using Newton-Euler approach, with the resultant force acting in each suspension system given as(1)Fij=Fkij+Fbij-Faij.The heave, pitch, and roll dynamics of the vehicle are [4]
(2)Msz¨=Ffr+Ffl+Frr+Frl,Iθθ¨=-Ffrlf-Ffllf+Frrlr+Frllr,Iϕϕ¨=tf[Ffr-Ffl+Frr-Frl].The dynamics of each wheel are determined as follows:(3)muijz¨tij=-Fkij-Fbij+Faij+Fktij+Fbtij.The suspension travel at each wheel is determined as(4)zfl=z-lfsinθ-tfsinϕ,zrl=z+lrsinθ-tfsinϕ,zfr=z-lfsinθ+tfsinϕ,zrr=z+lfsinθ+tfsinϕ.The suspension spring and damping forces are computed as(5)Fkij=ksijl(ztij-zij)+ksijnl(ztij-zij)3,Fbij=bsijl(z˙tij-z˙ij)+bsijnl|z˙tij-z˙ij|sgn(z˙tij-z˙ij)=bsijsym|z˙tij-z˙ij|with linear and nonlinear suspension spring coefficientsksijl and ksijnl, and linear, nonlinear, and symmetric damping coefficients bsijl, bsijnl, and bsijsym, respectively. The dynamics of the tyres are modelled as
(6)Fktij=ktij(wij-ztij),Fbtij=btij(wij-ztij),
where ktij and btij are the spring and damping coefficients of the wheel. The pressure developed by the servo-hydraulic actuators is modelled as [4, 32]
(7)P˙Lij=γijΦijxvij-βijPLij+αijAijz˙ij,
where
(8)αij=4βeijVtij,βij=αijCtpij,γij=CdijSij1ρij,Φij=ϕ1ij×ϕ2ij,ϕ1ij=sgn(Psij-sgn(xvij)PLij),ϕ2ij=|Psij-sgn(xvij)PLij|,wherexvij are the spool-valves’ displacements, Psij are the supply pressures, Aij are the pistons’ cross-sectional areas, Vtij are the volumes of the cylinders, βeij are the bulk moduli of the working fluids, Φij are the hydraulic load flows, Ctpij are the leakage coefficients.C
d
i
j are the discharge coefficients between the supply lines and the hydraulic cylinders and Sij are spool-valves’ area gradients. The valves’ positions are modelled with first-order lag as follows:
(9)x˙vij=1τij(Kvijui-xvij),whereKvij are the servo-valves’ gains and τij are the time constants. The disturbance is a deterministic road bump. The vehicle approaches a sinusoidal bump at 40 km/h. The heights of the bump on the right and left sides are 4 and 6 cm, respectively, such that roll motion is induced. The road profile at each wheel is given as
(10)wfr={a1(1-cos2π(V/λ)t)20.45≤t≤0.90otherwise,wfl={a2(1-cos2π(V/λ)t)20.45≤t≤0.90otherwise,wrr=wrl={a1(1-cos2π(V/λ)t)20.45+λV≤t≤0.9+λV0otherwise,
where a1 and a2 are the amplitudes of the bump, V is the horizontal speed, and λ is the wavelength of the bump. Values of the system parameters are given in Table 1.Table 1
Parameters of the full-car model.
Parameters
Value
Sprung mass (Ms)
1060 kg
Unsprung masses (muij)
40, 40, 35, 35 kg
Pitch moment of inertia (Iθ)
2200kgm2
Roll moment of inertia (Iϕ)
460 kgm2
Distance from vehicle front axle to its centre of gravity (lf)
1 m
Distance from vehicle rear axle to its centre of gravity (lr)
1.5 m
Vehicle width (tf)
1.5 m
Linear suspension stiffness at each wheel (ksijl)
2.35
×
1
0
4 N/m
Nonlinear suspension stiffness at each wheel (ksijnl)
2.35
×
1
0
6 N/m
Tyre stiffness at each wheel (ktij)
1.9
×
1
0
5 N/m
Linear suspension damping at each wheel (bsijl)
700 Ns/m
Nonlinear suspension damping at each wheel (bsijnl)
400 Ns/m
Asymmetric suspension damping at each wheel (bsijsym)
400 Ns/m
Tyre damping at each wheel (btij)
80, 70 Ns/m
Actuator parameter (αij)
4.515
×
1
0
13
Actuator parameter (βij)
1
Actuator parameter (γij)
1.545
×
1
0
9
Piston area (Aij)
3.35
×
1
0
-
4 m2
Supply pressure (Psij)
10,342,500Pa
Time constant (τij)
3.33
×
1
0
-
2 s
Servo-valve gains (Kvij)
0.001m/V
Bump amplitudes (a1, a2)
4, 6cm
Bump wavelength (λ)
5m
Vehicle speed(V)
40km/h


## 3. System Performance Specification and Evaluation

### 3.1. Performance Specifications

The following characteristics are required of the AVSS controller in a bid to meet the set performance objectives.(1)
Nominal stability: the closed loops should be nominally stable. Stability in the inner loop is enhanced through a force feedback loop. The enhanced stability of the actuator dynamics should improve the overall system stability.(2)
Disturbance rejection: the controller should demonstrate good low frequency disturbance attenuation.(3)
Good command following: the suspension travel response of the AVSS is examined in the presence of the deterministic road inputs shown. The controller should be able to keep the steady-state error as close to zero as possible.(4)
Suspension travel is constrained to physical limits to avoid damages due to topping and bottoming. Thus it is not to exceed ±0.1 m [33].(5)
The control voltage is also limited to ±10 V.(6)
The maximum actuator force must be less than the static weight of the vehicle; that is,Fhyd<msg.(7)
For good road holding, the dynamic load that is transmitted through the road should not be larger than the static weight of the vehicle.(8)
Ride comfort: this is quantified using the vehicle body acceleration in the vertical direction. The vertical acceleration of the vehicle body needs to be minimal for a good ride comfort, especially within the low frequency band of 0.1 to 10 Hz. The peak sprung mass acceleration isz¨<4.5m/s2 [34–38].


### 3.2. Optimal Trade-Off among Performance Criteria

The objective function employed in the optimisation algorithms is presented in (11). It is designed to minimise the suspension travels, yi, actuators’ forces, Faij, sprung mass acceleration, z¨, control voltages, ui, and wheels’ dynamic loads. This way, ride comfort and road holding are improved while control voltages and actuator forces are kept as small as possible.The objective of the proposed controller is to minimise the performance index:(11)J=J1+J2+J3+J4+J5,J1=1T∫0T‍[(z¨z¨max)2+(θ¨θ¨max)2+(ϕ¨ϕ¨max)2]dt,J2=∑j=fr‍∑i=rl‍1T∫0T‍[(Fktij+Fbtij(Fktij+Fbtij)max)2]dt,J3=∑i=14‍1T∫0T‍[(yiyimax)2]dt,J4=∑j=fr‍∑i=rl‍1T∫0T‍(FaijFaijmax)2dt,J5=∑i=14‍1T∫0T‍[(uiuimax)2]dt,whereJ is the performance index, J1, J2, J3, J4, and J5 relate to the vehicle ride comfort and vehicle handling, road holding properties, suspension travel, actuation force, and power consumption, respectively, z¨max, θ¨max, ϕ¨max, (Fktij+Fbtij)max, yimax, Faijmax, and uimax are the maximum allowable heave body acceleration, pitch acceleration, roll acceleration, tyre dynamic load, suspension travels, actuation forces, and control input voltages, respectively, T is the period over which the simulation runs. The controller gains are determined either using manual tuning or through DE.


## 3.1. Performance Specifications

The following characteristics are required of the AVSS controller in a bid to meet the set performance objectives.(1)
Nominal stability: the closed loops should be nominally stable. Stability in the inner loop is enhanced through a force feedback loop. The enhanced stability of the actuator dynamics should improve the overall system stability.(2)
Disturbance rejection: the controller should demonstrate good low frequency disturbance attenuation.(3)
Good command following: the suspension travel response of the AVSS is examined in the presence of the deterministic road inputs shown. The controller should be able to keep the steady-state error as close to zero as possible.(4)
Suspension travel is constrained to physical limits to avoid damages due to topping and bottoming. Thus it is not to exceed ±0.1 m [33].(5)
The control voltage is also limited to ±10 V.(6)
The maximum actuator force must be less than the static weight of the vehicle; that is,Fhyd<msg.(7)
For good road holding, the dynamic load that is transmitted through the road should not be larger than the static weight of the vehicle.(8)
Ride comfort: this is quantified using the vehicle body acceleration in the vertical direction. The vertical acceleration of the vehicle body needs to be minimal for a good ride comfort, especially within the low frequency band of 0.1 to 10 Hz. The peak sprung mass acceleration isz¨<4.5m/s2 [34–38].


## 3.2. Optimal Trade-Off among Performance Criteria

The objective function employed in the optimisation algorithms is presented in (11). It is designed to minimise the suspension travels, yi, actuators’ forces, Faij, sprung mass acceleration, z¨, control voltages, ui, and wheels’ dynamic loads. This way, ride comfort and road holding are improved while control voltages and actuator forces are kept as small as possible.The objective of the proposed controller is to minimise the performance index:(11)J=J1+J2+J3+J4+J5,J1=1T∫0T‍[(z¨z¨max)2+(θ¨θ¨max)2+(ϕ¨ϕ¨max)2]dt,J2=∑j=fr‍∑i=rl‍1T∫0T‍[(Fktij+Fbtij(Fktij+Fbtij)max)2]dt,J3=∑i=14‍1T∫0T‍[(yiyimax)2]dt,J4=∑j=fr‍∑i=rl‍1T∫0T‍(FaijFaijmax)2dt,J5=∑i=14‍1T∫0T‍[(uiuimax)2]dt,whereJ is the performance index, J1, J2, J3, J4, and J5 relate to the vehicle ride comfort and vehicle handling, road holding properties, suspension travel, actuation force, and power consumption, respectively, z¨max, θ¨max, ϕ¨max, (Fktij+Fbtij)max, yimax, Faijmax, and uimax are the maximum allowable heave body acceleration, pitch acceleration, roll acceleration, tyre dynamic load, suspension travels, actuation forces, and control input voltages, respectively, T is the period over which the simulation runs. The controller gains are determined either using manual tuning or through DE.


## 4. Controller Design

The control system is made up of multiple control loops: four outer control loops are used to regulate the controlled variables, while the four inner control loops are used to stabilise the hydraulic actuators [2]. Figure 3 shows the proposed control architecture. The controlled variables are the suspension travel at each of the wheels and force developed by the actuator is fed back in the inner control loops.Figure 3
Control architecture of the proposed multiloop PID controller.PID controllers are described as follows:(12)uij=KPijeij(t)+KDijdeij(t)dt+KIij∫0T‍eij(t)dt,whereKPij are the proportional gains, KDij are the derivative gains, KIij are the integral gains, and the error signals eij are given as
(13)eij=ydij(t)-yij(t)=Rd(t)-yij(t),
where ydij(t) are the desired suspension travels. Thus this investigation addresses the issue of regulation, and the desired setpoints are hence set to zero.


### 4.1. Differential Evolution Optimization Method

Differential evolution is another heuristic routine and operates by firstly generating a uniformly distributed populationS={x1,x2,…,xN} to search a feasible region, where each vector denotes a set of PID gains. For each iteration, every vector in S is targeted with the objective of being altered such that its fitness improves. This alteration consists of three steps. Firstly, three distinctive individuals apart from the targeted individual, xp(1), xp(2), and xp(3), are selected randomly and are thereafter utilized to produce a mutant vector using (14) [26]:
(14)x^i=xp(1)+F(xp(2)-xp(3)),whereF>0 is a scaling factor and xp(1) is known as the base vector. If the point x^i∉Ω, then the mutation operation is repeated. The trial point yi is found from its parents xi and x^i using the following crossover rule:
(15)yij={x^ijifRj≤CRorj=IixijifRj>CRandj≠Ii,whereIi denotes a randomly selected integer, which is an element of the set I=I1,I2,…,In, superscript j signifies the jth component of the respective vectors, and Rj∈(0,1) is a vector of random numbers drawn for each j. The ultimate aim of the crossover rule (15) is to obtain the trial vector yi with components coming from the components of target vector xi and mutated vector x^i. And this is ensured by introducing CR and the set I. Notice that for CR=1 the trial vector yi is the replica of the mutated vector x^i. The effect of CR has been studied in [26, 31] and it was found that CR=0.5 is a good choice. The targeting process continues until all members of S are considered. After all N trial points yi have been generated, acceptance is applied. In the acceptance phase, the function value at the trial point, f(yi), is compared to f(xi), the value at the target point. If f(yi)<f(xi), then yi replaces xi in S; otherwise, S retains the original xi. Reproduction (mutation and crossover) and acceptance continue until some stopping conditions are met.Algorithm: The DE AlgorithmStep 1.
Determine the initial setS={x1,x2,…,xN}, where the points xi,i=1,2,…,N, are sampled randomly in Ω; evaluate f(x) at each xi,i=1,2,…,N. Set iteration counter k=0.Step 2.
Stopping condition: if the stopping condition such ask>kmax is satisfied, then stop.Step 3.
Generate points to replace points inS for the next population (or iteration). For each xi∈S(i=1,2,…,N), determine yi by the following two operations. (i)
Mutation:x^i=xp(1)+F(xp(2)-xp(3)), where xp(1), xp(2), and xp(3) are three random vectors from S and F is a scaling factor, chosen randomly. The tournament selection is applied for each i. If the jth component x^ij∉Ω, then it is generated randomly.(ii)
Crossover: calculate the trial vectoryi corresponding to the target xi from xi and x^i using the crossover rule (15).Step 4.
Acceptance rule to replace points inS: select each trial vector yi for the k+1 iteration using the acceptance criterion: replace xi∈S with yi if f(yi)<f(xi), otherwise retain, xi. Set k:=k+1 and go to Step 2.We have usedCR=0.5, F=0.75, kmax=500, and N=100. The individual x0 of the population S was chosen as an initial condition to be used by the algorithm. It had a performance index of 2.4. The evolution of the performance index from this initial condition is plotted in Figure 4, and the controller gains selected for the manually tuned and DE-optimized cases are listed in Table 2.Table 2
PID tuning parameters.
Technique
Outer PID loop gains
Inner PID loop gains
K
p
K
i
K
d
k
p
k
i
k
d
Front right suspension system
Manual
1100
360
140
0.002
0.001
0
DE
5837.262
239.436
-
123.923
4.783
×
1
0
-
3
2.563
×
1
0
-
3
3.618
×
1
0
-
9
Front left suspension system
Manual
1050
170
220
0.002
0.001
0
DE
1692.418
267.7925
166.1722
0.003855
0.001056
3.124
×
1
0
-
9
Rear right suspension system
Manual
1200
340
150
0.002
0.001
0
DE
5051.295
269.968
702.321
0.003084
0.001891
3.117
×
1
0
-
9
Rear left suspension system
Manual
1000
200
200
0.002
0.001
0
DE
3762.704
-
54.3744
-
415.94
0.004145
0.004006
1.992
×
1
0
-
9Figure 4
Convergence history of performance index through the use of differential evolution algorithm.


## 4.1. Differential Evolution Optimization Method

Differential evolution is another heuristic routine and operates by firstly generating a uniformly distributed populationS={x1,x2,…,xN} to search a feasible region, where each vector denotes a set of PID gains. For each iteration, every vector in S is targeted with the objective of being altered such that its fitness improves. This alteration consists of three steps. Firstly, three distinctive individuals apart from the targeted individual, xp(1), xp(2), and xp(3), are selected randomly and are thereafter utilized to produce a mutant vector using (14) [26]:
(14)x^i=xp(1)+F(xp(2)-xp(3)),whereF>0 is a scaling factor and xp(1) is known as the base vector. If the point x^i∉Ω, then the mutation operation is repeated. The trial point yi is found from its parents xi and x^i using the following crossover rule:
(15)yij={x^ijifRj≤CRorj=IixijifRj>CRandj≠Ii,whereIi denotes a randomly selected integer, which is an element of the set I=I1,I2,…,In, superscript j signifies the jth component of the respective vectors, and Rj∈(0,1) is a vector of random numbers drawn for each j. The ultimate aim of the crossover rule (15) is to obtain the trial vector yi with components coming from the components of target vector xi and mutated vector x^i. And this is ensured by introducing CR and the set I. Notice that for CR=1 the trial vector yi is the replica of the mutated vector x^i. The effect of CR has been studied in [26, 31] and it was found that CR=0.5 is a good choice. The targeting process continues until all members of S are considered. After all N trial points yi have been generated, acceptance is applied. In the acceptance phase, the function value at the trial point, f(yi), is compared to f(xi), the value at the target point. If f(yi)<f(xi), then yi replaces xi in S; otherwise, S retains the original xi. Reproduction (mutation and crossover) and acceptance continue until some stopping conditions are met.Algorithm: The DE AlgorithmStep 1.
Determine the initial setS={x1,x2,…,xN}, where the points xi,i=1,2,…,N, are sampled randomly in Ω; evaluate f(x) at each xi,i=1,2,…,N. Set iteration counter k=0.Step 2.
Stopping condition: if the stopping condition such ask>kmax is satisfied, then stop.Step 3.
Generate points to replace points inS for the next population (or iteration). For each xi∈S(i=1,2,…,N), determine yi by the following two operations. (i)
Mutation:x^i=xp(1)+F(xp(2)-xp(3)), where xp(1), xp(2), and xp(3) are three random vectors from S and F is a scaling factor, chosen randomly. The tournament selection is applied for each i. If the jth component x^ij∉Ω, then it is generated randomly.(ii)
Crossover: calculate the trial vectoryi corresponding to the target xi from xi and x^i using the crossover rule (15).Step 4.
Acceptance rule to replace points inS: select each trial vector yi for the k+1 iteration using the acceptance criterion: replace xi∈S with yi if f(yi)<f(xi), otherwise retain, xi. Set k:=k+1 and go to Step 2.We have usedCR=0.5, F=0.75, kmax=500, and N=100. The individual x0 of the population S was chosen as an initial condition to be used by the algorithm. It had a performance index of 2.4. The evolution of the performance index from this initial condition is plotted in Figure 4, and the controller gains selected for the manually tuned and DE-optimized cases are listed in Table 2.Table 2
PID tuning parameters.
Technique
Outer PID loop gains
Inner PID loop gains
K
p
K
i
K
d
k
p
k
i
k
d
Front right suspension system
Manual
1100
360
140
0.002
0.001
0
DE
5837.262
239.436
-
123.923
4.783
×
1
0
-
3
2.563
×
1
0
-
3
3.618
×
1
0
-
9
Front left suspension system
Manual
1050
170
220
0.002
0.001
0
DE
1692.418
267.7925
166.1722
0.003855
0.001056
3.124
×
1
0
-
9
Rear right suspension system
Manual
1200
340
150
0.002
0.001
0
DE
5051.295
269.968
702.321
0.003084
0.001891
3.117
×
1
0
-
9
Rear left suspension system
Manual
1000
200
200
0.002
0.001
0
DE
3762.704
-
54.3744
-
415.94
0.004145
0.004006
1.992
×
1
0
-
9Figure 4
Convergence history of performance index through the use of differential evolution algorithm.


## 5. Simulation Results and Discussion

This paper investigates a novel tuning method, which uses DE optimization to compute gains with the objective of firstly improving the overall performance index of an AVSS and secondly finding the best compromise between these conflicting performance criteria. Simulations were executed in MATLAB/Simulink. Thus far, it was shown in Figure4 that the DE-optimized case was able to minimize the performance index, but this information alone cannot guarantee that the trade-offs between the various performance criterion have been resolved. Thus, it is still imperative that the AVSS system response in each performance criteria must be analysed to test the effectiveness of DE optimization algorithm in performing its objectives.The conflicting performance criteria results reported include the suspension travel at each wheel, road holding which is directly related to the force experienced by each of the tyres, ride comfort which is essentially analysed through the body-heave acceleration, vehicle handling which primarily is a function of roll and pitch accelerations, and power consumption which is the voltage utilised at each actuator and the resulting hydraulic force, each of which are plotted in Figures5–11, respectively. The suspension travel, road holding, and power consumption at the rear left wheel were only considered as it possessed the worst response in each of these facets. The root-mean-square (RMS) and peak values of all suspension performance criteria are listed in Table 3.Table 3
RMS and peak values of all the suspension performance criteria.
Technique
RMS
Peak
RMS
Peak
Front right suspension travel (m)
Front left suspension travel (m)
Passive
0.01313
0.04816
0.01079
0.03815
Manual
0.00915
0.03922
0.00667
0.02942
DE
0.00943
0.03915
0.00676
0.02516
Rear right suspension travel (m)
Rear Left suspension travel (m)
Passive
0.00798
0.03398
0.00766
0.03367
Manual
0.00866
0.03483
0.00557
0.02593
DE
0.00678
0.02863
0.00545
0.024698
Front right tyre dynamic load (N)
Front Left tyre dynamic load (N)
Passive
395.74
1594.05
321.46
1254.53
Manual
321.18
1531.90
263.39
1295.15
DE
189.83
984.83
132.095
713.38
Rear right tyre dynamic load (N)
Rear left tyre dynamic load (N)
Passive
248.95
1194.39
240.81
1158.58
Manual
248.77
1234.10
218.73
1096.45
DE
229.63
1196.07
185.24
930.12
Front right control input voltage (V)
Front left control input voltage (V)
Passive
N/A
N/A
N/A
N/A
Manual
0.3166
1.4946
0.2256
1.0739
DE
0.3189
1.6125
0.2255
1.1875
Rear right control input voltage (V)
Rear left control input voltage (V)
Passive
N/A
N/A
N/A
N/A
Manual
0.2659
1.3492
0.2129
1.0484
DE
0.2786
1.3710
0.2044
1.0156
Pitch acceleration (rad/s2)
Roll acceleration (rad/s2)
Passive
0.5163
2.3672
0.1308
0.5885
Manual
0.4631
2.2687
0.1169
0.5388
DE
0.3140
1.6223
0.06685
0.3289
Heave acceleration (m/s2)
Cumulative actuation force (N)
Passive
0.6488
2.4473
N/A
N/A
Manual
0.4298
1.5946
498.569
1849.75
DE
0.3070
1.7250
342.1241
1363.52Figure 5
Time response of the rear left suspension travel for the passive and the proposed tuning methods.Figure 6
Time response of the rear left tyre dynamic load for the passive and the proposed tuning methods.Figure 7
Time response of the rear left control input voltage for the proposed tuning methods.Figure 8
Time response of the pitch acceleration for the passive and the proposed tuning methods.Figure 9
Time response of the roll acceleration for the passive and the proposed tuning methods.Figure 10
Time response of the body-heave acceleration for the passive and the proposed tuning methods.Figure 11
Time response of the electrohydraulic actuator force for the proposed tuning methods.The suspension travel response shown in Figure5 indicated a considerable improvement in terms of peak and RMS values when using control. Furthermore, there was an improvement in transient behaviour where the controlled cases had one less oscillation. DE-optimized PID controller tuning produced marginally better peaks and RMS values than the manually tuned case and had a slightly quicker settling time (see Table 3). Such response of the DE-based PID-controlled case was expected as suspension travel was a major factor in the cost function of the algorithm.Figures6, 8, 9, and 10 showed a similar trend in the road holding, pitch, roll, and body-heave accelerations, respectively. This is for the very same reason that these aforementioned responses significantly affected the performance index. Moreover, in these criteria, the variations in peak and RMS values between the DE-optimized and manually tuned PID-controlled cases, as shown in Table 3, were much greater. This was because the weighting factors of these criteria with the performance index were slightly larger than that of the suspension travel.On the other hand, the control input voltages of the two tuning methods were fairly similar (see Figure7), whereas the control force of the DE-optimized case was considerably better in both peak and RMS values and had a transient behaviour (see Table 3 and Figure 11). In reality, the direct opposite would have been anticipated as a larger actuation force is often required to minimise the suspension travel, road holding, and vehicle handling which was previously reported. However, the inherent coupling as well as nonlinearities in this full-car system plays a significant role in altering such relationships. The fact that the DE-optimized case was able to produce such a solution infers that it is fully capable of dealing with nonlinearities and coupling.However, the controlled cases did possess a drawback in terms of chattering, with the DE-optimized case being more severe. Such a response could be due to the marginally high outer control loop gains used computed by the DE optimization algorithm. Such behaviour is undesirable as it will lead to component degradation.


### 5.1. Sensitivity to Parameter Variations Analysis

Each of the results plotted possessed a steady-state error that was in the order of magnitude of 0.001% of the peak values. Furthermore, a stability study was conducted to test the control systems sensitivity to parameter variations. The parameter variations apply to the mass and inertia of the vehicle, as the chassis will fluctuate in mass due to variations in fuel and passengers, tyre damping and stiffness, as the wheels will experience changes in pressure, and vehicle speed. Such parameters are expected to change within a range of ±20%, and hence a sensitivity investigation was conducted along these lines. Since the DE-optimized case produced the most favourable results, stability is conducted on it. Figures12 and 13 show the suspension travel response at the rear left wheel for −20% and +20% variations in the selected parameters, respectively.Figure 12
Suspension travel response of the DE-based PID-controlled case for a −20% variation in selected parameters.Figure 13
Suspension travel response of the DE-based PID-controlled case for a +20% variation in selected parameters.Figures12 and 13 show that the control system remains bounded-input-bounded-output (BIBO) stable within an acceptable steady-state error for all anticipated parameter variations. Apart from the +20% variation in vehicle speed, the sensitivity was excellent with good transient behaviour and peak values. Hence, it may be concluded that the DE-based PID controller has an acceptable sensitivity to parameter variations except for a 20% increase in vehicle speed.


### 5.2. Frequency Domain Analysis

This section presents the pseudo-frequency domain analysis (based on power spectrum density estimation) for the vehicle ride comfort and handling characteristics. While vehicle ride comfort is commonly associated with the level of vibrations perceived by the occupants of the vehicle, it is also necessary to minimise the pitch and roll accelerations to obtain good handling. Analysis for ride comfort is normally carried out in consonance with ISO 2631 [36, 38].Figures14–16 present a pseudo-frequency domain analysis carried out using the power-spectral-density-(PSD-) estimates based Welch algorithm in the MATLAB/Simulink signal processing toolbox. The Welch periodograms were computed using the following parameters: the windowing function—Hanning window function, the number of points used in forming each fast Fourier transform, NFFT = 1024, and length of the window, NWind = 256; and the sampling frequency of the windows was set at 80 Hz to accommodate the whole-body vibration (WBV) range.Figure 14
Frequency domain analysis of the body-heave acceleration.Figure 15
Frequency domain analysis of the pitch acceleration.Figure 16
Frequency domain analysis of the roll acceleration.Figure14 presents the frequency weighted vehicle body-heave acceleration mode. This mode of acceleration is most critical in the estimation of vehicle ride comfort. The analysis presented covers 0.5–80 Hz, which is the WBV frequency range. This range includes the critical frequency ranges that need to be avoided for both human and vehicle comfort. Apart from the PID-controlled AVSS whose magnitude reached 1dB at the onset of the analysis, all the signals showed attenuation until a sharp resonance occurred around 40 Hz. The attenuation trend stopped around 80 Hz.Figures15 and 16 show similar trend for the pitch and roll acceleration modes. The signals for the PID-controlled AVSS and the PVSS coincided all through the range of the analysis. All the signals showed attenuation of the signals at frequencies below 80 Hz, except at 40 Hz where there was a sharp spike of resonance.


## 5.1. Sensitivity to Parameter Variations Analysis

Each of the results plotted possessed a steady-state error that was in the order of magnitude of 0.001% of the peak values. Furthermore, a stability study was conducted to test the control systems sensitivity to parameter variations. The parameter variations apply to the mass and inertia of the vehicle, as the chassis will fluctuate in mass due to variations in fuel and passengers, tyre damping and stiffness, as the wheels will experience changes in pressure, and vehicle speed. Such parameters are expected to change within a range of ±20%, and hence a sensitivity investigation was conducted along these lines. Since the DE-optimized case produced the most favourable results, stability is conducted on it. Figures12 and 13 show the suspension travel response at the rear left wheel for −20% and +20% variations in the selected parameters, respectively.Figure 12
Suspension travel response of the DE-based PID-controlled case for a −20% variation in selected parameters.Figure 13
Suspension travel response of the DE-based PID-controlled case for a +20% variation in selected parameters.Figures12 and 13 show that the control system remains bounded-input-bounded-output (BIBO) stable within an acceptable steady-state error for all anticipated parameter variations. Apart from the +20% variation in vehicle speed, the sensitivity was excellent with good transient behaviour and peak values. Hence, it may be concluded that the DE-based PID controller has an acceptable sensitivity to parameter variations except for a 20% increase in vehicle speed.


## 5.2. Frequency Domain Analysis

This section presents the pseudo-frequency domain analysis (based on power spectrum density estimation) for the vehicle ride comfort and handling characteristics. While vehicle ride comfort is commonly associated with the level of vibrations perceived by the occupants of the vehicle, it is also necessary to minimise the pitch and roll accelerations to obtain good handling. Analysis for ride comfort is normally carried out in consonance with ISO 2631 [36, 38].Figures14–16 present a pseudo-frequency domain analysis carried out using the power-spectral-density-(PSD-) estimates based Welch algorithm in the MATLAB/Simulink signal processing toolbox. The Welch periodograms were computed using the following parameters: the windowing function—Hanning window function, the number of points used in forming each fast Fourier transform, NFFT = 1024, and length of the window, NWind = 256; and the sampling frequency of the windows was set at 80 Hz to accommodate the whole-body vibration (WBV) range.Figure 14
Frequency domain analysis of the body-heave acceleration.Figure 15
Frequency domain analysis of the pitch acceleration.Figure 16
Frequency domain analysis of the roll acceleration.Figure14 presents the frequency weighted vehicle body-heave acceleration mode. This mode of acceleration is most critical in the estimation of vehicle ride comfort. The analysis presented covers 0.5–80 Hz, which is the WBV frequency range. This range includes the critical frequency ranges that need to be avoided for both human and vehicle comfort. Apart from the PID-controlled AVSS whose magnitude reached 1dB at the onset of the analysis, all the signals showed attenuation until a sharp resonance occurred around 40 Hz. The attenuation trend stopped around 80 Hz.Figures15 and 16 show similar trend for the pitch and roll acceleration modes. The signals for the PID-controlled AVSS and the PVSS coincided all through the range of the analysis. All the signals showed attenuation of the signals at frequencies below 80 Hz, except at 40 Hz where there was a sharp spike of resonance.


## 6. Conclusion

DE-optimized multi-loop PID-based control involving force feedback produced significant improvements in comparison to the PVSS and manually tuned cases. The recommended DE optimization was successful both in improving the performance index of the AVSS and in finding a better compromise between the conflicting AVSS design requirements. It also attained superior RMS and peak values in suspension travel, road holding, vehicle handling, and cumulative actuation force as compared to the manually tuned PID-controlled case. The observed reduction in actuation force is counterintuitive as larger forces are often required to minimise the other performance criteria. However, when dealing with nonlinear and coupled systems such as the full-car AVSS, this is possible, which infers that the DE-based optimization algorithm was able to account for coupling and nonlinearities as well. Chattering was the only shortfall of the controlled systems as it would degrade system components.



---

*Source: 261582-2013-09-17.xml*