Ride comfort investigation of semi-active seat suspension integrated with quarter car model

. A method for parameter identi ﬁ cation of the magnetorheological damper (MRD) model with an improved ﬁ re ﬂ y algorithm (IFA) is proposed, and a semi-active seat control system with three-degree-of-freedom (3-DOF) is established by combining with a quarter car model to investigate the ride comfort. The dynamic characteristics of the MRD were analyzed by experimental method. Combined with the IFA, the parameters of the MRD phenomenon model were identi ﬁ ed, and the forward model of the MR damper was constructed. The semi-active control model of a 3-DOF seat suspension was established. The MRD controller and suspension system controller were designed. The passive control, PID control, and Fuzzy-PID control on the vibration reduction of the semi-active seat suspension were compared and analyzed, under different road excitation. The simulation results show that the semi-active seat suspension controlled by the PID and Fuzzy-PID can effectively reduce the seat acceleration and dynamic stroke, which signi ﬁ cantly improve the ride comfort and operation safety compared to the passive seat suspension.


Introduction
In the process of vehicle driving, the vibration caused by uneven road surfaces has a significant impact on the physical and mental health of drivers.It is easy to cause physical and psychological discomfort, resulting in waist and shoulder pain, muscle pain, fatigue, digestive problems, visual impairment, and other diseases [1][2][3].The seat is in direct contact with the driver, transmitting low-frequency vibrations to the human body filtered by the tire and car body vibration isolation structure.The seat suspension system with a good damping effect can improve the driver's comfort and operation efficiency and reduce the damage of vibration.Seat suspension systems are mainly divided into passive suspension system, semiactive suspension system, and active suspension system.Passive seat suspension is generally composed of spring and hydraulic damper.It has a simple structure and low price, which is most widely used.However, its suspension parameters cannot be adjusted, and it is only suitable for a specific road excitation [4].Active seat suspension is generally composed of springs, passive dampers, and actuators in parallel.It can be controlled according to different road excitation to obtain outstanding performances [5,6].However, the active seat suspension has the disadvantages of complexity, high cost, and highenergy consumption.Compared with active and passive suspension systems, the semi-active suspension system with adjustable damping has the characteristics of quick action response, wide dynamic range, good reversible performance, small volume, low energy consumption, continuous damping with adjustable range, and highfrequency response.It has attracted extensive attentions in vibration reduction and isolation [7,8].At present, the most studied semi-active suspension is composed of a magnetorheological damper (MRD) and a spring in parallel, which achieves the best trade-off between cost and performance.The MRD uses a smart material, magnetorheological fluids, which is a mixture of micron ferromagnetic particles and non-magnetic fluid carriers that can change according to the magnetic field to adjust the damping force of the MRD.When the MRD is not energized, no magnetic field is generated, and the MRD is the same as the ordinary hydraulic damper.When the MRD is energized, the magnetic field makes ferromagnetic particles reassemble into chain structure along the magnetic field direction, forming Bingham fluid effect [9].So far, many scholars have studied the control of semiactive seat suspension and put forward different control methods.As shown in Figure 1, two controllers need to be designed for semi-active seat suspension control: system controller and damper controller.According to the dynamic stroke of the seat suspension and the parameters of the seat body mass acceleration, the system controller controls the damping force required by the MRD output.The MRD controller converts the damping force into the voltage value required to drive the MRD, thus realizing the control of the seat suspension.The effectiveness of the MRD controller mainly depends on its ability to deal with nonlinearity and hysteresis.The existing MRD controllers include Heaviside function method [10,11], continuous state controller [12], sigmoid function method [13], symbolic function controller [14], neural network and other intelligent algorithm controllers [15,16].The MRD controller based on the Heaviside function method is a simple "on-off" controller.The voltage value is switched between the minimum value and the maximum value.It has the advantages of simple principles and fast response.The signum function control algorithm is an improvement of the Heaviside function method.It takes samples near a certain point and averages all the values to get the control voltage of the corresponding point.In this way, the control voltage can be discretized under the level below the maximum value, so that it is constantly close to the control voltage value that needs to be generated.The continuous state controller allows continuous voltage intervals between discontinuous pulses to be evaluated.The sigmoid function method is also one of the continuous control algorithms.Various intelligent algorithm controllers such as neural networks are non-parametric modelling and control.It has high requirements on hardware that the computer memory and CPU must meet high requirements for computing.Heaviside function algorithm is used to control the MRD in this study.
Theoretical and experimental results show that different control algorithms have different effects on semi-active suspension with the MRD [17].There are many control algorithms in the control of semi-active seat suspension systems.Each control algorithm has its advantages and disadvantages, among which PID control and fuzzy control are widely used [18].PID control is a classical control with mature control theory, simple algorithm, easy implementation, and good robustness [19].However, it is challenging to set PID control parameters and has poor adaptability to nonlinear systems.Fuzzy control is mainly used to solve the control problems of non-linear systems, difficult to establish mathematical models, and uncertain internal disturbances of multiple-input multiple-output systems (MIMO) [20].Singh and Aggarwal [21] compared and analyzed PID control, fuzzy control, and fuzzy-PID control strategies through the semi-active suspension system of the body.The simulation results show that the fuzzy-PID control is better than the other two control algorithms.Ning et al. [22] used the takagi-Sugeno fuzzy controller to study the control of an active seat suspension.They verified its control effect and stability through simulation and test.Chao et al. [23] introduced the fuzzy adaptive cuckoo search algorithm into the Fuzzy-PID control.They verified that the proposed algorithm could improve the ride comfort and adhesion ability through the simulation analysis of the suspension system.Alfadhli et al. [24] used a PSO optimized fuzzy control algorithm to analyze the active seat suspension system, and the simulation results showed the effectiveness of the control method on ride comfort.Phu et al. [25] proposed a new adaptive hybrid controller based on H ∞ , sliding mode control, and PID control to control the magnetorheological seat suspension system.They verified the effectiveness of the controller through experiments.The control algorithm is developing towards intelligent control and hybrid control.To explore the response characteristics of the semiactive seat suspension system integrated with a quarter car model under different road profiles, a three-degree-offreedom (3-DOF) seat suspension model was established.Dynamic stroke and acceleration of the seat suspension system are used as indicators to analyze and compare the influence of the PID controller and Fuzzy-PID controller on ride comfort.
2 System model and parameter identification 2.1 Quarter car model with MRD seat system Figure 2 shows the semi-active seat suspension system integrated with a quarter car model.The system is a 3-DOF model composed of a 1-DOF semi-active seat suspension and a 2-DOF quarter car model.m u is the unsprung mass; m s is the sprung mass; m seat is the mass of the seat and the driver; k t , k s , and k seat are the tires, quarter car suspension, and seat suspension stiffness, respectively; c s is the damping coefficient of the quarter car suspension; the tire damping is ignored; and seat suspension damper is an adjustable damper composed of the MRD.Z r is the displacement function of road irregularity, Z u , Z s , and Z seat are the vertical displacement coordinates of the wheel axle, body, and seat cushion upper surface, respectively, and the coordinate origin is at their respective equilibrium positions.F MRD is the damping force generated by the MRD.Parameters of the 3-DOF suspension system are obtained from the literatures [26,27], and the specific parameters are listed in Table 1.
According to Newton's second law, the dynamic differential equation of the 3-DOF semi-active seat suspension system is built as follows: Take the state vector as , and convert the motion differential equation of the semiactive suspension system into the system state equation and output equation: where, The parameters of the 3-DOF seat suspension system are described by the phenomenon model [28], and the F MRD equations of the seat adjustable damping force are as follows: where, v is the voltage required to control the driving current of the MRD; u is the first-order filter input control voltage; Z STD represents the seat travel deformation; Z 0 represents the initial displacement of the spring k 1 , which is  related to the standard damping force generated by the compensator; Z ID is the internal displacement of the MRD, which is a fictitious variable; z is the evolution variable of the hysteresis effect produced by the damper; k 1 represents the stiffness coefficient of the compensator; c 1 is used to generate the damping coefficient at low speeds in the forcevelocity relationship; a, b, g, and A are used to control the size and shape of the nonlinear hysteresis cycles; c 0 and k 0 are the viscous damping coefficient and stiffness coefficient at high speed, respectively; h is the filter gain, which determines the change rate of the magnetic field to the equilibrium state of the magnetorheological fluid.

