Matlab work
Dynamic Model of an Automobile Suspension
Jake Bailey, CSU Fullerton
Spring 2016
Contents
1 Introduction 1
2 Modeling 2 2.1 Total System Energy . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2
2.1.1 Potential Energy . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2 2.1.2 Kinetic Energy . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 2.1.3 Dissipative Terms . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4
2.2 Lagrangian and the Equations of Motion . . . . . . . . . . . . . . . . . . . . . . . . . 4 2.3 Linearizing a Nonlinear Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
2.3.1 Small Angles . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7
3 Conclusion 7
1 Introduction
In this document, we shall use the methods of Lagrangian Mechanics to develop a complete set of equations of motion for a simplified version of an automobile suspension. The suspension we will model includes one half of a complete vehicle, with half a chassis suspended upon two individual wheels. Figure 1 shows this graphically.
Figure 1: A simplified dynamic model of half of a car’s suspension, including tire stiffness and viscoelastic damping. For our purposes, we shall ignore the effects of the tires and treat them as rigid bodies.
As can be seen in the figure, we have two spring stiffnesses and two damping coefficients to deal with, along with the two ”unsprung” masses and the mass of the vehicle. These are the parameters on which we shall base the dynamic model of the vehicle.
1
Before we begin modeling, we shall make a few assumptions about the performance of the system.
• Both wheels maintain contact with the ground throughout the time period of interest
• Motion of the vehicle is constrained to the plane (i.e. there are only 3 degrees of freedom)
• Neither wheel is allowed to slip along the ground as the vehicle moves in the x direction
2 Modeling
2.1 Total System Energy
2.1.1 Potential Energy
To use the Lagrangian approach to arrive at our equations of motion, we must first define the Lagrangian. As you’ll recall, the Lagrangian for systems with conservative potential energies can be written as
L = T − U
where T represents the total kinetic energy in the system in terms of the generalized coordinates, and U represents the total potential energy. We will first concern ourselves with the potential energy of the system. By examination of the system, we see that there are only three sources of potential energy in the model:
1. Energy stored via compression of the front spring
2. Energy stored via compression of the rear spring
3. Gravitational potential energy due to the height of the chassis’ mass
You’ll notice that we ignore the potential energy due to the height of the “unsprung” masses of the wheels, due to our assumption of constant contact with the ground. While there may be some potential energy stored by these masses, it will ultimately be included in our analysis of the base excitation of the system, which is due to the unevenness of the ground.
Items 1 and 2 of the potential energy can be expressed in similar form:
U1 = 1
2 kfy
2 f
U2 = 1
2 kry
2 r
The third term in the total energy is simply the familiar gravitational potential energy,
U3 = Mg(yf + yr)
2 ,
where we take the height of the chassis’ mass as the average of the front and rear suspension heights. Finally, we can write the total potential energy of the system, Utotal, as
Utotal = U1 + U2 + U3
Utotal = 1
2 kfy
2 f +
1
2 kry
2 r +
Mg(yf + yr)
2 (1)
2
2.1.2 Kinetic Energy
Moving on to the kinetic energy, we will see that again we have three possible sources of kinetic energy, although their equations are slightly more complicated. These sources are
1. Kinetic energy of the front wheel
2. Kinetic energy of the rear wheel
3. Kinetic energy of the chassis mass
Since we are considering full plane motion, we will have energy terms for these three masses in all three coordinates x, yf, and yr. However, a few of these will simplify out.
We can write the total kinetic energy of the front wheel as
Tf = 1
2 mfẋ
2 + 1
2 mf ẏ
2 f +
1
2 mf ẏ
2 r
However, we know that, by definition in our model, the front wheel has no motion due to yr, and thus the third term in equation X is zero. Thus, the front wheel’s kinetic energy simplifies to
Tf = 1
2 mfẋ
2 + 1
2 mf ẏ
2 f
Similarly, for the rear wheel,
Tr = 1
2 mrẋ
2 + 1
2 mrẏ
2 r
Next, we turn our attention to the chassis mass of the vehicle. While we could find the kinetic energies of the body due to our current coordinates, it is slightly more intuitive to employ a trans- formation of coordinates into a more natural set for the chassis. We note that our coordinates of interest for this simulation are the pitch and bounce of the vehicle, as discussed in the introduction. Further, we note that these two can easily be obtained from our current coordinates via the simple transformations
ybounce = (yf + yr)
2
Θ = sin−1( yr − yf
L )
The final coordinate for the chassis mass, x, remains unchanged. We will also need the velocities of these transformed coordinates, which you can verify are as shown below.
ẏbounce = (ẏf + ẏr)
2
Θ̇ = (ẏr − ẏf )
L √
1 − ( (yr−yf ) L
)2
With the coordinate transformation in hand, all that remains is to write the kinetic energies of the chassis mass in the typical fashion,
TM = 1
2 Mẋ2 +
1
2 Mẏ2bounce +
1
2 JM Θ̇
2
3
noting that JM represents the mass moment of inertia of the chassis. Substituting our coordinate transformations for the output coordinates, we obtain the final kinetic energy for the chassis mass.
TM = 1
2 Mẋ2 +
1
2 M
(ẏf + ẏr)
4
2
+ 1
2 JM
ẏ2r − 2ẏrẏf + ẏ2f (L2 + 2yfyr − y2f − y2r )
We combine the kinetic energies of the three bodies to arrive at the system’s total kinetic energy.
Ttotal = 1
2 mfẋ
2 + 1
2 mf ẏ
2 f +
1
2 mrẋ
2 + 1
2 mrẏ
2 r +
1
2 Mẋ2 +
1
2 M
(ẏf + ẏr)
4
2
+ 1
2 JM
ẏ2r − 2ẏrẏf + ẏ2f (L2 + 2yfyr − y2f − y2r )
(2)
2.1.3 Dissipative Terms
With the system’s kinetic and potential energies in hand, we now turn our attention to the energy dissipation due to the front and rear dampers. As you’ll recall from our readings, we can construct a Rayleigh dissipation function, R, that includes the effects of all the dampers in the system, and include it in our Lagrangian model as a part of the general forces applied to the system. Recall that for low velocity, viscous damping, R takes the form
R = − 1
2
n∑ i=1
n∑ j=1
cijq̇iq̇j
Since we only have two dampers, our dissipation function is relatively simple:
R = − 1
2 cf ẏ
2 f −
1
2 crẏ
2 r
Recalling that the generalized force due to the Rayleigh dissipation function is
Fi = ∂R
∂q̇i
We find that the forces due to our dampers are simply the familiar linear terms in the respective velocities.
Ff = − 1
2 cf ẏ
2 f
Fr = − 1
2 crẏ
2 r
2.2 Lagrangian and the Equations of Motion
We now have the requisite material at hand to produce the Lagrangian of the system, and thus the equations of motion for the entire system. As noted in the previous section, the Lagrangian is simply L = T − U. Bringing together our complete terms for the potential and kinetic energies of the system, we arrive at the system’s Lagrangian:
L = Ttotal − Utotal
4
(3) L =
1
2 mfẋ
2 + 1
2 mf ẏ
2 f +
1
2 mrẋ
2 + 1
2 mrẏ
2 r +
1
2 Mẋ2 +
1
2 M
(ẏf + ẏr)
4
2
+ 1
2 JM
ẏ2r − 2ẏrẏf + ẏ2f (L2 + 2yfyr − y2f − y2r )
− ( 1
2 kfy
2 f +
1
2 kry
2 r +
Mg(yf + yr)
2 )
We can simplify this expression for the Lagrangian significantly by combining like terms, and introducing a few new ones. First, we can see that the term (L2 +2yfyr −y2f −y
2 r ) appears frequently
in the Lagrangian. For convenience, we rename this term as Le, a sort of effective length between the wheels. Since we will have need of it shortly in deriving the equations of motion, we also find that the time rate of change of the effective length can be written as L̇e = 2 (yrẏf + ẏryf )−2yrẏr−2yf ẏf . Substituting this new length and gathering like terms, we arrive at the simplified Lagrangian:
(4) L =
(M + mf + mr)
2 ẋ2 + (
M
8 +
JM 2Le
+ mf 2
)ẏ2f + ( M
8 +
JM 2Le
+ mr 2
)ẏ2r
+ ( M
4 −
JM Le
)ẏrẏf − 1
2 kfy
2 f −
1
2 kry
2 r −
Mg(yf + yr)
2
Next, we begin to derive the equations of motion from the Lagrangian. Starting with the equation of motion for the x coordinate, we take the appropriate partial derivatives as
d
dt
∂L
∂ẋ −
∂L
∂x =
∑ Fx
∂L
∂ẋ = (M + mf + mr)ẋ
d
dt
∂L
∂ẋ = (M + mf + mr)ẍ
∂L
∂x = 0
and arrive at the simplified equation of motion:
(5)(M + mf + mr)ẍ = ∑
Fx
This unsurprising result tells us simply that the acceleration of the total mass of the vehicle, and thus the vehicle itself, depends solely on the forces applied in the direction of the x coordinate. The equation of motion in the yf and yr coordinates are a bit more involved, as we shall see.
Again taking the appropriate partial derivatives,
∂L
∂ẏf = (
M
4 +
JM Le
+ mf )ẏf + ( M
4 −
JM Le
)ẏr
(6) d
dt
∂L
∂ẏf = (
M
4 +
JM Le
+ mf )ÿf + ( M
4 −
JM Le
)ÿr − JM L̇e(ẏr − ẏf )
L2e
Noting here that the third term in this derivative appears due to the time variant nature of our effective length, Le.
Continuing our work with partial differentiation, we move on to those with respect to yf
∂L
∂yf = −kfyf −
Mgyr 2
5
Before completing the equation of motion, we recall that our forces due to damping were included as part of the generalized force acting in each of the coordinate directions. Thus, we may write the generalized force in the yf direction as a sum of the known damping force and any external forces as such: ∑
Fyf = ∑
Fexternal − cf ẏf Now, we may write the full equation of motion for the yf direction.
d
dt
∂L
∂ẏf −
∂L
∂yf =
∑ Fyf
(7) ( M
4 +
JM Le
+ mf )ÿf + ( M
4 −
JM Le
)ÿr
− JM L̇e(ẏr − ẏf )
L2e − (−kfyf −
Mgyr 2
) = ∑
Fexternal − cf ẏf
Which, after a bit of rearrangement, gives us our final (non-linear) form of the equation of motion in the yf direction.
(8) ( M
4 +
JM Le
+ mf )ÿf + ( M
4 −
JM Le
)ÿr
− JM L̇e(ẏr − ẏf )
L2e + cf ẏf + kfyf +
Mgyr 2
= ∑
Ffexternal
By a nearly identical set of arguments, we may find the equation of motion for the yr direction. Since it would be highly repetitive to demonstrate the derivation of said equation, we will not labor it here. We will, however, call your attention to the fact that the cross-coupled term in the velocities (i.e. that term involving the time derivative of the effective length, L̇e) has had its sign switched in the yr equation.
(9) ( M
4 +
JM Le
+ mr)ÿr + ( M
4 −
JM Le
)ÿf
+ JM L̇e(ẏr − ẏf )
L2e + crẏr + kryr +
Mgyf 2
= ∑
Frexternal
Thus, we note that equations 5, 8, and 9 represent a complete set of equations of motion for our three degree of freedom model of an automobile suspension. The astute reader will notice that these equations are both highly coupled, and non-linear. While non-linear equations such as these pose little trouble to solution via modern numerical methods, they do pose a problem to us mere mortals in that we cannot use any of the methods we have developed for investigating and interrogating linear, multiple degree of freedom systems.
2.3 Linearizing a Nonlinear Model
To aid us in our system inspection, design, and selection of parameters for functional improvements, we shall linearize the previous equations of motion about a chosen operating point. There are many methods which would allow such a linearization, such as the Jacobian technique. However, the nature of our problem allows a much simpler method to be employed by adding an extra assumption to our initial set.
6
2.3.1 Small Angles
By restricting our attention to the performance of the system within a small range about the “origin,” or our equilibrium point, we can further simplify our model. First we note that our equilibrium point is specified completely as
x = 0,
ybounce = 0,
Θ = 0.
By restricting our attention to small angles about this equilibrium point, we can take advantage of the fact that at small angles,
sinΘ ≈ Θ.
Substituting this into our coordinate transform equation, we arrive at the following result:
Θ = (yr − yf )
L ,
Θ̇ = (ẏr − ẏf )
L ,
Θ̇2 = (ẏ2r − 2ẏrẏf + ẏ2f )
L2
This leads us to the result of an effective length, Le, which has a constant value of L 2. Most
importantly, we see that L̇e = 0, which shall simplify our model significantly. Substituting these results into our previous equations of motion lets us write:
(10)(M + mf + mr)ẍ = ∑
Fx,
(11) ( M
4 +
JM L2
+ mf )ÿf + ( M
4 −
JM L2
)ÿr
+ cf ẏf + kfyf + Mgyr
2 =
∑ Ffexternal,
(12) ( M
4 +
JM L2
+ mr)ÿr + ( M
4 −
JM L2
)ÿf
+ crẏr + kryr + Mgyf
2 =
∑ Frexternal,
which we can see are a set of three linear, coupled differential equations.
3 Conclusion
This linear set of equations can be used to interrogate the basic properties of our system, such as eigenvalues (natural frequencies), eigenvectors (mode shapes), and stability. The linearized equa- tions can also be used for lightweight simulation. In particular, MATLAB and Simulink offer very straightforward procedures for time domain simulation of systems represented in state space form, which is an excellent place to start as you work toward the design goals of your project. Note, however, that the full nonlinear model developed in the previous section will still be necessary as you look to validate your design results on that model.
7