# Commercial Vehicle Ride Comfort Optimization Based on Intelligent Algorithms and Nonlinear Damping

Authors: Shuilong He; Keren Chen; Enyong Xu; Wei Wang; Zhansi Jiang
Journal: Shock and Vibration
 (2019)

---

## Abstract

The method chosen to conduct vehicle dynamic modeling has a significant impact on the evaluation and optimization of ride comfort. This paper summarizes the current modeling methods of ride comfort and their limitations. Then, models based on nonlinear damping and equivalent damping and the multibody dynamic model are developed and simulated in Matlab/Simulink and Adams/Car. The driver seat responses from these models are compared, showing that the accuracy of the ride comfort model based on nonlinear damping is higher than the one based on equivalent damping. To improve the reliability of ride comfort optimization and analysis, a ride comfort optimization method based on nonlinear damping and intelligent algorithms is proposed. The sum of the frequency-weighted RMS of the driver seat acceleration, the RMS of dynamic tyre load, and suspension working space is taken as the objective function in this article, using nonlinear damping coefficients and stiffness of suspension as design variables. By applying the particle swarm optimization (PSO), cuckoo search (CS), dividing rectangles (DIRECT), and genetic algorithm (GA), a set of optimal solutions are obtained. The method efficiency is verified through a comparison between frequency-weighted RMS before and after optimization. Results show that the frequency-weighted RMS of driver seat acceleration, RMS values of the suspension working space of the front and rear axles, and RMS values of the dynamic tyre load of front and rear wheels are decreased by an average of 27.4%, 21.6%, 25.0%, 19.3%, and 22.3%, respectively. The developed model is studied in a pilot commercial vehicle, and the results show that the optimization method proposed in this paper is more practical and features improvement over previous models.


---

## Body

## 1. Introduction

Commercial vehicles are an essential part of the modern transport network and are responsible for the bulk of freight transport around the world. The demand for better ride comfort and safety in these vehicles is increasing, as such commercial vehicle ride comfort optimization has been an active area of research. Ride comfort in commercial vehicles provides the following advantages: (1) increased driver’s comfort and driving safety, ensuring good ride ability, and reducing the incidence of traffic accidents; (2) assurance goods arrive in better condition, increasing the utility value of commercial vehicles; and (3) improved service life of commercial vehicle parts, as vehicles are subjected to reduced impact and vibration forces constantly [1–4].Typically, software used for vehicle ride comfort simulation includes Adams and Matlab. The former adopts the structure-oriented modeling method, which is more accurate, and the simulation is reliable. However, the modeling procedures need to set characteristic parameters and curves for a large number of parts, resulting in a complex and time-consuming modeling solution. The latter establishes a ride comfort model from the perspective of mechanics and mathematics with lower accuracy but simpler modeling, a shorter simulation time, and convenient parameter setting and adjustment, with a powerful ability to perform algorithm-based optimization.Many studies have been devoted to vibration control and ride comfort optimization of vehicles. Ding et al. [5] established an eleven degree-of-freedom (DOF) ride comfort model for a three-axle truck that assumes the body is slightly vibrating near the equilibrium position. Li et al. [6] studied the synthetic mechanical model of plane rigid-frame heavy-duty vehicles, considering only the microvibrations caused by road excitation. A dynamic model is developed that aims to reduce the road excitation transmitted to the vehicle body and improve the ride comfort by using optimization algorithms. The spring stiffness and damping coefficient are regarded as a fixed value [7–9]. The ride model is used to study the relationship between ride comfort and the in-wheel motor. There is a linear relationship between damping and damping force [10, 11]. Zhu et al. [12] established a ten-DOF ride model for analyzing vehicle ride comfort and driving safety, and the damping force is considered a linear parameter in this analysis.Most of the research on simulation and optimization for vehicle ride comfort make the following assumptions on the actual situation: (1) the microvibration caused by road excitation is only considered and (2) the damping force is a linear function of its velocity and is regarded a constant value (equivalent damping). The research on ride comfort mainly focuses on equivalent damping when the model consists of passive suspension, with work based on nonlinear damping poorly established. However, in terms of the analysis for the ride comfort of commercial vehicles, the above assumptions are different from the actual situation. Commercial vehicle suspension damping is mostly a nonlinear damping curve that varies with extension or compression velocity. In addition, the driving conditions of commercial vehicles is far from ideal, with additional factors contributing to this, which goes against the assumption that only microvibration should be considered. Furthermore, this work is mainly focused on ride comfort and disregards the suspension working space and dynamic tyre load, which affect handling stability. The reduction of suspension stiffness has a negative effect on the suspension working space and dynamic tyre load.In this paper, the ride comfort models of a commercial vehicle based on nonlinear damping and equivalent damping are developed. The responses of theoretical models are verified by simulation performed in Adams/Car on a multibody dynamic model of the target vehicle. To improve ride comfort, more performance criteria are considered and optimized. To achieve this, a ride comfort optimization method based on nonlinear damping of suspension and intelligent algorithms is proposed. The objective functions are optimized by algorithms, and the optimal design variables that consider nonlinear damping are obtained. Lastly, the theoretical model and optimization are validated by experiment. In short, the method is more practical and accurate for vehicle ride comfort optimization compared to those based on equivalent damping.


## 2. Ride Comfort Model

In this section, the equations are serially developed for the vibration system, road excitation model, nonlinear damping, and equivalent damping model. The effectiveness of the ride comfort model using nonlinear damping is demonstrated and compared with equivalent damping and multibody dynamic models.


### 2.1. Vibration Model