Parameter identification of the phenomenon model
According to equations ( 6)-( 13), 14 parameters need to be determined for the phenomenon model of the MRD.In this study, the experimental data were obtained by the experimental method, and the dynamic model of the MRD was established by the improved firefly algorithm (IFA).Wu et al. [29] proposed an improved firefly algorithm in 2020 based on Xin-She Yang's classic firefly algorithm [30] and the Lévy-flight firefly algorithm [31].This The intensity of light emitted by fireflies is inversely proportional to the square of the distance and decreases with increasing distance.The mathematical expression is as follows:  where I i is the absolute brightness of the firefly i, g is the light absorption coefficient, I ij is the relative brightness of the firefly i to the firefly j, r ij is the distance between firefly i and firefly j.
Fireflies with high absolute brightness attract fireflies with low absolute brightness nearby and constantly gather around them.The mathematical expressions of the firefly attractiveness b, distance r ij , and updated position x i are established as shown in formulas ( 15)- (17).
where, b 0 is the maximum attractiveness (r = 0); d is dimension; d is a random number between [0,1]; ⊗is Hadamard product; h is the uniform random vector in [-1,1]; g is the switch parameter of global search and local search, which is a random number in [0,1]; sign ( • ) is a symbolic function; b is a constant that defines the shape of the logarithmic spiral; x i,t and x j,t are the positions of the firefly i and the firefly j at the time t, respectively; Levy is the Lévy-flight random number, which is calculated by equations ( 18) and (19).
where m 1 and m 2 conform the standard normal distributions, b = 1.5, and ' is calculated as follows, whereG (<>)is the standard Gamma function.The 14 parameters to be identified of the phenomenon model are expressed in vector form as follows: Phenomenal model parameter identification is a constrained optimization problem.In this study, root mean square error (RMSE) is taken as the objective function to express the error between the experimental data of the MRD and the simulation value of the phenomenal model.Under the boundary constraints of the parameter vector Q, the objective function is minimized using the IFA.The mathematical description of the constrained optimization problem is as follows: where, f is the objective function, F Exp i is the experimental data of the damping force output by the MRD, F Sim i is the damping force data obtained through simulation of the phenomenal model, N is the number of the comparing data, Q Lower and Q Upper are vectors consisted of the lower and upper ranges of the parameter vector Q, respectively, and the detailed parameter boundaries are listed in Table 2.
The dynamic characteristics of the MRD were obtained by the experimental method.Its test bench is shown in Figure 4.A RD 8040-1 MRD for the seat suspension system was made by LORD Company.The hydraulic servo shock absorber test bench made by the MTS Company was taken.The current of the MRD was regulated by the DC power source.Sinusoidal excitation at different frequencies (0.53 Hz, 0.85 Hz, 1.06 Hz, 1.60 Hz, 3.18 Hz, and 6.37 Hz) with amplitude of 15 mm was set, and different driving currents (0 A, 0.2 A, 0.4 A 0.6 A, and 0.8 A) were applied for the experiment.The force-displacement characteristic data of different frequencies and currents were obtained.The corresponding speed value can be calculated using the five-point center method [32].
To estimate the parameters of the MRD phenomenon model, MATLAB (R2020b) was used to program to minimize the RMSE between the experimental data and the simulation value.The IFA parameter settings: b 0 = 1.0, g = 1.0, a = 0.2, b = 1.0, the population size was 20, the number of iterations was 800, the identification parameters were taken by the average value of 10 simulations for every group to eliminate any spurious data.The parameter identification results of the phenomenon model are listed in Table 3.The fitting results of the predicted values and the experimental data are depicted in Figure 5. Dashed lines represent experimental data and solid lines forecast data.It can be seen from Figure 5 that the simulation results of the phenomenon model identified by the IFA are in good agreement with the experimental data, which can better simulate the nonlinear hysteretic characteristics of the MRD.

