Development of Matlab Simulink model for dynamics analysis of passive suspension system for lightweight vehicle

In designing suitable isolators to reduce unwanted vibration in vehicles, the response from a mathematical model which characterizes the transmissibility ratio of the input and output of the vehicle is required. In this study, a Matlab Simulink model is developed to study the dynamic behaviour performance of passive suspension system for a lightweight electric vehicle. The Simulink model is based on the two degrees of freedom system quarter car model. The model is compared to the theoretical plots of the transmissibility ratios between the amplitudes of the displacements and accelerations of the sprung and unsprung masses to the amplitudes of the ground, against the frequencies at different damping values. It was found that the frequency responses obtained from the theoretical calculations and from the Simulink simulation is comparable to each other. Hence, the model may be extended to a full vehicle model.


Introduction
The suspension system of a vehicle is primarily designed to provide good handling ability and to prevent occupants from discomfort, fatigue and health problems due to various road profiles and irregularities. However, it is difficult to achieve an optimal condition of ride comfort and handling simultaneously when using a passive suspension system [1]- [3]. The passive type suspension consists of spring and damper as the main components where no energy is added into the system compared to an active suspension system. An active suspension is capable of improving the ride and handling condition but it requires external power to function which results to added complexity in the design, cost and weight of the system. As a result, a passive type suspension system is preferred for a lightweight electric vehicle due to its simplicity and less cost involved [4]- [5].
Noise, vibration and harshness (NVH) engineering has been practiced in the automotive industry since 1980s due to customers demand for improved vehicle refinement, ride and quality. This discipline has become a crossover point for engineering practice, ergonomics and psychology where prediction of comfort levels in production vehicle is one of the objectives. This has demanded for the development of predictive models to evaluate the ride comfort level without using extensive prototyping [6]. Most of the numerical analysis studies related to vehicle vibration were based on mass-spring-damper systems. Generally, the studies were done by adopting parameters of ready-made systems and necessary modification is done to existing systems to overcome the problems. This approach is limited because of the complexity of the system. However, the use of simulation software to analyse vibration characteristics for structures has enabled designers to analyse and determine the components and suitable parameters that are able to overcome the vibration problems. The purpose of this study is to investigate the vehicle passive suspension performance and its ability to provide good driving condition by comparing the transmissibility ratios for different damping values. A mathematical model for the vehicle suspension system is developed in order to simulate its dynamic performance. This model will help to find the best suspension system parameters in providing an optimum driving comfort level.

Mathematical modeling
A vehicle development most likely commences with a mathematical model in order to optimise its dynamic and seating responses before a prototype is constructed. Furthermore, the model is used as an evaluation method that considers the sensitivity of human body to different frequencies, directions, and durations of vibration, and to obtain subjective evaluations that can contribute to design improvement [7].
A complex vehicle body supported by several elastic suspensions can be modelled as a simplified lumped parameter system. In this chapter, the mathematical equations governing the behaviour of the lumped parameters system for a lightweight solar vehicle are derived and discussed. The rigid body of the vehicle or the sprung mass is connected to other masses such as the wheels through massless springs and dampers that represent the suspension system. The wheels or the unsprung masses are connected to the ground through massless springs and dampers simulating the tyres. This simplified model is an approximation since the suspensions and the tyres actually have their own mass and natural frequencies. Apart from that, the body and the linkages can be considered to be non-rigid bodies [8].
The equations of motion of a vibrating system are usually presented in the form of a set of ordinary differential equations for a discrete system and partial differential equations for a continuous system. The fundamental law from which most vehicle dynamics analyses begin is the Newton's second law that can be applied to translational and rotational systems [9]- [10]. The law states that, for a translational systems, the sum of the external forces acting on a body in a given direction is equal to the product of its mass and acceleration in a particular direction, as represented by the following equation, where F i is the forces in the i-direction, M is the mass of the body and a i is the acceleration in the idirection.

Development of model for electric vehicle
Various components inside a solar vehicle are assembled together with the body to become a unit vehicle and moving together which makes it possible to represent them as one lumped mass at the centre of gravity (CG) of a system with appropriate total mass and inertia properties. Fundamentally, all road vehicles have common ride isolation properties that use primary suspension system at each wheel to support a sprung mass. This system dynamic behaviour is the first level of isolation from the unevenness of the road. An isolation system can be passive or active depending on whether there is external power to function as an isolator. Normally, a passive isolator consists of a resilient member (stiffness) such as metal spring, elastomer spring or pneumatic spring, and energy dissipator (damping). An active control isolator applies actuator such as hydraulic piston, piezoelectric device or electric motor to provide a force to the system whose vibration properties are to be changed [11]- [12]. Even though an active control system performs better than a passive isolator, it is costly and less reliable than the latter.