Generally, there are three types of vehicle ride comfort models, including the 1/4 vehicle model, half vehicle model, and full vehicle model. When a vehicle is symmetrical to its longitudinal axis, the main vibrations affecting ride comfort are vertical vibration along thez-axis and pitch vibration around the y-axis of the vehicle coordinate system. Therefore, a half commercial vehicle model that consists of five basic subsystems of the tires, suspension systems, frame, cab, and seat is selected in this paper, as shown in Figure 1. For this model, the vibration differential equation (1) is generated by using Newton’s second law of motion:where xf and xr are the displacement along the z-axis of the front wheel and the rear wheel, respectively; xfu, xru, xb, xc, and xp represent the displacement of the front axle unsprung mass, rear axle unsprung mass, sprung mass, mass of cab, and mass of the seat, respectively; θc, θb, and θ are the pitch angle displacements of the cab, sprung mass, and rear axle balance suspension around the y-axis, respectively. Moreover, the fixed parameters are presented in Table 1, and the nonlinear damping curves are shown in Figure 2(b).Figure 1
Seven-DOF vibration model.(1)
m
p
x
¨
p
=
−
k
p
x
p
−
x
c
−
l
1
θ
c
−
c
p
x
˙
p
−
x
˙
c
−
l
1
θ
˙
c
,
m
c
x
¨
c
=
k
p
x
p
−
x
c
−
l
1
θ
c
+
c
p
x
˙
p
−
x
˙
c
−
l
1
θ
˙
c
−
k
fc
x
c
−
x
b
−
l
2
θ
c
+
l
2
+
l
4
θ
b
−
c
fc
x
˙
c
−
x
˙
b
−
l
2
θ
˙
c
+
l
2
+
l
4
θ
˙
b
−
k
rc
x
c
−
x
b
+
l
3
θ
c
+
l
4
−
l
3
θ
b
−
c
rc
x
˙
c
−
x
˙
b
+
l
3
θ
˙
c
+
l
4
−
l
3
θ
˙
b
,
J
c
θ
¨
c
=
k
p
x
p
−
x
c
−
l
1
θ
c
+
c
p
x
˙
p
−
x
˙
c
−
l
1
θ
˙
c
l
1
+
k
fc
x
c
−
x
b
−
l
2
θ
c
+
l
2
+
l
4
θ
b
+
c
fc
x
˙
c
−
x
˙
b
−
l
2
θ
˙
c
+
l
2
+
l
4
θ
˙
b
l
2
−
k
rc
x
c
−
x
b
+
l
3
θ
c
+
l
4
−
l
3
θ
b
+
c
rc
x
˙
c
−
x
˙
b
+
l
3
θ
˙
c
+
l
4
−
l
3
θ
˙
b
l
3
,
m
b
x
¨
b
=
k
fc
x
c
−
x
b
−
l
2
θ
c
+
l
2
+
l
4
θ
b
+
c
fc
x
˙
c
−
x
˙
b
−
l
2
θ
˙
c
+
l
2
+
l
4
θ
˙
b
−
k
fs
x
b
−
x
fu
−
l
5
θ
b
−
c
fs
x
˙
b
−
x
˙
fu
−
l
5
θ
˙
b
+
k
rc
x
c
−
x
b
+
l
3
θ
c
+
l
4
−
l
3
θ
b
+
c
rc
x
˙
c
−
x
˙
b
+
l
3
θ
˙
c
+
l
4
−
l
3
θ
˙
b
−
k
rs
x
b
−
x
ru
+
l
6
θ
b
−
c
rs
x
˙
b
−
x
˙
ru
+
l
6
θ
˙
b
,
J
b
θ
¨
b
=
k
fs
x
b
−
x
fu
−
l
5
θ
b
+
c
fs
x
˙
b
−
x
˙
fu
−
l
5
θ
˙
b
l
5
−
k
rs
x
b
−
x
ru
+
l
6
θ
b
+
c
rs
x
˙
b
−
x
˙
ru
+
l
6
θ
˙
b
l
6
−
k
fc
x
c
−
x
b
−
l
2
θ
c
+
l
2
+
l
4
θ
b
+
c
fc
x
˙
c
−
x
˙
b
−
l
2
θ
˙
c
+
l
2
+
l
4
θ
˙
b
l
2
+
l
4
−
k
rc
x
c
−
x
b
+
l
3
θ
c
+
l
4
−
l
3
θ
b
+
c
rc
x
˙
c
−
x
˙
b
+
l
3
θ
˙
c
+
l
4
−
l
3
θ
˙
b
l
4
−
l
3
,
m
fu
x
¨
fu
=
k
fs
x
b
−
x
fu
−
l
5
θ
b
+
c
fs
x
˙
b
−
x
˙
fu
−
l
5
θ
˙
b
−
k
ft
x
fu
−
x
f
,
m
ru
x
¨
ru
=
k
rs
x
b
−
x
ru
+
l
6
θ
b
+
c
rs
x
˙
b
−
x
˙
ru
+
l
6
θ
˙
b
−
k
rt
x
ru
−
x
r
,Table 1
Fixed parameters of the vibration model.
Symbol
Description
Value
Unit
m
p
Driver and seat mass
100
kg
m
c
Cab mass
850
kg
m
b
Sprung mass
14238
kg
m
fu
Unsprung mass of front axle
607
kg
m
ru
Unsprung mass of rear axle
1054
kg
J
c
Rotational inertia of cab mass aroundy-axis
560
kg·m2
J
b
Rotational inertia of sprung mass aroundy-axis
115000
kg·m2
J
Rotational inertia of rear axle balance suspension mass aroundy-axis
615
kg·m2
k
p
Seat stiffness
20000
N/m
k
fc
Front suspended stiffness
25000
N/m
k
rc
Rear suspended stiffness
25000
N/m
k
fs
Front suspension stiffness
410000
N/m
k
rs
Rear suspension stiffness
1476000
N/m
k
ft
Front wheel stiffness
1800000
N/m
k
rt
Rear wheel stiffness
3600000
N/m
c
p
Equivalent damping coefficient of driver seat
800
N·s/m
c
fs
Equivalent damping coefficient of front suspended
5000
N·s/m
c
rs
Equivalent damping coefficient of rear suspended
5000
N·s/m
c
ft
Equivalent damping coefficient of front suspension
12000
N·s/m
c
rt
Equivalent damping coefficient of rear suspension
15000
N·s/m
l
1
Distance between the center and seat of the cab
0.2
m
l
2
Distance between front suspended and the center of the cab
0.76
m
l
3
Distance between rear suspended and the center of the cab
0.87
m
l
4
Distance between the center of the cab and the center of target commercial vehicle
2.47
m
l
5
Distance between the front suspension and the center of target commercial vehicle
2.1
m
l
6
Distance between the rear suspension and the center of target commercial vehicle
3
mFigure 2
Damping curves. (a) Damping ratio. (b) Damping coefficient. (c) Damping force.
(a)
(b)
(c)


### 2.2. Road Excitation Model

As the most important excitation source when driving, the acquisition of road excitation is highly significant for ride comfort analysis. The road level and vehicle speed are primary factors affecting vehicle vibration. Power spectral density (PSD) is usually used for statistical characteristics of road excitation [13, 14]. The PSD of road excitation satisfies(2)Gqn=Gqn0nn0−w,where n is the spatial frequency; n0 is the reference spatial frequency, n0 = 0.1 m−1; Gqn0 is the road surface PSD corresponding to the reference spatial frequency n0, or the road excitation coefficient; and w is the frequency index that determines the frequency structure of the road surface PSD.The time-domain model of road excitation based on filtered white noise is as follows:(3)q˙t=−2πunqqt+2πn0Gqn0uωt,where u is the vehicle speed (m/s); nq is the cutoff space frequency, nq = 0.011 m−1; qt is the displacement excitation of the road surface; n0 represents the spatial reference frequency and is equal to 0.1 m−1; Gqn0 is the road excitation coefficient, Gqn0 = 256 × 10−6 m3; and ωt is the Gauss white noise with zero mean.Figure3 shows a road profile when the vehicle travels on the C-level road at 60 km/h. The profile PSD and standard road PSD according to ISO 8608 are presented in Figure 4 and show that a simulation model of road excitation is feasible.Figure 3
Road excitation time-domain profile.Figure 4
Road excitation frequency-domain profile.


### 2.3. Nonlinear Damping Model

To capture the damping nonlinear characteristic, a set of damping ratio curves are provided by a commercial vehicle manufacturer and shown in Figure2(a). The damping force is calculated by the following equation:(4)F=cv=2ξkmv,where F and c are the damping force and damping coefficient, respectively; ξ represents the damping ratio; k represents the stiffness; m refers to the sprung mass; and v is the compression or extension velocity.The damping force lines presented in Figure2(c) can be expressed as follows:(5)F=c1+η⋅sgnvvn,where η is the asymmetry coefficient and indicates the difference between the compression and extension coefficient and n is the damping characteristic index.


### 2.4. Equivalent Damping Model

Usually, the equivalent damping coefficient is calculated through the energy dissipated during one vibration period, which is equal to the energy dissipated by the equivalent damping. When a system is forced to vibrate, the viscous damping force can be expressed as the following equation:(6)x=Xsinωt−φ,x˙=Xωcosωt−φ,Fe=cx˙=cXωcosωt−φ.Energy consumption by equivalent damping in a vibration cycle is(7)We=∫0TFex˙dt=∫02π/ωceX2ω2cos2ωt−φ=πceωX2.When the damping of the system is nonlinear, it can be replaced by equivalent damping:(8)πceωX2=Wd,ce=WdωX2π,where Wd is the energy consumed by the equivalent damping during one vibration period and ce is the equivalent damping coefficient.