MRD controller
Heaviside step function is used in the MRD control algorithm to realize real-time tracking and semi-active control.Dyke et al. [10] used the Heaviside step function in the MRD control algorithm for the first time to solve the clipped-optimal control algorithm of the MRD for semiactive control.The mathematical expression of the Heaviside step function to control the MRD voltage v is as follows:     where v max is the maximum control voltage of the MRD, H( ) is the Heaviside step function, F c is the desired damping force estimated by the system controller, F d is the actual damping force provided by the MRD.The working principle of the MRD is described as follows: (1) When F c = F d , the desired damping force is equal to the actual damping force of the MRD.The voltage applied to the MRD remains unchanged.(2) When F c > F d , that is, the desired damping force is greater than the actual damping force of MRD, the maximum voltage is applied to satisfy the desired value.(3) In cases other than those mentioned above, the voltage value of the MRD controller is zero, that is, the MR damper has no field effect.
4 System controller 4.1 PID controller PID controller consists of proportional, integral, and derivative.Proportional enters the control state through the system and corrects in case of deviation.Integral improves the accuracy of system control by calibration when the system is static.The purpose of the derivative is to improve the accuracy and sensitivity and improve the stability of the control by correcting the deviations that occur in the system.The principle of the PID control is shown in Figure 6.
The mathematical expression of the PID control algorithm is as follows: where F c is the desired damping force of the MRD controller, e is the difference between the preset value and the actual value of the controlled object, K p is the proportional coefficient, K i is the integral coefficient, K d is the derivative coefficient.
In the digital simulation control system, The PID controller needs to be discretized.The expression of the PID controller in discrete form is as follows: where T is the sampling period, k is the sampling sequence, k = 0, 1, 2, ….
The objective of this PID control research is to control the vertical acceleration of the driver-seat.The difference between the preset value and the actual value of the controlled object is used as the input.It is controlled by changing the adjustable damping force F c .
To achieve the ideal control effect, it is necessary to set the controller parameters of K p , K i , and K d .In this study, the trial-and-error method was used to set parameters in the order of proportion, integral, and differential.According to the adjustment characteristics in Table 4, the parameters were determined by following the transformation of the control effect after changing the PID parameters.When adjusting K p , increasing the value of K p from very small to very large, and then gradually reducing it from very large to very small, and observe the change of the control effect.The optimal K p value is selected by following the transformation of the control effect curve.The adjustment of the integral coefficient K i is the same as that of the proportional coefficient K p .K i is also increased from small to large by following the transformation of control effect.If the control effect is not ideal, the influence of the derivative coefficient K d in the control strategy should be considered.The adjustment method of K d is the  same as that of K p and K i .At the same time, increasing K d from minor to large, fine-tune the K p and the Ki, and repeatedly adjust the values of the three parameters to obtain the ideal control effect.Finally, the parameters of the PID are selected as K p = 25, K i = 600, and K d = 0.2.