Quarter vehicle model
The essential dynamics can be represented by two DOF quarter vehicle models as shown in figure 1. The sprung mass is supported by a primary suspension that is placed on the unsprung mass of the axle. The equations of motion for this model represent the motions of unsprung mass, M u and sprung mass, M s according to their free body diagram, and solving the forces acting on the bodies at equilibrium. The motion of each mass is represented by the displacement, velocity and acceleration, in the vertical direction, represented by Z , Z and Z respectively from the static equilibrium position. The input displacement due to road profile is given by Z r while Z u and Z s are displacement of M u and M s respectively.  Firstly, the equation of motion for each mass is developed when the system is at equilibrium. The free body diagram for mass, M u as in figure 2 is used to develop the equation by assuming that the system is linear and Z r is greater than Z u , and the latter is greater than Z s . The stiffness and damping constants of the shock absorber are K s and C s respectively while K t and C t represents the stiffness and damping constants of the tyre respectively. The equation of motion for the M u is given by where F u is the force acting on M u . Similarly, the equation of motion for M s is developed based on the free body diagram as shown in figure 3.
where F s is the force acting on M s . Therefore, the equations of motion for the two DOF model as in figure 1 can be written in the form are vectors representing accelerations, velocities, displacements and external forces respectively. The matrix equation can be simplified by neglecting some parameters such as the weight of the masses due to static equilibrium position and also the damping constant of the tyre which is usually very small [10]. The simplified equation is given by  There are different types of road excitation that can be applied to the system such as step, sine wave and speed bump signals. The MATLAB Function as in figure 4 allows a sinusoidal excitation to be developed and used as input to the model. The mathematical expression of the excitation is given by: for 0 ≤ u < 10 y = 0.01sin(λ*u) for 10 ≤ u ≤ 30 (9) y = 0 for u > 30

Results and discussion
The frequency responses related to both M s and M u are plotted as shown in figure 7 and figure 8 for a system with K t = 10 kN/m, K s = 2.5 kN/m (K t = 4K s ) and M s = 150 kg, M u = 10 kg (M s = 15M u ). The values of the masses, M and stiffness constants, K are based on an actual three-wheeled solar vehicle previously developed. The plots show the amplitude ratio of the output and input against the nondimensional frequency λ * = λ(M s /K s ) 0.5 obtained using different values of C s . In figure 7, two natural frequencies can be observed for C s = 0 withhigher peak being at approximately 4 Hz. It is also noticed thatthe peak is highest for the case of C s = ∞ where the whole system is considered to be rigid. All of the peaks occurred at frequencies below 10 Hz except for the second natural frequency of Cs = 0 which has its highest peak at 20 Hz. The lowest peak of the ratio happens when the optimum damping value of C s-opt was used.  In figure 8, it is observed that the highest peaks occurred when C s = 0 and C s = ∞. All of the peaks are approximately at 8 Hz except for the case of C s = ∞ in which the recorded peak is at 20 Hz. The ratio of amplitudes Z u0 to Z r0 (for M s ) is higher than the ratio of amplitudes Z s0 to Z r0 (for M u ) especially at frequencies above 10 Hz for C s = C s-opt , 2C s-opt and 0.   In figure 9, the ratio of the acceleration s Z to Z r0 increases with decreasing C s for frequencies below 6 Hz and between 16 Hz and 25 Hz while the ratio increases with increasing C s for frequencies between 6 Hz and 16 Hz and from 25 Hz upwards. In figure 10, the ratio of the acceleration u Z to Z r0 increases from 4 Hz up until 12 Hz with increasing C s but the ratio decreases with increasing C s beyond that. The optimum value of C s can be found by trying to keep the acceleration as low as possible in a large field of the frequencies. In the above cases, the value of C s-opt is obtained from  Figure 9 also shows that this value of damping is effective in maintaining low acceleration in a wide field of frequencies. Figure 11 and figure 12 show the frequency responses of the ratios between amplitudes of displacement M s , Z s and displacement M u , Z u to the amplitude of the displacement of the input signal or the ground, Z r for different values of C s using MATLAB Simulink. Both set of graphs exhibit similar output patterns with figure 7 and figure 8. In figure 11, the ratio between Z s and Z r increases with decreasing C s at 4 Hz with one significant peak at 8 Hz when C s = ∞. Similarly as in the previous cases, two peaks are obtained if C s = 0. However, the magnitude of the ratios are different  figure 12, the peaks for C s = C s-opt , 2C s-opt and 0.5C s-opt are observed to be not so significant compared to responses in figure 8. However, the amplitude of the ratios are of the same order of magnitude after 12 Hz when C s = 0.5C s-opt > C s = C s-opt > C s = 2C s-opt .     Figure 13 and figure 14 show the frequency responses of the ratios between acceleration s Z to Z r and acceleration u Z to Z r for different values of C s using MATLAB Simulink. The order of amplitude for the ratios of s Z to Z r, obtained below 4 Hz, between 4 Hz and 24 Hz, and above the 24 Hz are almost similar, as shown in figure 13 with the responses in figure 9. Similarly in figure 14, the ratios of u Z to Z r increase with increase in frequency except for C s = ∞.

Conclusion
A quarter-vehicle model for an electric lightweight vehicle was developed in order to analyse the performance of the vehicle's passive suspension system in isolating disturbances from the roughness of the road. Differential equations were used to plot the frequency responses of the system in terms of the ratios of the displacements M s and M u to the displacement of input Z r representing the road disturbance. The ratios of the accelerations for M s and M u to the input Z r are also plotted against the non-dimensional frequency. Simulink model developed from the differential equations were also used to analyse the frequency responses of the system with different values of the damping constant, C s . It can be seen that the most suitable value of the damping is C s = C s-opt since in all analyses the responses are minimum at most level of frequency range as shown. It is also important to show that the frequency responses obtained from the calculation and simulation using Simulink are comparable and agree with the frequency responses given by Genta, 1997 [8]. However, the differences that can be observed from the graphs might be due to the inaccuracies in the characteristics of the input signal applied to the simulation.

Future scope
A full-vehicle model for a lightweight electric vehicle will be developed using Simulink to analyse the vibration responses using suspension system parameters similar to the off-the-shelf components. Furthermore, the input signal can be easily modified based on actual road surface roughness and the vibration responses from the involvement of more masses can be obtained easily.