### 2.5. Verification of Ride Comfort Model

Equations (1), (3), and (5) have been solved using Matlab/Simulink and validated through a multibody dynamic model in Adams/Car. The multibody dynamic model is shown in Figure 5. The response of the driver seat, obtained from the Simulink models based on nonlinear damping, equivalent damping, and multibody dynamics model, is compared through time-domain curves and frequency-weighted RMS to verify the accuracy of the nonlinear damping ride model.Figure 5
Multibody dynamic model.The time-domain curves are shown in Figure6. Overall, the Matlab-based model has the same level of acceleration amplitude as the Adams-based model. The acceleration amplitude of the model based on nonlinear damping is greater than the one based on equivalent damping, leading to a higher frequency-weighted RMS, and the situation is consistent with the RMS in Figure 7. The difference in the curves is typically caused by the different methods of road modeling in Adams and Matlab and the bushings in the Adams model.Figure 6
Time-domain responses of the driver seat. (a) 40 km/h. (b) 80 km/h.
(a)
(b)Figure 7
Frequency-weighted RMS.As shown in Figure7, there is a similar tendency that ride comfort increases with speed among the three models. The results show that simulation results of the model based on the nonlinear damping are more consistent with the Adams model than the model based on equivalent damping. Therefore, the ride comfort model based on nonlinear damping has improved accuracy compared to equivalent damping models.


## 2.1. Vibration Model

Generally, there are three types of vehicle ride comfort models, including the 1/4 vehicle model, half vehicle model, and full vehicle model. When a vehicle is symmetrical to its longitudinal axis, the main vibrations affecting ride comfort are vertical vibration along thez-axis and pitch vibration around the y-axis of the vehicle coordinate system. Therefore, a half commercial vehicle model that consists of five basic subsystems of the tires, suspension systems, frame, cab, and seat is selected in this paper, as shown in Figure 1. For this model, the vibration differential equation (1) is generated by using Newton’s second law of motion:where xf and xr are the displacement along the z-axis of the front wheel and the rear wheel, respectively; xfu, xru, xb, xc, and xp represent the displacement of the front axle unsprung mass, rear axle unsprung mass, sprung mass, mass of cab, and mass of the seat, respectively; θc, θb, and θ are the pitch angle displacements of the cab, sprung mass, and rear axle balance suspension around the y-axis, respectively. Moreover, the fixed parameters are presented in Table 1, and the nonlinear damping curves are shown in Figure 2(b).Figure 1
Seven-DOF vibration model.(1)
m
p
x
¨
p
=
−
k
p
x
p
−
x
c
−
l
1
θ
c
−
c
p
x
˙
p
−
x
˙
c
−
l
1
θ
˙
c
,
m
c
x
¨
c
=
k
p
x
p
−
x
c
−
l
1
θ
c
+
c
p
x
˙
p
−
x
˙
c
−
l
1
θ
˙
c
−
k
fc
x
c
−
x
b
−
l
2
θ
c
+
l
2
+
l
4
θ
b
−
c
fc
x
˙
c
−
x
˙
b
−
l
2
θ
˙
c
+
l
2
+
l
4
θ
˙
b
−
k
rc
x
c
−
x
b
+
l
3
θ
c
+
l
4
−
l
3
θ
b
−
c
rc
x
˙
c
−
x
˙
b
+
l
3
θ
˙
c
+
l
4
−
l
3
θ
˙
b
,
J
c
θ
¨
c
=
k
p
x
p
−
x
c
−
l
1
θ
c
+
c
p
x
˙
p
−
x
˙
c
−
l
1
θ
˙
c
l
1
+
k
fc
x
c
−
x
b
−
l
2
θ
c
+
l
2
+
l
4
θ
b
+
c
fc
x
˙
c
−
x
˙
b
−
l
2
θ
˙
c
+
l
2
+
l
4
θ
˙
b
l
2
−
k
rc
x
c
−
x
b
+
l
3
θ
c
+
l
4
−
l
3
θ
b
+
c
rc
x
˙
c
−
x
˙
b
+
l
3
θ
˙
c
+
l
4
−
l
3
θ
˙
b
l
3
,
m
b
x
¨
b
=
k
fc
x
c
−
x
b
−
l
2
θ
c
+
l
2
+
l
4
θ
b
+
c
fc
x
˙
c
−
x
˙
b
−
l
2
θ
˙
c
+
l
2
+
l
4
θ
˙
b
−
k
fs
x
b
−
x
fu
−
l
5
θ
b
−
c
fs
x
˙
b
−
x
˙
fu
−
l
5
θ
˙
b
+
k
rc
x
c
−
x
b
+
l
3
θ
c
+
l
4
−
l
3
θ
b
+
c
rc
x
˙
c
−
x
˙
b
+
l
3
θ
˙
c
+
l
4
−
l
3
θ
˙
b
−
k
rs
x
b
−
x
ru
+
l
6
θ
b
−
c
rs
x
˙
b
−
x
˙
ru
+
l
6
θ
˙
b
,
J
b
θ
¨
b
=
k
fs
x
b
−
x
fu
−
l
5
θ
b
+
c
fs
x
˙
b
−
x
˙
fu
−
l
5
θ
˙
b
l
5
−
k
rs
x
b
−
x
ru
+
l
6
θ
b
+
c
rs
x
˙
b
−
x
˙
ru
+
l
6
θ
˙
b
l
6
−
k
fc
x
c
−
x
b
−
l
2
θ
c
+
l
2
+
l
4
θ
b
+
c
fc
x
˙
c
−
x
˙
b
−
l
2
θ
˙
c
+
l
2
+
l
4
θ
˙
b
l
2
+
l
4
−
k
rc
x
c
−
x
b
+
l
3
θ
c
+
l
4
−
l
3
θ
b
+
c
rc
x
˙
c
−
x
˙
b
+
l
3
θ
˙
c
+
l
4
−
l
3
θ
˙
b
l
4
−
l
3
,
m
fu
x
¨
fu
=
k
fs
x
b
−
x
fu
−
l
5
θ
b
+
c
fs
x
˙
b
−
x
˙
fu
−
l
5
θ
˙
b
−
k
ft
x
fu
−
x
f
,
m
ru
x
¨
ru
=
k
rs
x
b
−
x
ru
+
l
6
θ
b
+
c
rs
x
˙
b
−
x
˙
ru
+
l
6
θ
˙
b
−
k
rt
x
ru
−
x
r
,Table 1
Fixed parameters of the vibration model.
Symbol
Description
Value
Unit
m
p
Driver and seat mass
100
kg
m
c
Cab mass
850
kg
m
b
Sprung mass
14238
kg
m
fu
Unsprung mass of front axle
607
kg
m
ru
Unsprung mass of rear axle
1054
kg
J
c
Rotational inertia of cab mass aroundy-axis
560
kg·m2
J
b
Rotational inertia of sprung mass aroundy-axis
115000
kg·m2
J
Rotational inertia of rear axle balance suspension mass aroundy-axis
615
kg·m2
k
p
Seat stiffness
20000
N/m
k
fc
Front suspended stiffness
25000
N/m
k
rc
Rear suspended stiffness
25000
N/m
k
fs
Front suspension stiffness
410000
N/m
k
rs
Rear suspension stiffness
1476000
N/m
k
ft
Front wheel stiffness
1800000
N/m
k
rt
Rear wheel stiffness
3600000
N/m
c
p
Equivalent damping coefficient of driver seat
800
N·s/m
c
fs
Equivalent damping coefficient of front suspended
5000
N·s/m
c
rs
Equivalent damping coefficient of rear suspended
5000
N·s/m
c
ft
Equivalent damping coefficient of front suspension
12000
N·s/m
c
rt
Equivalent damping coefficient of rear suspension
15000
N·s/m
l
1
Distance between the center and seat of the cab
0.2
m
l
2
Distance between front suspended and the center of the cab
0.76
m
l
3
Distance between rear suspended and the center of the cab
0.87
m
l
4
Distance between the center of the cab and the center of target commercial vehicle
2.47
m
l
5
Distance between the front suspension and the center of target commercial vehicle
2.1
m
l
6
Distance between the rear suspension and the center of target commercial vehicle
3
mFigure 2
Damping curves. (a) Damping ratio. (b) Damping coefficient. (c) Damping force.
(a)
(b)
(c)