Fuzzy-PID controller
The seat suspension system is a nonlinear time-varying system, and the PID controller has the disadvantage that the parameters cannot be dynamically adjusted.Therefore, PID control is often combined with other intelligent control systems to form a composite controller to control the system.This research combines PID control and fuzzy control to create a Fuzzy-PID controller, the principle of which is shown in Figure 7.The fuzzy control is used to adjust the PID control parameters in real-time to meet the nonlinear characteristics of the semi-active seat suspension system to achieve better control effects [34].
In the Fuzzy-PID controller, the vertical acceleration (e) of the seat and its change rate (ec) are selected as the input of the fuzzy controller, and the input value is set to zero.e is the vertical acceleration of the seat (Z :: seat ).The output value of the fuzzy controller is the adjustment bulks of Dk p , Dk i and Dk d for the PID controller.Through the PID controller calculation, the best k p , k i and k d values are obtained, and the final output is the damping force required by the MRD.The relationships between k p , k i , k d , Dk p , Dk i and Dk d are as follows: where The initial parameters of k p0 , k i0 , and k d0 are 25, 600, and 0.2, respectively.The two inputs and three outputs of the fuzzy logic controller are described by seven fuzzy language variables, namely, negative-big (NB), negative-medium (NM), negative-small (NS), zero (ZO), positive-small (PS), positive-median (PM) and positive-big (PB).Membership functions were selected in terms of the PID setting principle, the knowledge and experience of the experts, and the literatures [21,[35][36][37][38].The membership function of triangle distribution was taken as fuzzy input and output variables.The Mamdani method was chosen for fuzzy reasoning.A total of 49 control rules were established.The output was defuzzied by the centroid method.MATLAB software was used to establish the fuzzy controller, as shown in Figure 8. Table 5 lists the fuzzy control rule table of the change value of the output Dk p , Dk i and Dk d of the fuzzy controller with the change of the acceleration error e and the acceleration change rate error ec.

Numerical simulation and results discussion
To evaluate the performance of the seat suspension system, the acceleration and the dynamic stroke of the seat suspension system were analyzed.In this paper, two different road profiles were used to simulate the forced vibration of the seat suspension system.

Random road profile
The random white noise excitation function [39] is adopted, and its time-domain expression is as follows: where f 0 is the lower cut-off frequency, which is 0.1 Hz; n 0 is the spatial reference frequency, which is 0.1 m À1 ; Z r (t) is the displacement of road irregularity, whose unit is m; G 0 is road irregularity coefficient, whose unit is m 3 /cycle.According to the ISO 8608-2016 standard, the road surface is divided into eight grades from A to H. D-grade road profile is adopted in this research, and the value of G 0 is 1024 Â 10 À6 m 3 /cycle.v 0 is the vehicle speed, with a value of 36 km/h.w(t) is the band-limited white noise signal.The random road profile of the D-grade is shown in Figure 9. Comparing the passive seat suspension system, the PID controlled MRD semi-active seat suspension system, and the Fuzzy-PID controlled MRD semi-active seat suspension system, the time-domain characteristics of the seat acceleration and dynamic stroke of the seat suspension system under random road profile are shown in Figure 10.The RMS values of the seat travel deformation (STD) and seat acceleration are shown in Table 6 (the value in parentheses is the reduction of fuzzy PID compared with PID).
From Figure 10 and Table 6, it can be seen that both PID control and Fuzzy-PID control can significantly improve the ride comfort and operational safety, and the Fuzzy-PID control is better than the PID control.Compared with the passive seat suspension system, the Fuzzy-PID control system can increase the ride comfort by 32.64%, which is less than 2% compared with the PID control.In terms of operational safety, the Fuzzy-PID control system reduces the seat suspension dynamic stroke by 9.09% compared with the passive seat suspension system and 1.64% compared with the PID control system.The peak-to-peak values of the passive suspension system, PID control system, and Fuzzy-PID control system for the seat acceleration and STD are listed in Table 7.It can be seen from Table 7 that the Fuzzy-PID control significantly affects the vibration reduction performance of the semiactive suspension system of the 3-DOF model.Compared with the passive suspension system, the fuzzy PID control system reduces the peak-to-peak value of the seat vertical vibration acceleration by 20.54%, the PID control system reduces the peak-to-peak value of the seat vertical vibration acceleration by 18.48%, the fuzzy PID control system reduces the peak-to-peak value of the seat suspension dynamic travel by 13.18%, and the PID control system reduces the peak-to-peak value of the seat suspension dynamic travel by 12.35%.Compared with the PID control system, the Fuzzy-PID control system has no significant improvement on seat acceleration and dynamic stroke of the seat suspension, which thoroughly verifies the rationality of parameter setting of the standard PID control system.
To further verify the influence of the Fuzzy-PID controller and the PID controller on seat comfort and operational safety, the energy method was used to analyze the frequency response characteristics of the suspension system of the 3-DOF system.The frequency-domain characteristic curves of the seat acceleration and STD are shown in Figure 11.It can be seen from Figure 11 that the semi-active seat suspension system controlled by PID and Fuzzy-PID has a significant damping effect at the natural frequency of the seat (0.89 Hz) and the natural frequency of the sprung mass (1.49Hz).It also shows that the Fuzzy-PID control has a slightly higher effect than the PID control.PID control and Fuzzy-PID control are not as effective as a passive suspension in suppressing vibration acceleration in the high-frequency range.Still, they have a better suppressing effect on the seat suspension dynamic stroke.

Bump road profile
The bump road profile was used to simulate the transient response characteristics of the vehicle seat suspension system traveling through the barrier.The displacement function is defined as follows [22]: where a is the height of the bump, which is 0.07 m; l is the width of the bump, which is 0.8 m; v 0 is the velocity of the vehicle, which is taken as 3.6 km/h.The stroke response curve of the STD and the seat acceleration response curve under the bump road profile are shown in Figure 12.The RMS values and the peak-topeak values of the STD and seat acceleration under bump road profile are listed in Tables 8 and 9, respectively.It can be seen from Figure 12, Table 8, and Table 9 that Fuzzy-PID control and PID control can still show superior vibration reduction performance under the excitation of the bump road profile.The Fuzzy-PID control shows better damping characteristics than the PID control, which can significantly reduce the acceleration and the STD of the seat suspension.Compared with the passive seat suspension, the Fuzzy-PID control semi-active seat suspension has a 35.03% reduction in acceleration RMS value and an 13.78% reduction in the STD.Compared with the PID control semi-active suspension system, the peak-to-peak acceleration of the Fuzzy-PID control system is reduced by 12.25%, and the STD reduced by 11.74%.Fuzzy-PID control and PID control semi-active suspension system can quickly enter a stable state than passive suspension system.Under the bump road profile, the passive seat suspension is prone to end-stop impacts, which seriously affecting the ride comfort and handling safety.Fuzzy PID control and PID control can significantly reduce the occurrence of the end-stop impacts by adjusting the control parameters to improve the ride comfort.

Conclusions
The semi-active seat suspension was integrated with a quarter vehicle model to analyze the ride comfort.The experimental method was used to analyze the dynamic characteristics of the MRD, and the IFA was used to identify the parameters of the phenomenon model to build the mathematical model of the MRD.The Heaviside step function was used as the controller of the MRD to track and feedback the desired damping force and the actual damping force of the MRD in real-time to obtain the command voltage of the MRD and realize the semi-active control of the seat suspension.Random road profile and bump road profile were used to simulate the 3-DOF semi-active seat suspension, and the seat acceleration and dynamic travel stroke of the seat suspension were compared between Fuzzy-PID control, PID control, and passive seat suspension, so as to verify the influence of semi-active control algorithm on ride comfort and operation safety.
The simulation results show that the seat suspension system of the MRD can significantly reduce the vibration response under different road profiles to the human body and reduces the impact of vibration on the physical and  mental health of drivers and passengers.Under the D-grade random road profile, Fuzzy-PID control reduce the RMS of the seat acceleration by 32.64% and the peak-to-peak value of suspension motion by 18.48%, as compared to that of the passive seat suspension, thus improving the ride comfort and reducing the probability of end-stop collision.Under the bump road excitation, compared with passive suspension, the acceleration RMS value of the semi-active seat suspension controlled by Fuzzy-PID decreases by 35.03%, and the STD of the seat suspension decreases by 13.78%.In summary, the 3-DOF seat suspension semi-active control system using Fuzzy-PID controller and Heaviside step function MRD controller can significantly attenuate seat vibration and reduce the end-stop collision, which has the characteristics of simple structure, good stability, and high robustness.It can be used in various vehicle seat suspension systems to provide passengers with a safe and comfortable driving environment, which has broad application prospects.

Fig. 5 .
Fig. 5. Comparison between the experimental data and the simulation results of the phenomenon model identified by IFA: (a) Force vs. displacement; (b) Force vs. velocity.
k p , k i and k d are the final PID controller parameters; k p0 ,k i0 and k d0 are the initial parameters before tuning; Dk p , Dk i and Dk d are the adjustment values of the proportional, integral, and derivative coefficients, respectively.According to the simulation results of the passive seat suspension system, the basic domain of the seat vertical acceleration e is [-2.5, 2.5] m/s 2 , the basic domain of the seat vertical acceleration change rate ec is [-50, 50] m/s 3 , the basic domain of k p is [20, 40], the basic domain of k i is [300, 900], the basic domain of k d is [0.15, 0.25].The input and output fuzzy domains are normalized to the interval [-1,1], then the input quantization factors are k e = 0.4 and k ec = 0.02, and the output scale factors are U p = 20, U i = 600, and U d = 0.1.

Fig. 10 .
Fig. 10.Time domain characteristics of the seat suspension system under random excitation: (a) seat acceleration; (b) seat travel stroke.

Fig. 11 .
Fig. 11.Frequency domain characteristics of the seat suspension system under random excitation: (a) seat acceleration; (b) seat travel stroke.

Fig. 12 .
Fig. 12.Time domain characteristics of the seat suspension system under bump road profile: (a) seat acceleration; (b) seat travel stroke.

Table 2 .
The constraint values of the parameter vector Q.

Table 3 .
Parameter identification results of the Phenomenon model by IFA.

Table 6 .
RMS values under random road profile.

Table 7 .
Peak-to-peak values under random road profile.

Table 8 .
RMS values under bump road profile.