## 2.2. Road Excitation Model

As the most important excitation source when driving, the acquisition of road excitation is highly significant for ride comfort analysis. The road level and vehicle speed are primary factors affecting vehicle vibration. Power spectral density (PSD) is usually used for statistical characteristics of road excitation [13, 14]. The PSD of road excitation satisfies(2)Gqn=Gqn0nn0−w,where n is the spatial frequency; n0 is the reference spatial frequency, n0 = 0.1 m−1; Gqn0 is the road surface PSD corresponding to the reference spatial frequency n0, or the road excitation coefficient; and w is the frequency index that determines the frequency structure of the road surface PSD.The time-domain model of road excitation based on filtered white noise is as follows:(3)q˙t=−2πunqqt+2πn0Gqn0uωt,where u is the vehicle speed (m/s); nq is the cutoff space frequency, nq = 0.011 m−1; qt is the displacement excitation of the road surface; n0 represents the spatial reference frequency and is equal to 0.1 m−1; Gqn0 is the road excitation coefficient, Gqn0 = 256 × 10−6 m3; and ωt is the Gauss white noise with zero mean.Figure3 shows a road profile when the vehicle travels on the C-level road at 60 km/h. The profile PSD and standard road PSD according to ISO 8608 are presented in Figure 4 and show that a simulation model of road excitation is feasible.Figure 3
Road excitation time-domain profile.Figure 4
Road excitation frequency-domain profile.


## 2.3. Nonlinear Damping Model

To capture the damping nonlinear characteristic, a set of damping ratio curves are provided by a commercial vehicle manufacturer and shown in Figure2(a). The damping force is calculated by the following equation:(4)F=cv=2ξkmv,where F and c are the damping force and damping coefficient, respectively; ξ represents the damping ratio; k represents the stiffness; m refers to the sprung mass; and v is the compression or extension velocity.The damping force lines presented in Figure2(c) can be expressed as follows:(5)F=c1+η⋅sgnvvn,where η is the asymmetry coefficient and indicates the difference between the compression and extension coefficient and n is the damping characteristic index.


## 2.4. Equivalent Damping Model

Usually, the equivalent damping coefficient is calculated through the energy dissipated during one vibration period, which is equal to the energy dissipated by the equivalent damping. When a system is forced to vibrate, the viscous damping force can be expressed as the following equation:(6)x=Xsinωt−φ,x˙=Xωcosωt−φ,Fe=cx˙=cXωcosωt−φ.Energy consumption by equivalent damping in a vibration cycle is(7)We=∫0TFex˙dt=∫02π/ωceX2ω2cos2ωt−φ=πceωX2.When the damping of the system is nonlinear, it can be replaced by equivalent damping:(8)πceωX2=Wd,ce=WdωX2π,where Wd is the energy consumed by the equivalent damping during one vibration period and ce is the equivalent damping coefficient.


## 2.5. Verification of Ride Comfort Model

Equations (1), (3), and (5) have been solved using Matlab/Simulink and validated through a multibody dynamic model in Adams/Car. The multibody dynamic model is shown in Figure 5. The response of the driver seat, obtained from the Simulink models based on nonlinear damping, equivalent damping, and multibody dynamics model, is compared through time-domain curves and frequency-weighted RMS to verify the accuracy of the nonlinear damping ride model.Figure 5
Multibody dynamic model.The time-domain curves are shown in Figure6. Overall, the Matlab-based model has the same level of acceleration amplitude as the Adams-based model. The acceleration amplitude of the model based on nonlinear damping is greater than the one based on equivalent damping, leading to a higher frequency-weighted RMS, and the situation is consistent with the RMS in Figure 7. The difference in the curves is typically caused by the different methods of road modeling in Adams and Matlab and the bushings in the Adams model.Figure 6
Time-domain responses of the driver seat. (a) 40 km/h. (b) 80 km/h.
(a)
(b)Figure 7
Frequency-weighted RMS.As shown in Figure7, there is a similar tendency that ride comfort increases with speed among the three models. The results show that simulation results of the model based on the nonlinear damping are more consistent with the Adams model than the model based on equivalent damping. Therefore, the ride comfort model based on nonlinear damping has improved accuracy compared to equivalent damping models.


## 3. Multiobjective Optimization

In this section, the optimization of suspension systems based on nonlinear damping is presented. The efficiency and robustness of the optimization method affect the results. PSO, CS, DIRECT, and GA are chosen as the optimization methods to optimize performance with respect to comfort, dynamic tyre load, and suspension working space.


### 3.1. Intelligent Algorithms

#### 3.1.1. Particle Swarm Optimization

The PSO is an iterative algorithm formed by a family of particles in which each particle keeps track of its coordinates and shares them with the other particles. The particles fly in theN-dimensional problem space by learning from the historical information of all the particles [15, 16]. The particle i at iteration t + 1 should be updated according to the following equation:(9)Vik+1=Vik+c1×rand1×Pi,bestk−Xik+c2×rand2×Gi,bestk−Xik,Xik+1=Xik+Vik+1,where V and X are the velocity and position of the particles, respectively; k is the iteration; c1 and c2 are the constants called cognitive factor and social factor; Pi,best is the best position in the neighborhood of particle i; Gi,best is the best position found in the whole population of particles; and rand1 and rand2 represent random value between 0 and 1.


#### 3.1.2. Cuckoo Search

CS is a nature-inspired algorithm developed based on the evolution of cuckoo birds [17, 18]. Cuckoos typically lay their fertilized eggs in other cuckoos’ nests, hoping that their offspring will be raised by other cuckoos. There are times when the cuckoos discover that the eggs in their nests do not belong to them, and in these cases, the foreign eggs are either thrown out of the nest or the whole nests are abandoned. The three main idealized rules on which the algorithm is based are as follows:(1)
Each cuckoo has only one egg and selects the nest position randomly(2)
The best nests with high-quality eggs will be carried to the next generation(3)
The host nest number is constant during iteration, and the probability of the host finding eggs ispa ∈ [0, 1]The particlei update in the CS is performed according to following equationfd10:(10)Xig+1=Xig+α0Xbestg−Xig⊕Levyλ,where Xi represents the nest i position; α0 is a constant; Xbest is the current best nest position; ⊕ indicates entrywise multiplications; and Levy (λ) is a Levy distribution function describing randomly walked steps.A simple version of Levy distribution can be mathematically defined by(11)Levyλ∼ϕuv1/λ,where u and v are both standard normal distributions; λ = 1.5; and ϕ is calculated by the following equation:(12)ϕ=Γ1+λsinπλ/2Γ0.5λ1+λ2λ−1/21/λ,where Γ is known as a standard Gamma function.


#### 3.1.3. Dividing Rectangles

The DIRECT algorithm is composed of potentially optimal hyper-rectangles and the dividing strategy for hypercubes [19].We deﬁne the potentially optimal hyper-rectangle as follows: letε > 0 be a positive constant and let fmin be the current best function value. A hyper-rectangle j is then said to be potentially optimal if there exists some K > 0 such that(13)fmj−Kdj≤fmi−Kdi,∀i,fmj−Kdj≤fmin−fmin,where mj is the center of hyper-rectangle j and dj is a measure for this hyper-rectangle.Furthermore,m is the center point of a hypercube with the smallest function value si. One dimension is divided into thirds so that m − δei, m, and m + δei are the new centers for hyper-rectangles. The function is then evaluated at the points m ± δei, where δ is 1/3 of the side length of the hypercube and ei is the i-th Euclidean base vector. si is defined as(14)si=minfm+δei,fm−δei.Then,si is the new center of a hyper-rectangle with the smallest function value. This pattern is repeated for all dimensions on the center of the hyper-rectangle, choosing the next dimension by determining the next smallest si [20].


#### 3.1.4. Genetic Algorithm

The GA is an optimization method based on the principles of natural genetics and natural selection. The basic elements of GA are performed by three operations as reproduction, crossover, and mutation. Solutions are evaluated with respect to their fitness value that indicates how well the individual will solve the problem [21].


### 3.2. Objectives

#### 3.2.1. Ride Comfort

The vibration perception of the human body is related to the frequency, response of the driver seat, and vibration acceleration and has been weighted according to the ISO 2631:1997 standard. This standard evaluates human exposure to whole-body vibration and is usually based on one objective function. The frequency-weighted RMS valueaw is then defined as follows [22]:(15)aw=∫0.580W2fGafdf1/2,where Wf and Gaf can be expressed by the following equations, respectively:(16)Wf=0.5,0.5Hz<f<2Hz,f4,2Hz<f<4Hz,1,4Hz<f<12.5Hz,12.5f,12.5Hz<f<80Hz,(17)Gaf=∫−∞+∞atexp−i2πωtdt,where at is the time-domain signal of the driver seat.


#### 3.2.2. Dynamic Tyre Load

The dynamic tyre load is defined as the RMS value of the change in tyre load relative to the static equilibrium position. The mechanism of tyre adhesion loss caused by the tyre load is described as the following: when a tyre is required to generate lateral or longitudinal force, the contact portion between the tyre and the ground must be deformed before the force is sufficiently generated. However, this deformation requires a roll distance, so there is a time delay before the tyre force is fully obtained. When the tyre load fluctuates with the suspension motion, the effective lateral or longitudinal force available is reduced because of the tyre dynamic mechanism. Therefore, if a stable normal load of the tyre can be maintained, a large tyre force can be obtained; if the dynamic load fluctuation of the tyre increases, the tyre grip ability weakens as the tyre jumps. The dynamic tyre load can be expressed as(18)DTLf=kftxfu−xfmfu+m1g,DTLr=krtxru−xrmru+m2g,(19)RMSDTLf=1T∫0TDTLf2tdt,RMSDTLr=1T∫0TDTLr2tdt.


#### 3.2.3. Suspension Working Space

The suspension working space is defined as the difference between the wheel and the body displacement. According to the Gaussian distribution of random roads, for linear systems, the response should have a Gaussian property and can be described by a normal distribution. Therefore, the suspension working space can be considered the probability that the relative displacement of the wheel and the body remains within the SWS, 2SWS, and 3SWS under static equilibrium position conditions 68.3%, 95.4%, and 99.7% of the time, respectively. According to the RMS value of the suspension working space, the suspension working space required by the vehicle under certain road input conditions can be determined:(20)SWSf=xb−xfu−l5θb,SWSr=xb−xru−l6θb,(21)RMSSWSf=1T∫0TSWSf2tdt,RMSSWSr=1T∫0TSWSr2tdt.In summary, there are the following objective functions:(22)fX=minf1X,f2X,f3X,f4X,f5X=aw+RMSSWSf+RMSSWSr+RMSDTLf+RMSDTLr.


### 3.3. Design Variables

Figure1 and equation (1) show that there is a complex nonlinear relationship between the mechanical parameters, quality parameters, and geometric parameters of commercial vehicles. By considering all factors, it is not conducive to optimize ride comfort. In addition, the suspension system is the key mechanical structure affecting the ride performance of vehicles. The parameter analysis and optimization of suspension have long been an important research topic [23, 24]. The stiffness and damping of suspension have the greatest influence and are coupled with each other to affect ride comfort. Therefore, the stiffness and damping coefficients of suspension of each force section are considered design variables:(23)X=x1,x2,x3,x4,x5,x6,x7,x8T=kfs,krs,cf1,cf2,cf3,cr1,cr2,cr3T.The range of variables should not be too large or too small. When the variation ranges are too large, other performances of the vehicle may be drastically reduced after optimization, and the matching and installation of the entire vehicle will be affected. In contrast, the effect of optimization is not obvious. Therefore, variables are altered up and down by 10% as the optimization interval, and the bounds are defined in Table2.Table 2
Bounds of design variables.
Variables
Original
Lower
Upper
kfs (N/m)
410000
369000
451000
krs (N/m)
1476000
1328400
1623600
c
f1
36.54
32.886
40.194
c
f2
1.487
1.338
1.636
c
f3
0.5996
0.540
0.660
c
r1
49.66
44.694
54.626
c
r2
1.487
1.338
1.636
c
r3
0.5996
0.540
0.660


### 3.4. Constraints

The static deflection of the suspension is the ratio of the sprung mass to the stiffness under full loads. Typically, the static deflection of the suspension is 50–110 mm. At the same time, having a static deflection of the rear suspensionSr smaller than the static deflection of the front suspension Sf helps prevent large longitudinal angular vibration of the vehicle body (the longitudinal angular vibration is smaller):(24)0.05≤Sf=mfg2Kf≤0.11,0.05≤Sr=mrg2Kr≤0.11,0.6mfmr≤kfskrs≤0.8mfmr,where mf and mr are the sprung mass of the front and rear axles, respectively, and Kf and Kr represent the one-side stiffness of the front and rear suspensions, respectively.The natural frequency of suspension usually ranges from 1.5 to 2.2. Thus,(25)1.5≤12πkfsmfs≤2.0,1.7≤12πkrsmrs≤2.2.The probability that the wheel jumps off the ground is 0.15%, and the RMS value of the dynamic tyre load must satisfy the following equation:(26)DTLf≤13,DTLr≤13.


### 3.5. Nonlinear Ride Optimization Method

A ride potimization method based on nonlinear damping and intelligent lagorithms is presented, and the procedure is shown in Figure8.Figure 8
Ride comfort optimization procedure based on nonlinear damping.


## 3.1. Intelligent Algorithms

### 3.1.1. Particle Swarm Optimization

The PSO is an iterative algorithm formed by a family of particles in which each particle keeps track of its coordinates and shares them with the other particles. The particles fly in theN-dimensional problem space by learning from the historical information of all the particles [15, 16]. The particle i at iteration t + 1 should be updated according to the following equation:(9)Vik+1=Vik+c1×rand1×Pi,bestk−Xik+c2×rand2×Gi,bestk−Xik,Xik+1=Xik+Vik+1,where V and X are the velocity and position of the particles, respectively; k is the iteration; c1 and c2 are the constants called cognitive factor and social factor; Pi,best is the best position in the neighborhood of particle i; Gi,best is the best position found in the whole population of particles; and rand1 and rand2 represent random value between 0 and 1.


### 3.1.2. Cuckoo Search

CS is a nature-inspired algorithm developed based on the evolution of cuckoo birds [17, 18]. Cuckoos typically lay their fertilized eggs in other cuckoos’ nests, hoping that their offspring will be raised by other cuckoos. There are times when the cuckoos discover that the eggs in their nests do not belong to them, and in these cases, the foreign eggs are either thrown out of the nest or the whole nests are abandoned. The three main idealized rules on which the algorithm is based are as follows:(1)
Each cuckoo has only one egg and selects the nest position randomly(2)
The best nests with high-quality eggs will be carried to the next generation(3)
The host nest number is constant during iteration, and the probability of the host finding eggs ispa ∈ [0, 1]The particlei update in the CS is performed according to following equationfd10:(10)Xig+1=Xig+α0Xbestg−Xig⊕Levyλ,where Xi represents the nest i position; α0 is a constant; Xbest is the current best nest position; ⊕ indicates entrywise multiplications; and Levy (λ) is a Levy distribution function describing randomly walked steps.A simple version of Levy distribution can be mathematically defined by(11)Levyλ∼ϕuv1/λ,where u and v are both standard normal distributions; λ = 1.5; and ϕ is calculated by the following equation:(12)ϕ=Γ1+λsinπλ/2Γ0.5λ1+λ2λ−1/21/λ,where Γ is known as a standard Gamma function.


### 3.1.3. Dividing Rectangles

The DIRECT algorithm is composed of potentially optimal hyper-rectangles and the dividing strategy for hypercubes [19].We deﬁne the potentially optimal hyper-rectangle as follows: letε > 0 be a positive constant and let fmin be the current best function value. A hyper-rectangle j is then said to be potentially optimal if there exists some K > 0 such that(13)fmj−Kdj≤fmi−Kdi,∀i,fmj−Kdj≤fmin−fmin,where mj is the center of hyper-rectangle j and dj is a measure for this hyper-rectangle.Furthermore,m is the center point of a hypercube with the smallest function value si. One dimension is divided into thirds so that m − δei, m, and m + δei are the new centers for hyper-rectangles. The function is then evaluated at the points m ± δei, where δ is 1/3 of the side length of the hypercube and ei is the i-th Euclidean base vector. si is defined as(14)si=minfm+δei,fm−δei.Then,si is the new center of a hyper-rectangle with the smallest function value. This pattern is repeated for all dimensions on the center of the hyper-rectangle, choosing the next dimension by determining the next smallest si [20].


### 3.1.4. Genetic Algorithm

The GA is an optimization method based on the principles of natural genetics and natural selection. The basic elements of GA are performed by three operations as reproduction, crossover, and mutation. Solutions are evaluated with respect to their fitness value that indicates how well the individual will solve the problem [21].


## 3.1.1. Particle Swarm Optimization

The PSO is an iterative algorithm formed by a family of particles in which each particle keeps track of its coordinates and shares them with the other particles. The particles fly in theN-dimensional problem space by learning from the historical information of all the particles [15, 16]. The particle i at iteration t + 1 should be updated according to the following equation:(9)Vik+1=Vik+c1×rand1×Pi,bestk−Xik+c2×rand2×Gi,bestk−Xik,Xik+1=Xik+Vik+1,where V and X are the velocity and position of the particles, respectively; k is the iteration; c1 and c2 are the constants called cognitive factor and social factor; Pi,best is the best position in the neighborhood of particle i; Gi,best is the best position found in the whole population of particles; and rand1 and rand2 represent random value between 0 and 1.


## 3.1.2. Cuckoo Search

CS is a nature-inspired algorithm developed based on the evolution of cuckoo birds [17, 18]. Cuckoos typically lay their fertilized eggs in other cuckoos’ nests, hoping that their offspring will be raised by other cuckoos. There are times when the cuckoos discover that the eggs in their nests do not belong to them, and in these cases, the foreign eggs are either thrown out of the nest or the whole nests are abandoned. The three main idealized rules on which the algorithm is based are as follows:(1)
Each cuckoo has only one egg and selects the nest position randomly(2)
The best nests with high-quality eggs will be carried to the next generation(3)
The host nest number is constant during iteration, and the probability of the host finding eggs ispa ∈ [0, 1]The particlei update in the CS is performed according to following equationfd10:(10)Xig+1=Xig+α0Xbestg−Xig⊕Levyλ,where Xi represents the nest i position; α0 is a constant; Xbest is the current best nest position; ⊕ indicates entrywise multiplications; and Levy (λ) is a Levy distribution function describing randomly walked steps.A simple version of Levy distribution can be mathematically defined by(11)Levyλ∼ϕuv1/λ,where u and v are both standard normal distributions; λ = 1.5; and ϕ is calculated by the following equation:(12)ϕ=Γ1+λsinπλ/2Γ0.5λ1+λ2λ−1/21/λ,where Γ is known as a standard Gamma function.


## 3.1.3. Dividing Rectangles

The DIRECT algorithm is composed of potentially optimal hyper-rectangles and the dividing strategy for hypercubes [19].We deﬁne the potentially optimal hyper-rectangle as follows: letε > 0 be a positive constant and let fmin be the current best function value. A hyper-rectangle j is then said to be potentially optimal if there exists some K > 0 such that(13)fmj−Kdj≤fmi−Kdi,∀i,fmj−Kdj≤fmin−fmin,where mj is the center of hyper-rectangle j and dj is a measure for this hyper-rectangle.Furthermore,m is the center point of a hypercube with the smallest function value si. One dimension is divided into thirds so that m − δei, m, and m + δei are the new centers for hyper-rectangles. The function is then evaluated at the points m ± δei, where δ is 1/3 of the side length of the hypercube and ei is the i-th Euclidean base vector. si is defined as(14)si=minfm+δei,fm−δei.Then,si is the new center of a hyper-rectangle with the smallest function value. This pattern is repeated for all dimensions on the center of the hyper-rectangle, choosing the next dimension by determining the next smallest si [20].


## 3.1.4. Genetic Algorithm

The GA is an optimization method based on the principles of natural genetics and natural selection. The basic elements of GA are performed by three operations as reproduction, crossover, and mutation. Solutions are evaluated with respect to their fitness value that indicates how well the individual will solve the problem [21].


## 3.2. Objectives

### 3.2.1. Ride Comfort

The vibration perception of the human body is related to the frequency, response of the driver seat, and vibration acceleration and has been weighted according to the ISO 2631:1997 standard. This standard evaluates human exposure to whole-body vibration and is usually based on one objective function. The frequency-weighted RMS valueaw is then defined as follows [22]:(15)aw=∫0.580W2fGafdf1/2,where Wf and Gaf can be expressed by the following equations, respectively:(16)Wf=0.5,0.5Hz<f<2Hz,f4,2Hz<f<4Hz,1,4Hz<f<12.5Hz,12.5f,12.5Hz<f<80Hz,(17)Gaf=∫−∞+∞atexp−i2πωtdt,where at is the time-domain signal of the driver seat.


### 3.2.2. Dynamic Tyre Load

The dynamic tyre load is defined as the RMS value of the change in tyre load relative to the static equilibrium position. The mechanism of tyre adhesion loss caused by the tyre load is described as the following: when a tyre is required to generate lateral or longitudinal force, the contact portion between the tyre and the ground must be deformed before the force is sufficiently generated. However, this deformation requires a roll distance, so there is a time delay before the tyre force is fully obtained. When the tyre load fluctuates with the suspension motion, the effective lateral or longitudinal force available is reduced because of the tyre dynamic mechanism. Therefore, if a stable normal load of the tyre can be maintained, a large tyre force can be obtained; if the dynamic load fluctuation of the tyre increases, the tyre grip ability weakens as the tyre jumps. The dynamic tyre load can be expressed as(18)DTLf=kftxfu−xfmfu+m1g,DTLr=krtxru−xrmru+m2g,(19)RMSDTLf=1T∫0TDTLf2tdt,RMSDTLr=1T∫0TDTLr2tdt.


### 3.2.3. Suspension Working Space

The suspension working space is defined as the difference between the wheel and the body displacement. According to the Gaussian distribution of random roads, for linear systems, the response should have a Gaussian property and can be described by a normal distribution. Therefore, the suspension working space can be considered the probability that the relative displacement of the wheel and the body remains within the SWS, 2SWS, and 3SWS under static equilibrium position conditions 68.3%, 95.4%, and 99.7% of the time, respectively. According to the RMS value of the suspension working space, the suspension working space required by the vehicle under certain road input conditions can be determined:(20)SWSf=xb−xfu−l5θb,SWSr=xb−xru−l6θb,(21)RMSSWSf=1T∫0TSWSf2tdt,RMSSWSr=1T∫0TSWSr2tdt.In summary, there are the following objective functions:(22)fX=minf1X,f2X,f3X,f4X,f5X=aw+RMSSWSf+RMSSWSr+RMSDTLf+RMSDTLr.


## 3.2.1. Ride Comfort

The vibration perception of the human body is related to the frequency, response of the driver seat, and vibration acceleration and has been weighted according to the ISO 2631:1997 standard. This standard evaluates human exposure to whole-body vibration and is usually based on one objective function. The frequency-weighted RMS valueaw is then defined as follows [22]:(15)aw=∫0.580W2fGafdf1/2,where Wf and Gaf can be expressed by the following equations, respectively:(16)Wf=0.5,0.5Hz<f<2Hz,f4,2Hz<f<4Hz,1,4Hz<f<12.5Hz,12.5f,12.5Hz<f<80Hz,(17)Gaf=∫−∞+∞atexp−i2πωtdt,where at is the time-domain signal of the driver seat.


## 3.2.2. Dynamic Tyre Load

The dynamic tyre load is defined as the RMS value of the change in tyre load relative to the static equilibrium position. The mechanism of tyre adhesion loss caused by the tyre load is described as the following: when a tyre is required to generate lateral or longitudinal force, the contact portion between the tyre and the ground must be deformed before the force is sufficiently generated. However, this deformation requires a roll distance, so there is a time delay before the tyre force is fully obtained. When the tyre load fluctuates with the suspension motion, the effective lateral or longitudinal force available is reduced because of the tyre dynamic mechanism. Therefore, if a stable normal load of the tyre can be maintained, a large tyre force can be obtained; if the dynamic load fluctuation of the tyre increases, the tyre grip ability weakens as the tyre jumps. The dynamic tyre load can be expressed as(18)DTLf=kftxfu−xfmfu+m1g,DTLr=krtxru−xrmru+m2g,(19)RMSDTLf=1T∫0TDTLf2tdt,RMSDTLr=1T∫0TDTLr2tdt.


## 3.2.3. Suspension Working Space

The suspension working space is defined as the difference between the wheel and the body displacement. According to the Gaussian distribution of random roads, for linear systems, the response should have a Gaussian property and can be described by a normal distribution. Therefore, the suspension working space can be considered the probability that the relative displacement of the wheel and the body remains within the SWS, 2SWS, and 3SWS under static equilibrium position conditions 68.3%, 95.4%, and 99.7% of the time, respectively. According to the RMS value of the suspension working space, the suspension working space required by the vehicle under certain road input conditions can be determined:(20)SWSf=xb−xfu−l5θb,SWSr=xb−xru−l6θb,(21)RMSSWSf=1T∫0TSWSf2tdt,RMSSWSr=1T∫0TSWSr2tdt.In summary, there are the following objective functions:(22)fX=minf1X,f2X,f3X,f4X,f5X=aw+RMSSWSf+RMSSWSr+RMSDTLf+RMSDTLr.


## 3.3. Design Variables

Figure1 and equation (1) show that there is a complex nonlinear relationship between the mechanical parameters, quality parameters, and geometric parameters of commercial vehicles. By considering all factors, it is not conducive to optimize ride comfort. In addition, the suspension system is the key mechanical structure affecting the ride performance of vehicles. The parameter analysis and optimization of suspension have long been an important research topic [23, 24]. The stiffness and damping of suspension have the greatest influence and are coupled with each other to affect ride comfort. Therefore, the stiffness and damping coefficients of suspension of each force section are considered design variables:(23)X=x1,x2,x3,x4,x5,x6,x7,x8T=kfs,krs,cf1,cf2,cf3,cr1,cr2,cr3T.The range of variables should not be too large or too small. When the variation ranges are too large, other performances of the vehicle may be drastically reduced after optimization, and the matching and installation of the entire vehicle will be affected. In contrast, the effect of optimization is not obvious. Therefore, variables are altered up and down by 10% as the optimization interval, and the bounds are defined in Table2.Table 2
Bounds of design variables.
Variables
Original
Lower
Upper
kfs (N/m)
410000
369000
451000
krs (N/m)
1476000
1328400
1623600
c
f1
36.54
32.886
40.194
c
f2
1.487
1.338
1.636
c
f3
0.5996
0.540
0.660
c
r1
49.66
44.694
54.626
c
r2
1.487
1.338
1.636
c
r3
0.5996
0.540
0.660


## 3.4. Constraints

The static deflection of the suspension is the ratio of the sprung mass to the stiffness under full loads. Typically, the static deflection of the suspension is 50–110 mm. At the same time, having a static deflection of the rear suspensionSr smaller than the static deflection of the front suspension Sf helps prevent large longitudinal angular vibration of the vehicle body (the longitudinal angular vibration is smaller):(24)0.05≤Sf=mfg2Kf≤0.11,0.05≤Sr=mrg2Kr≤0.11,0.6mfmr≤kfskrs≤0.8mfmr,where mf and mr are the sprung mass of the front and rear axles, respectively, and Kf and Kr represent the one-side stiffness of the front and rear suspensions, respectively.The natural frequency of suspension usually ranges from 1.5 to 2.2. Thus,(25)1.5≤12πkfsmfs≤2.0,1.7≤12πkrsmrs≤2.2.The probability that the wheel jumps off the ground is 0.15%, and the RMS value of the dynamic tyre load must satisfy the following equation:(26)DTLf≤13,DTLr≤13.


## 3.5. Nonlinear Ride Optimization Method

A ride potimization method based on nonlinear damping and intelligent lagorithms is presented, and the procedure is shown in Figure8.Figure 8
Ride comfort optimization procedure based on nonlinear damping.


## 4. Optimization Results

In this section, the results of ride optimization are discussed. First, the iterations of algorithms are shown in Figure9, and the optimal design variables are listed in Table 3. The responses, which represent the time-domain signals of objective functions based on original design variables and optimized ones, are shown in Figure 9. The indexes for objectives are calculated and presented in Table 4.Figure 9
Objective function iteration.Table 3
Optimization variables.
Variables
kfs (N/m)
krs (N/m)
c
f1
c
f2
c
f3
c
r1
c
r2
c
r3
Original
410000
1476000
36.54
1.487
0.5996
49.66
1.487
0.5996
CS
369000
1330000
40.18
1.635
0.6599
54.63
1.636
0.6600
DIRECT
371000
1345000
40.06
1.631
0.6578
52.97
1.608
0.6578
GA
381000
1340000
40.19
1.636
0.6600
54.62
1.636
0.6600
PSO
369000
1328000
40.19
1.636
0.6600
54.63
1.636
0.6600Table 4
Evaluation indexes of objective function before and after optimization.
Speed (km/h)
30
40
50
60
70
80
90
100
f
1
Original
0.4107
0.5190
0.4430
0.6077
0.6988
0.6914
0.8738
0.7832
Optimized
0.3067
0.3531
0.3527
0.4318
0.4742
0.5284
0.5786
0.6012
Fall ratio (%)
25.3
32.0
20.4
28.9
32.1
23.6
33.8
23.3
Average (%)
27.4
f
2
Original
0.2442
0.2903
0.3118
0.3740
0.4085
0.4042
0.4731
0.4755
Optimized
0.1922
0.2263
0.2557
0.2908
0.3122
0.3275
0.3584
0.3701
Fall ratio (%)
21.3
22.0
18.0
22.2
23.6
19.0
24.2
22.2
Average (%)
21.6
f
3
Original
0.1413
0.1672
0.2050
0.2405
0.2577
0.2674
0.3045
0.3021
Optimized
0.1075
0.1314
0.1567
0.1756
0.1875
0.2014
0.2196
0.2282
Fall ratio (%)
23.9
21.4
23.6
27.0
27.2
24.7
27.9
24.5
Average (%)
25.0
f
4
Original
0.0142
0.0166
0.0171
0.0215
0.0233
0.0218
0.0271
0.0261
Optimized
0.0113
0.0132
0.0148
0.0172
0.0183
0.0186
0.0207
0.0209
Fall ratio (%)
20.4
20.5
13.5
20.0
21.5
14.7
23.6
19.9
Average (%)
19.3
f
5
Original
0.0086
0.0101
0.0125
0.0150
0.0159
0.0161
0.0188
0.0183
Optimized
0.0068
0.0082
0.0099
0.0114
0.0119
0.0127
0.0139
0.0144
Fall ratio (%)
20.9
18.8
20.8
24.0
25.2
21.1
26.1
21.3
Average (%)
22.3The iterations of ride comfort optimization are shown in Figure9. Obviously, the CS, DIRECT, and PSO get the similar and ideal results besides GA, and the iteration of PSO is better than others, which is usually caused by the following reasons: (1) different algorithms have different characteristic for various problems and (2) the generation of individuals and population is random during the optimization processing. Consequently, it is necessary to compare the results among different algorithms that can ensure a relatively better result. Then, we listed the original variables and optimized ones in Table 3.To verify the feasibility and effectiveness of the ride comfort optimization method, nonlinear damping based on intelligent algorithms is proposed in this paper. A comparative analysis of ride comfort based on the original and optimized design variables is performed. Design variables are substituted into the ride comfort model based on nonlinear damping. Then, the optimized responses of the time history of the acceleration of driver seat, dynamic tyre load, and suspension working space at 40 km/h and 80 km/h are obtained. Data are presented in Figure10 in comparison with the original response. The evaluation indexes for objective functions are calculated from 30 to 100 km/h and shown in Table 4 according to equations (15), (19), and (21)Figure 10
Comparison between original and optimized objective functions. (a, b) Acceleration of the driver seat at 40 and 80 km/h, respectively. (c, d) Dynamic tyre load of the front wheel at 40 and 80 km/h, respectively. (e, f) Dynamic tyre load of the rear wheel at 40 and 80 km/h, respectively. (g, h) Suspension working space of the front axle at 40 and 80 km/h, respectively. (i, j) Suspension working space of the rear axle at 40 and 80 km/h, respectively.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
(j)As shown in Table4, the target commercial vehicle ride comfort is significantly improved after optimization. The proposed optimizer with intelligent algorithms and nonlinear damping successfully optimized the design variables with an average fall ratio of 27.4%, 21.6%, 25.0%, 19.3%, and 22.3% for f1, f2, f3, f4, and f5. Additionally, the evaluation indexes for objective increase with increased vehicle speed; that is, ride comfort deteriorates totally as speed increases.


## 5. Experiment

To further verify the feasibility and efficiency of the ride optimization method in this paper, the ride comfort experiment of target commercial vehicle traveling on theC-level road with 30–100 km/h was conducted in accordance with GB/T 4970-2009 “Method of running test—Automotive ride comfort.” The experimental vehicle is shown in Figure 11(a). A section of the cement road is selected as the ride comfort experiment road, and its surface is shown in Figure 11(b). The sensor is used to measure the vertical vibration acceleration of the seat. The installation position is shown in Figure 11(c). The CPCI signal acquisition instrument is used for data acquisition, as shown in Figure 11(d).Figure 11
Target commercial vehicle ride comfort experiment.The ride comfort experiment test requirements are as follows: (1) The road surface is even and straight, and there is an acceleration section initially. (2) The target commercial vehicle is fully loaded with standard driving equipment. The components of the vehicle are in good condition, and the tyre pressure complies with regulations. (3) When the test truck speed reaches the speed to be analyzed, the tester starts timing for 60 s. (4) The speed of vehicle ride comfort experiment is 30–100 km/h, and the speed is increased by 5 km/h each test. The target commercial vehicle drives evenly at the specified speed during the experiment.Figure12 shows that the frequency-weighted RMS of optimized seat vibration acceleration is significantly lower than the original situation. The results in Figure 12 are similar to those in Figure 7, demonstrating that the ride comfort model based on nonlinear damping constructed in this paper is feasible and the optimization method is effective. In addition, the fluctuations of the experimental result after 70 km/h may be affected by engine vibration because of increased rotating speed.Figure 12
Frequency-weighted RMS of experiment.


## 6. Conclusions

In this paper, we summarized the situation and shortcoming for the ride comfort modeling method, currently. Furthermore, the vibration model of commercial vehicle, road excitation model, nonlinear damping model, and equivalent damping model are successively developed. The accuracy of a ride comfort model using nonlinear damping is demonstrated and compared with other models. Then, the ride comfort model based on the nonlinear damping is optimized using four intelligent algorithms. Lastly, the effectiveness of the optimization method is tested according to the performances before and after optimization. The main conclusions are as follows:(1)
The frequency-weighted RMS values of seat vibration acceleration based on models with nonlinear damping and equivalent damping are obtained and compared with the response of a multibody dynamic model. The result indicates that the ride comfort model based on nonlinear damping has higher accuracy than a model based on equivalent damping.(2)
The comparison of driver seat vibration acceleration and relevant responses before and after optimization is performed, and the evaluation indexes at different speeds are obtained. The frequency-weighted RMS of driver seat acceleration, RMS of suspension working space of front and rear axle, and RMS of dynamic tyre load of front and rear wheel are decreased by an average of 27.4%, 21.6%, 25.0%, 19.3%, and 22.3%, respectively, through the optimization. Thus, the ride comfort of the commercial vehicle is improved. In addition, the optimization method developed in this paper has theoretical and practical significance for similar problems.



---

*Source: 2973190-2019-08-19.xml*