ME513VehicleDynamics LectureSlides 2021

ME 513
VEHICLE
DYNAMICS
Prof. Dr. Y. Samim ÜNLÜSOY

Mechanical Engineering
Department, Middle East
Technical University, Ankara
Prof. Dr. Y. Samim
ME 513 Vehicle Dynamics 1
Ünlüsoy
 CH I-INTRODUCTION CONTENTS
 CH II-PERFORMANCE (ME 466)
 CH III-VEHICLE HANDLING
 Pneumatic Tire Characteristics
 Plane Motions of Vehicles
 Articulated Vehicles
 Four Wheel Steering
 Roll Centers and Roll Axis
 Higher Order Models
 CH IV-VEHICLE RIDE
 Vibrational Characteristics & Modeling
 Road Modeling
 Response of Linear Models to Random
Inputs
 Passive Suspension Optimization
 Active Suspensions
Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 4

CH I - INTRODUCTION
 High speed vehicles are the main
concern of the subject of Vehicle
Dynamics.
 There exists many different
classifications for high speed
vehicles which are capable of
travelling at maximum speeds of
over 100 kph.

Classification
 A general classification will include :
 Ground Vehicles
- Guided and
- Nonguided
 Fluid Vehicles
- Marine,
- Air, and
- Space-crafts.

Road Vehicles
 In this course, dynamics of the road
vehicle which is in the group of nonguided
ground vehicles, will be considered.
 Road vehicles move on road surfaces
and, therefore, their motions are affected
by
 road irregularities,
 road and environmental conditions,
and
 commands from the driver.

Definition
 Vehicle Dynamics is formally
defined as
“The study,
measurement,
analysis, and
prediction
of the motion of vehicles.”

Areas
 With respect to road vehicles, the
subject of vehicle dynamics is
customarily subdivided into the three
areas of
 performance,
 handling, and
 ride.

Performance
 The performance of a vehicle refers to
its ability to
Subject of  accelerate, decelerate, or travel at a

ME 466 constant velocity on level or
sloping road surfaces in response
to driver commands.
 It includes the study and reduction
of fuel consumption of the vehicle
for a specified driving cycle.

Handling
 Handling of a vehicle refers to its
directional control and stability
behavior in response to
 steering commands, or
 environmental inputs.

Ride
 Ride (or Ride Comfort) is related to

the vibratory activity of the vehicle,
which is excited mainly by the road
irregularities among other inputs.

Approach
 In real life the vehicle is simultaneously
subjected to
 road irregularities,
 control inputs, and
 A vehicle, however, is a very complex
dynamic system and it is difficult to
model and then interpret its response to
these simultaneous inputs.

Decoupling Areas
 Thus, in most of the studies, it is common
to completely decouple performance,
ride, and handling, e.g.,
 study the ride motion in the absence
of any control inputs (but not the
environmental inputs).
 similarly, handling is easier to study
if it is assumed that the vehicle is
driven on a flat road with no road
irregularities.

ME 513
VEHICLE
DYNAMICS

Prof. Dr. Y. Samim
Ünlüsoy
CH II
VEHICLE HANDLING
Prof. Dr. Y. Samim Ünlüsoy ME 465 Automotive Engineering 2

VEHICLE HANDLING
 Vehicle handling is an area of vehicle

dynamics related to the motion of a vehicle
under the action of
 directional commands and/or
 A vehicle system consists of two main parts;
 driver and
 vehicle.

Man-vehicle System
Road Surface
Environment
Steering wheel Position
Accelerator pedal Velocity

DRIVER VEHICLE
Brake Pedal Acceleration
Feedback
Feedback

Vehicle System
Road Surface
Environment
Steering wheel Position
Accelerator pedal Velocity

VEHICLE
Brake Pedal Acceleration

LATERAL ACCELERATION

Fc
Lateral Acceleration
U
 When a vehicle is cornering, it will be
ac experiencing a centripetal acceleration,
where
2  U : forward speed of the vehicle

U
ac   and
R
 R : radius of curvature of the turn.
Fc  m ac
 The cornering vehicle will then be
U2 subjected to a centrifugal force, Fc,
m
R acting at its center of gravity, away from
the center of curvature.

WaL U
aL
 It is customary, in vehicle dynamics to
use the nondimensional lateral
ac 1 U 2 acceleration which is defined as the
aL  
g g R ratio of centripetal acceleration to
gravitational acceleration
Fc  m ac  Thus, centrifugal force can be written
 mgaL in terms of weight and lateral
 WaL acceleration.

LOAD TRANSFER

Wa L
Load
Transfer
h
Fco
W
F
center of During
curvature
t/2
ci
x Cornering
t/2
W Wi
o
t
 During cornering, centrifugal force acts at the center of

gravity, away from the center of curvature. Thus, there
will be a load transfer from inside wheels (closer to the
center of curvature) to outside wheels. This results in
outer wheels carrying larger normal (vertical) loads
compared with those of inner wheels.

Wa L
Load
W h
Transfer
Fco F
ci
center of
curvature
During
x
W
t/2 t/2 Cornering
o Wi
t
On straight road W
W  Wi  Wo Wi   W
2 Thus, the outer wheel
On curves :
W load is higher than
t Wo   W
W  hWaL  tWi 2 the inner wheel load
2 during cornering.
h
t
W  hWaL  tWo W  WaL
2 t
Wa L
Load
W h Transfer
center of
Fco F
ci curvature During
x
W
t/2 t/2
Wi
Cornering
o
t
Lateral acceleration at A static Stability Safety Factor (SSF)

which the vehicle may be used to assess basic rollover
starts to turn over. Wi = 0 stability.
 h
ΔW = WaL =
W
t 2 SSF = a to = t
At impending turn L 2h

over, inner wheel load
a to = t
becomes zero. Thus L 2h

Wa L
Load
W h Transfer
center of
Fco F
ci curvature During
x
W
t/2 t/2
Wi
Cornering
o
t
Estimate lateral acceleration at which the vehicle

starts to turn over, i.e. static stability factor.
to 1340
SSF = aL =
h=580 mm 2  580 
t=1340 mm
= 1.2

Load Transfer During Cornering
 Note that this analysis is based on rigid
suspensions and assumes equal loads on
front and rear wheels.
 In reality these conditions are almost
never satisfied. Tire and suspension
deformations result in lowered roll over
thresholds. Thus, impending turn over
will start with one wheel leaving the road
first and this may occur at considerably
lower lateral accelerations than the result
obtained in the example.

 Realistic values considering suspension and tire

deformations are :
 0.8 to 1.0 [g] for automobiles, and
 0.4 to 0.5 [g] for commercial vehicles.
Since, on dry roads the road-tire adhesion is ≥ 0.8, the
automobiles are more likely to slide laterally rather than
roll over. Commercial vehicles, on the other hand, are
more likely to overturn rather than slide sideways.
Thus, vehicle dynamics tests involving severe
maneuvers are performed on wet or icy surfaces.


TIRES
(Supplementary Information)

Nominal Section
Width (NSW)
Tire Dimensions
Nominal
 Three main dimensions of
Section pneumatic tires are
Height
(NSH)
indicated on the figure.
 Note that the NSW is
always specified in tire
Nominal Rim
designations, while NSH
Diameter (NRD) is never specified.
 Thus, NSH can be calculated if

NSH
aspect ratio, which is usually 
expressed as a percentage, is NSW
known.
Tire Designation
235 / 55 R 17 80 W Speed Rating [kph]

(See Table)
NSW [mm]
. NRD ["]
165 (x100) Load Rating [kgf]
175 . 13
14 (See Table)
185 70 R for radial
195 65 D for diagonal 15
. 60 B for bias belted
.

TIRE
CONSTRUCTION
 The pneumatic tire
consists of three primary
components.
- Tread
- Carcass
- Bead & bead wires

Tread
Breaker belts
Carcass plies

Circumferential
OBSERVATIONS
 When the direction of the cords

approaches the circumferential
direction,
 directional control and
stability improves, and
 ride gets harsher.

Circumferential
Radial OBSERVATIONS
 If, on the other hand, the cords

are placed radially (at right
angles to the circumference)
 softest ride is obtained,
but
 directional control and
stability are deteriorated.

CROSS PLY TIRES
 As a result of these experiments,
cross ply tires were introduced
in the early 1920s.
 In this type of construction, the
carcass is made up of two or
more plies (layers) of fabric
which has all its cords in any
one ply running in one direction
only.
 The adjacent layers are placed at
opposite cord directions.

Bias (crown) angle CROSS PLY TIRES
 To obtain the compromise
between the two extremes,
layers are placed such that the
cords make an angle (bias or
crown angle) of 45 degrees with
the circumferential line of the
tire.
 The bias angle has now been
reduced to about 40 degrees. In
racing car tires bias angles
down to 25 degrees are
encountered.

RADIAL PLY TIRE
 The idea of dividing the tasks of
providing soft ride and good
directional control and stability
between the two parts of the carcass
led to the development of the radial
ply tire.
 In this construction the cords in the
plies run from bead to bead across
the circumferential line resulting in a
very flexible sidewall. A very soft ride
is thus assured.

RADIAL PLY TIRE
 Directional stability and

control are provided by a
number (a minimum of two) of
belts called breaker belts
running round the
circumference, but not
extending to the side sides,
beneath the tread.

RADIAL PLY TIRE
 The cords of
the breaker
belts are
slightly
diagonal with
bias angles of
about 20
degrees.
Body plies
Breaker belt

RADIAL PLY TIRE
 The zone beneath the tread is firmly

restricted by the breaker belts and the
lateral stiffness is increased while the
tire is soft in the vertical direction.
 The overall result is the reduced
distortion as the tire passes through
the contact area with the road surface.
There is less friction between the cords
of the plies and less heat is generated.

Cornering Force Characteristics
of Pneumatic Tires

Tire Cornering Force
Characteristics
 The motion of a vehicle under the
driver control and aerodynamic inputs
depends mainly on the forces applied
through the tire-road contact.
 The study of vehicle handling,
therefore, requires a thorough
understanding of the force and
moment characteristics of rolling
tires.

Tire Cornering Force
Characteristics
 The dynamic behavior of the tire,
however, is so complex and
nonlinear that no adequate general
mathematical model has so far
been obtained.
 Experimental methods have been,
therefore, used extensively in the
study of tire behavior.

Finite Element Models
Models for dynamic

analysis (i.e.
cornering forces, SAT,
etc.) are still in their
infancy and the
simulations results
are yet neither reliable
nor accurate enough.
FEM Stress analysis

Cornering
 If a moving vehicle is steered
away from the straight-ahead
position, it will still tend to go
straight on because of its
inertia.
 The tires, on the other hand
will try to pull the vehicle
towards the direction they are
steered.

Direction of Motion
Cornering
 Therefore, during cornering, there
will always be a difference between
 the direction of motion of the
vehicle and
 the direction to which the
wheels are steered.
 This can only be possible if the tires
are deformed such that the part in
contact with the ground is forced to
follow the direction of motion.

Direction of
Motion
Tire Definitions
Plane
Slip  Contact patch (tire footprint) is
Angle
the part of the tire in contact with
the road surface.
 The angle between the direction
Contact
Patch
of motion and the plane of the
Fc wheel is called the slip angle.
 The resultant of the elastic
forces due to the deformation of
the tire in the contact patch, at
right angles to the plane of the
tire, is called the cornering force
(or lateral force).

Direction of
Motion Definitions
Slip
 The distance between the line of
Angle action of the lateral force and the
Pneumatic projection of the wheel axis of
Trail, t
rotation on the road is called the
pneumatic trail.
Fc
 The moment exerted by the
lateral force about the center
Self point of the contact patch
Aligning tending to reduce the slip angle
Torque and thus turning the wheel
towards the straight-ahead
position is called the self-
aligning torque.

Rolling Tire
 Underinflated tire rolling

under load over a
flywheel at a slip angle*.
 Note the deformation of
the part of the tire in
contact with the flywheel.
* Tire testing apparatus, Automotive

Laboratory, Mechanical Engineering
Department, Middle East Technical
University.

Choice of Reference Systems
rear z
y rear
front
x
x
z front y
SAE reference system ISO reference system
Note that all rotations are clockwise when viewed from

the origin in the positive directions of the axes !
SAE
reference
frame Tire Forces
& Moments
Fx : Longitudinal
(drive/brake force)
Fy : Lateral (cornering
force)
Fz : Tire load
My : Pitch moment
Mz : Yaw moment
 : Camber angle
w : angular speed

Direction Direction
of motion a of motion a
Fc
Self-
Fc Aligning
Torque
Self
Aligning
Torque
Note that the sign of SAT will be positive in the configuration shown,
i.e., trying to reduce a, which is the usual case at lower speeds. At
higher speeds, Fc may move in front of the center of tire contact
patch center, and thus SAT becomes negative.

Cornering Force Characteristics
 The cornering behavior of a tire depends on
many factors. The three main parameters,
however, are :
 magnitude of the cornering force,
 vertical (normal) tire load, and
 slip angle.
As there are three parameters of interest,
there are different ways of presenting the
tire cornering force characteristics.
Note : In all these presentations, parameters other
than slip angle and normal tire load (camber, traction,
etc.) are assumed to be zero.

35
Tire : 12 R 22.5 inflated to 8.0 bar
30 W=52.5 kN
W=38.5 kN
Cornering Force (kN)
25
20
W=28.5 kN
15
W=14.5 kN
10
5 W=9 kN
0
0 2 4 6 8 10 12
Slip Angle (o)
Cornering (Lateral) Force vs. Slip Angle

Normal Load as a Parameter
Carpet Plot
4000
155 R 13 4500 N Vertical Tire Load

165 kPa 4000 N
3500 N
3000
3000 N
Cornering Force [N]
2500 N
2000 N
2000
1500 N
7° 8°
6°
1000 5°
4°
3° Slip Angle
2°
1°
0

30
Fc kN
20 Cornering Force
15
vs. Normal Load
Slip Angle as a
a [o ] Parameter
10
5
Now obselete !
W kN
0
0 2 4 6 8

3-D Plot

Cornering Stiffness
35 Observe that at
W=38.5 kN small slip angles
30
(below about 4o)*
cornering force
25
is almost a linear
20 function of the
slip angle.
15
Fc
10
* 4o slip angle, as
a the limit of linear
5
tire characteristics,
seems to apply to a
0
0 2 4 6 8 10 12 large majority of
o
Slip Angle ( ) pneumatic tires!
20
Cornering Stiffness
W=38.5 kN  This linear relationship can
15
be described by a constant
10
Fc coefficient relating cornering
force to slip angle. This is the
a
5 linear tire model.
0
0 2
o
Slip Angle ( )
4 6 Fc = Cs .α
The constant, Cs , is defined as
the derivative of the cornering ∂Fc ΔFc
force with respect to slip angle Cs = α=0 
∂α Δα
evaluated at zero slip angle and
is called the cornering stiffness.

Direction Direction
of motion of motion
Sign
a a
+Fc -Fc
Convention
Steady state
Top view
right hand turn
For a steady state right hand turn: Slip angle is

positive for a counterclockwise rotation of the central
plane of the tire, negative otherwise.

Sign Convention
+Fc
A steady-state right-hand turn represents
positive centripetal acceleration (to the
ac right) which requires a positive lateral
force from the tires. Thus :
A negative slip angle is required for a
positive lateral force.
As a result :
Cs < 0 Cornering stiffness is by definition a

negative quantity. In vehicle handling
simulations, a negative cornering
stiffness value must be entered .

Sign Convention
Note that, with the

sign convention,
the actual variation
of the lateral force
versus slip angle
should be as in the
figure, i.e., positive
Fy for negative a.

Sign Convention
You may come

across some
publications,
however, that take
this into account
in the vehicle
equations. Thus,
they take positive
Fy for positive a.

Secondary Parameters

 The cornering behavior of tires are
influenced by various parameters and the
cornering characteristics must be modified
to take their effects into consideration.
 Here, the variation of basic tire cornering

characteristics with varying
 tire inflation pressure,
 camber angle, and
 tractive or braking forces
will be examined.

Effect of Inflation
Pressure
In the usual working range
of tire pressures (around
recommended value),
increasing pressure
improves the cornering
force by increasing lateral
stiffness of the tire.
This improvement continues at higher pressures,

but at a lower rate because of the reduction in the
size of the contact patch.

g g
Camber Angle
 The deviation of the

central plane of a tire
from the vertical plane
is known as camber
and the angle between
the two planes is called
Center of curvature the camber angle.
 The angle is positive when the tire leans away from the
vehicle, and negative when it leans towards the vehicle.

Camber Angle
 Camber angle given to the wheels

during assembly is called static
camber.
 Dynamic camber, on the other hand,
may result from the roll motion of the
body during cornering or from wheel
travel.

Camber
Angle
The 1960 Miliken MX1 Camber Car showing a large negative camber

To achieve highest possible side forces all modern racing cars run
with negative camber at all wheels. Formula 1 cars work with a
negative camber of over -4 degree at the front. The wider rear tires
can withstand less camber (appr. -2 degree as a maximum)

g g
Effect of
Camber Angle
When the camber angle is
positive, cornering force
is reduced compared to
that of a vertical tire;
whereas a negative
F
ca
F
c
F
c
F
ca camber angle results in
increased cornering
Center of curvature force, in general.
Cambered
wheels result in
a rolling cone
like action.

+g +g
Fca Fc Fc Fca
Center of curvature Effect of
-g -g Camber
Angle
Fc Fca Fca Fc
Center of curvature
 It is clear that the effect of camber angle on cornering

force is different on inner and outer tires.
 However, load transfer during cornering from inner to
outer tires must also be taken into account.

Effect of
Camber
Angle
Carpet Plot
Reduction for
positive camber
angle is
observed.

Camber Stiffness
 The effect of camber angle on cornering
force can be represented by defining a
camber stiffness, i.e. assuming a linear
variation of camber force by camber
angle.
 Camber stiffness is defined, in a similar
manner to cornering stiffness, as the
derivative of the camber force with
respect to camber angle, evaluated at
zero camber angle.

Camber Stiffness & Lateral Force
Fca
Cca = γ=0
γ
Fca
g
Fca = Cca . γ
FL = Cs .α - Cca . γ
Camber angle

Decreased inner Initially equal left Increased outer
wheel load. and right tire loads. wheel load.
-W +W
Effect of load
transfer on
+Fc Cornering
-Fc Force
Due to increased
slope at lower wheel
loads, cornering
force at a given slip
angle decreases!

Effects of Longitudinal Forces
 The cornering (lateral) force developed at any fixed
slip angle is reduced by the application of
longitudinal (tractive or braking) forces. This
basically follows from :
FT2  FC2  mW
where
FT = longitudinal (tractive or braking) force,
FC = cornering (lateral) force,
m = available coefficient of adhesion,
W = vertical tire load.

 Therefore tractive effort reduces the local adhesion
available in lateral direction.
 The same argument applies to the case of braking.
 An approximate method of determining

the lateral or longitudinal tire force,
available during the presence of
tractive and braking effort is possible:
 ‘If it is assumed that the tire may reach
its limiting force condition in any
direction, but that the maximum force
may not exceed a given value in either
the lateral or longitudinal direction’.

Friction Circle
The so called friction circle
concept is based on the
F
c condition that the limiting
mW
force is the same for
cornering and braking
F forces.
T
Braking Traction
FT2  FC2  mW 

2

Friction Ellipse
max
Fc
FcT  In the development of

some tires, however,
either traction or
FTc F max
T braking effort is
Braking Traction
( = m W) emphasized. Thus, the
force limit is different
for traction and
2 2 braking.
 FcT   FTc 
 max  +  max  ≤ 1, FTmax = μW In this case, the
F  
 c   FT  friction ellipse takes
F 
2 the place of friction
FcT = Fcmax 1 -  Tc  circle.
 μW 

Self Aligning Torque

300
250 Fz=2 [kN]

Fz=4 [kN]
200
Self
Fz=6 [kN]
Fz=8 [kN]
150
SAT [Nm]
100
50
Aligning
0
0 2 4 6 8 10 12 14
Torque
-50
-100
Slip Angle [o]
 The variation of self-aligning torque (SAT) may be

presented as a function of slip angle similar to the
presentation of the cornering force characteristics.
 It is, however, more meaningful to plot SAT against
cornering force.

Factors
Affecting
SAT
The application of
tractive effort
increases SAT,
whereas
the application of
a braking force
on the tire has the
reverse effect.

Factors
Affecting SAT
 SAT increases initially
with increasing slip
angle, but after a certain
slip angle it starts
decreasing with
increasing slip angles.
 In the extreme cases of
slip angles around 15o,
SAT may even become
negative.

Factors Affecting SAT
 Another factor affecting SAT is the normal load on
the tire: the greater the load, the greater is the
deflection of the tire, and correspondingly the
length of the contact patch is greater.
From this, it follows that the amount of pneumatic
trail will be increased, and hence the SAT will also
rise.
 Similarly, increasing the inflation pressure for a
given load will reduce the pneumatic trail by the
same means and thus diminish the self-aligning
torque.

 Self aligning torque is one of the factors
related to the straight ahead stability of
the vehicle. Its effect, however, is not as
strong as those of the suspension
parameters such as caster and steering
axis inclination.
 Therefore, SAT is not normally included
in basic vehicle handling models.
 On the other hand, SAT is functional in
giving the driver the so called ‘steering
feel’.


 The cornering behavior of tires are
influenced by various parameters and the
cornering characteristics must be modified
to take their effects into consideration.
 Here, the variation of basic tire cornering

characteristics with varying
 tire inflation pressure,
 camber angle, and
 tractive or braking forces
will be examined.

Effect of Inflation
Pressure
In the usual working range
of tire pressures (around
recommended value),
increasing pressure
improves the cornering
force by increasing lateral
stiffness of the tire.
This improvement continues at higher pressures,

but at a lower rate because of the reduction in the
size of the contact patch.

+g -g
Camber Angle
 The deviation of the

central plane of a tire
from the vertical plane
is known as camber
and the angle between
the two planes is called
Center of curvature the camber angle.
 The angle is positive when the tire leans away from the
vehicle, and negative when it leans towards the vehicle.

Camber Angle
 Camber angle given to the wheels

during assembly is called static
camber.
 Dynamic camber, on the other hand,
may result from the roll motion of the
body during cornering or from wheel
travel.

Camber
Angle
The 1960 Miliken MX1 Camber Car showing a large negative camber

To achieve highest possible side forces all modern racing cars run
with negative camber at all wheels. Formula 1 cars work with a
negative camber of over -4 degree at the front. The wider rear tires
can withstand less camber (appr. -2 degree as a maximum)

+g -g
Effect of
Camber Angle
When the camber angle is
positive, cornering force
is reduced compared to
that of a vertical tire;
whereas a negative
F
ca
F
c
F
c
F
ca camber angle results in
increased cornering
Center of curvature force, in general.
Cambered
wheels result in
a rolling cone
like action.

+g +g
Fca Fc Fc Fca
Center of curvature Effect of
-g -g Camber
Angle
Fc Fca Fca Fc
Center of curvature
 It is clear that the effect of camber angle on cornering

force is different on inner and outer tires.
 However, load transfer during cornering from inner to
outer tires must also be taken into account.

Effect of
Camber
Angle
Carpet Plot
Reduction for
positive camber
angle is
observed.

Camber Stiffness
 The effect of camber angle on cornering
force can be represented by defining a
camber stiffness, i.e. assuming a linear
variation of camber force by camber
angle.
 Camber stiffness is defined, in a similar
manner to cornering stiffness, as the
derivative of the camber force with
respect to camber angle, evaluated at
zero camber angle.

Camber Stiffness & Lateral Force
Fca
Cca = γ=0
γ
DFca
Dg
Fca = Cca . γ
FL = Cs .α - Cca . γ
Camber angle

-DW +DW
Effect of load
transfer on
+DFc Cornering
-DFc Force
Due to increased
slope at lower wheel
loads, cornering
force at a given slip
angle decreases!

 The cornering (lateral) force developed at any fixed
slip angle is reduced by the application of
longitudinal (tractive or braking) forces. This
basically follows from :
FT2 + FC2  mW
where
FT = longitudinal (tractive or braking) force,
FC = cornering (lateral) force,
m = available coefficient of adhesion,
W = vertical tire load.

 Therefore tractive effort reduces the local adhesion
available in lateral direction.
 The same argument applies to the case of braking.
 An approximate method of determining

the lateral or longitudinal tire force,
available during the presence of
tractive and braking effort is possible:
 ‘If it is assumed that the tire may reach
its limiting force condition in any
direction, but that the maximum force
may not exceed a given value in either
the lateral or longitudinal direction’.

Friction Circle
The so called friction circle
concept is based on the
F
c condition that the limiting
mW
force is the same for
cornering and braking
F forces.
T
Braking Traction
FT2 + FC2  mW 

2

Friction Ellipse
max
Fc
FcT  In the development of

some tires, however,
either traction or
FTc F max
T braking effort is
Braking Traction
( = m W) emphasized. Thus, the
force limit is different
for traction and
2 2 braking.
 FcT   FTc 
 max  +  max  ≤ 1, FTmax = μW In this case, the
F  
 c   FT  friction ellipse takes
F 
2 the place of friction
FcT = Fcmax 1 -  Tc  circle.
 μW 


300
250 Fz=2 [kN]

Fz=4 [kN]
200
Self
Fz=6 [kN]
Fz=8 [kN]
150
SAT [Nm]
100
50
Aligning
0
0 2 4 6 8 10 12 14
Torque
-50
-100
Slip Angle [o]
 The variation of self-aligning torque (SAT) may be

presented as a function of slip angle similar to the
presentation of the cornering force characteristics.
 It is, however, more meaningful to plot SAT against
cornering force.

Factors
Affecting
SAT
The application of
tractive effort
increases SAT,
whereas
the application of
a braking force
on the tire has the
reverse effect.

Factors
Affecting SAT
 SAT increases initially
with increasing slip
angle, but after a certain
slip angle it starts
decreasing with
increasing slip angles.
 In the extreme cases of
slip angles around 15o,
SAT may even become
negative.

Factors Affecting SAT
 Another factor affecting SAT is the normal load on
the tire: the greater the load, the greater is the
deflection of the tire, and correspondingly the
length of the contact patch is greater.
From this, it follows that the amount of pneumatic
trail will be increased, and hence the SAT will also
rise.
 Similarly, increasing the inflation pressure for a
given load will reduce the pneumatic trail by the
same means and thus diminish the self-aligning
torque.

 Self aligning torque is one of the factors
related to the straight ahead stability of
the vehicle. Its effect, however, is not as
strong as those of the suspension
parameters such as caster and steering
axis inclination.
 Therefore, SAT is not normally included
in basic vehicle handling models.
 On the other hand, SAT is functional in
giving the driver the so called ‘steering
feel’.

Tire Models

Tire Models
 Mathematical models of tires are particularly
needed
 for the study of complex tire behavior in
an attempt to better understand the
relation between tire construction and
materials and the resulting tire
characteristics,
 for use in complex vehicle dynamic
models in which tires are only one set of
many different components.

Tire Models
 Many different models of various

levels of complexity have been
proposed in the past.
 Efforts to develop new models to
provide more accurate results in
relatively short computer solution
times are still continuing.

Tire Models
 Tire models found in literature may be broadly
classified as follows:
 Experimental (empirical) tire models which are
obtained by fitting a set of curves to experimental
data.
 Analytical tire models, obtained by oversimplifying
the physical behavior of tires so that governing
equations can be solved by relative ease.
 Physical tire models, obtained by including the
actual construction and material properties of tires,
which require advanced numerical techniques for
solution.

Tire Models
 Experimental and analytic tire models are relatively
simple and thus commonly used to simulate tire
behavior in complex vehicle dynamic models.
 Physical tire models, on the other hand, are mainly
developed for the study of the behavior of tires in detail.
 Thus, these models include the actual
construction and material characteristics of tires.
 Due to large computer resources needed and
longer computational time required for solution,
when compared to the other two models, they are
not suitable for use in vehicle dynamic models.

Tire Kinematics
 For a tire of radius rw rotating at a speed of w
vithout slip, the velocity of the center of tire
contact patch in the longitudinal direction is :
u' = ωrw
 Thus, for a tire rotating with tractive or braking
longitudinal slip velocities of Vsx, the velocity of
the center of tire contact patch in the longitudinal
direction becomes :
u = u' ± Vsx
Note: “-” for driving and “+” for braking.

x

 Vs Freely Rolling Tire
u 
V 0
(No slip)
sx
  
V Vs  Vsy j
Vs  u tan 
 
y u  rw w i
rw
w
One full rotation
L=2prw
Braking Torque
Velocity of the center


of tire contact patch

u  rw w + Vsx
rw
w
One full rotation
L=2prw
Translational
Slip
Distance travelled is the sum of distance
covered by rotation plus translational slip.

Driving Torque
Velocity of the center of


tire contact patch :
u  rw w - Vsx
rw
w
One full rotation
Rotational Slip
L=2prw
A portion of the displacement for one full

rotation is lost because of tire deformation!

Braking or Driving Torque
x Vsx x

V
Vsx  Vsy Vs
sx
Vs rw w
 Vsy
 u
rw w u
 
V V
y y
    
Velocity of the u   rw w  Vsx  i V  u + Vsy j
center of tire     
Vs  Vsx i + Vsy j   rw w  Vsx  i + Vsy j
contact patch

Definitions
V u - rw w r w
 Longitudinal Slip S  sx   1 - w (braking : u  rw w)
u u u
V r w-u u
S  sx  w  1- (traction : u  rw w)
rw w rw w rw w
Vsy
 Lateral Slip Sy   tan 
u
Case1 : u  rw w then S  0
w
rw (freely rolling tire!)
Case 2 :  w  0 and u  0  then S  1
(tire sliding without rotation!)
Case 3 : u  0 and w  0  then S  1
(tire slipping while rotating!)
Slip Velocity
 Traction
Vs2  Vsx
2 2
+ Vsy
Slip velocity can be
 

2 2
written in terms of Vs   uS  + uS y
the longitudinal and
lateral slips.  u S 2 + S 2y
2
 u  2
Vs  u  1 -  + tan 
 rw w 

Slip Velocity
 Braking
Vs2  Vsx
2 2
+ Vsy
Slip velocity can be
 
2 2
written in terms of Vs   uS  + uS y
the longitudinal and
lateral slips.  u S 2 + S 2y
2
 rw w  2
Vs  u  1 -  + tan 
 u 

Friction Coefficient
The dependence of tire-

road friction coefficient, m,
to slip velocity, Vs is μ = μ o 1 - A s Vs 
usually approximated by a
linear equation.
where
mo : peak coefficient of friction
As : friction reduction factor.
Vs : slip velocity.

Linear Tire Model
We have already seen that a linear

representation of tire cornering force
can be used for slip angles of below
about 4o.
Fc = Cs 
Here, the cornering stiffness, Cs, is the
slope of the cornering force versus slip
angle characteristic.

Linear Tire Model
 Note that the linear tire model is valid

in a region of tire slip angles between
the origin and up to about 4o.
 If the slip angle of the tire during any
maneuver exceeds this value, then the
nonlinear cornering force
characteristics of tire behavior should
be taken into consideration.

Exercise
Discuss the effects on simulation

results, in case you use a linear
tire model but the slip angles go
into the nonlinear range, i.e.,
above say 4o, for the maneuvers
simulated.

Nonlinear Tire Models
 In the following section a couple of
nonlinear tire models will be introduced.
 Dugoff model – for simplified
simulations when qualitative
representation of tire lateral and
longitudinal forces is sufficient.
 Magic Formula – used for simulations
that require accurate tire cornering and
longitudinal force characteristics.

Dugoff Model
 Dugoff model1 is based on the simplified
analytical tire mechanics analysis, which is
based on the classic theoretical study of the
freely rolling tire due to Fiala2.
 It is mathematically simple and provides a
qualitatively valid nonlinear tire model, which
describes the complex tire behavior in terms
of a limited number of basic parameters.
1 Dugoff, H., P. S. Fancher, and L., “Tire performance characteristics affecting

vehicle response to steering and braking control inputs” Final Technical
Report, Highway Safety Research Institute, Ann Arbor, Michigan, 1969.
2 E. Fiala, "Seitenkraefte am Rollenden Luftreifen", VDI Zeitschrift, vol.96, 1954.

Dugoff Model
 Although there are quantitative
discrepancies between the curves
generated by the equations derived by
Dugoff and the ones derived from the
experimental data, there is considerable
qualitative agreement between theory and
measurement.
 Further, Dugoff model provides
longitudinal and lateral tire forces in
combined maneuvers (traction and
cornering) directly.

Formulation
Vs = u S 2 + tan 2 α
 slip angle [rad] m  mo 1 - As Vs 
S longitudinal slip μ Fz 1 - S 
Vs slip velocity [m/s] z=
 Cl S  +  Cs tan α 
2 2
u velocity component in wheel plane [m/s] 2
mo static (peak) tire/road friction coefficient
As friction reduction factor [s/m] z 2 - z  for z  1
f (z )  
Cs cornering stiffness [N/rad]  1 for z  1
Cl longitudinal stiffness [N/unit slip]
Fx driving/braking force [N] Cl S f(z)
Fx =
Fy cornering force [N] 1-S
Fz normal load on tire [N]
Cs tan α f(z)
Fy =
1- S

Typical Data
Normal tire load, Fz : 2 [kN]
Cornering stiffness, Cs : 22.9 [kN/rad]
Longitudinal Stiffness, Cl : 80 [kN/unit slip],
Static tire/road friction coefficient, mo : 0.83,
Friction reduction factor, As : 0.0115 [s/m],
Velocity component in wheel plane, u : 15 [m/s]
Note that:
- a single cornering stiffness is given
in the basic data. Cornering stiffness
should be specified for each tire load. Fx
- the longitudinal stiffness is the initial
slope of the longitudinal force Fx
versus longitudinal slip S
characteristic. 00 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
Longitudinal slip, S

1400
= 2 [o ]
1200 = 4 [o ]
F z = 2 kN = 6 [o ]
1000
800
600
400
200
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
Longitudinal slip, S
Lateral tire force versus longitudinal slip – Dugoff Model

Fx [N]
Longitudinal tire force versus longitudinal slip – Dugoff Model

Fy [N]
Lateral force versus longitudinal force – Dugoff Model

Fy [N]
Lateral force versus slip angle – Dugoff Model

A single
cornering
stiffness is
F y [N]
used for all

tire loads.
This is not
realistic !
Here different
cornering
stiffness is
F y [N]
used for
each tire
load.
Realistic !

Dugoff Model
The accuracy of the Dugoff model is poor at higher
longitudinal slips.
Exp. Exp.

Dugoff Model
Better accuracy can be obtained if the
road/tire friction coefficient is defined
as follows:
μ x = Kμ max
0.3 - S x
K = 1+ 2
S x + 3.237μ max - 1.456μmax + 0.7

MAGIC FORMULA

Magic Formula
 Of the many experimental (empirical) tire
models, a trigonometric model first
suggested by Bakker, Nyberg, and Pacejka1
has gained wide usage. This model is now
called the Magic Formula because of its
power to provide many important tire
functions accurately.
1 Bakker, E., Nyborg, L., and Pacejka, H. B., “Tyre Modelling

for Use in Vehicle Dynamics Studies”, SAE Technical Paper
870432 (1987).

Magic
Formula
Examples of
curves that can be
generated by the
Magic Formula

Magic Formula
 In its full form, the Magic Formula is given as :

Y = S y +Dsin C tan-1  B  X - Sx 1- E  +E tan-1 B  X - Sx  
  
where
 Y : any tire response quantity (lateral force,
longitudinal force, self aligning torque, etc.).
 X : lateral slip angle, longitudinal slip, or
camber angle.
 The six constants Sx, Sy, B, C, D, and E must
be obtained from experimental data by some
curve fitting procedure.

Magic Formula

Y = Sy +Dsin C tan-1  B  X - Sx  1- E  +E tan-1 B  X - Sx  
  
 Here the constants Sx and Sy represent the offset
of the variable of interest (Fy, Fx, Mz, etc.) at zero
independent variable (lateral or longitudinal slip,
or camber angle, etc.) from the origin. These are
usually assumed to be zero for simplicity of
analysis.
 Note: the constant D gives the peak value, and the
product BCD gives the initial slope of the curve.
 Warning: keep in mind that Magic Formula is
strictly dimensional – Note the units used!

P 205/60 VR 15 6J  
Y  S y + D sin C tan - 1 B X - S x 1 - E  + E tan - 1 B X - S x  
a0 1.6929
Y = Fy X=α Lateral Force Fy [N]
a1 -55.2084
Slip Angle  [o]
a2 1271.28 C = a0 Normal Load Fz [kN]
a3 1601.8
D =  a1Fz + a2  Fz Camber Angle g [o]
a4 6.4946
a5 4.80E-03
 -1  Fz 
a6 -0.3875 B  a 3 sin  2 tan     1 - a5 g  / CD
a7 1   a4  
a8 -4.54E-02
D = max(Fy)
a9 4.28E-03 BCD = tire cornering stiffness [N/o]
a10 8.65E-02 a2 / 1000 = mym max. lateral friction coefficient
a111 -7.9730
a112 -0.2231 E = a6 Fz + a7 S x = a8 γ + a9 Fz + a10
a12 7.6680
Sy = a11γFz + a12Fz + a13 a11 = a111Fz + a112
a13 45.8764

Camber angle = 0 [o ]
7000
Fz  1 to 8 [kN]
6000
5000
4000
3000
2000
1000
-1000 P 205/60 VR 15 6J
0 5 10 15
Slip angle, [ o]
Lateral force versus slip angle – Magic Formula


Y = S y + Dsin Ctan-1 B  X - S x 1 - E  + Etan-1 B  X - Sx  
  
P 205/60 VR 15 6J
Fx : Longitudinal Force [N]
b0 1.6500 Y = Fx
S : Longitudinal Slip [%]
b1 -7.6118
b2 1122.6
X=S Fz : Normal Tire Load [kN]
mxm : Maximum Longitudinal
b3 -7.3600E-03 C = b0 Friction Coefficient
b4 144.8200
b5 -7.6614E-03 D =  b1Fz + b2  Fz D = max(Fx)
b6 -3.8600E-03
b7
b8
8.5055E-02
7.5719E-02
 
B  b 3 Fz2 + b 4 Fz e- b5 Fz / CD
BCD = tire
longitudinal
stiffness
b9 2.3655E-02
E = b 6 Fz2 + b7 Fz + b8 [N / % slip]
b10 2.3655E-02
b2 / 1000 = mxm
S x = b9 Fz + b10 S y = 0 max. longitudinal
friction coefficient
Fxmax
Fx
Fxmin
Sx

Fz  1 to 8 [kN]
Longitudinal force, F [N]
x
P 205/60 VR 15 6J
Longitudinal force versus longitudinal slip – Magic Formula

P 205/60 VR 15 6J
c0 2.2264

Y = S y + Dsin Ctan-1  B  X - S x 1 - E  + Etan-1 B  X - S x  
  
c1 -3.0428
c2 -9.2284
Mz : Self Aligning Torque [Nm]
 : Slip Angle [o]
c3 0.5000088 Y = Mz X = 
c4 5.56696 g : Camber Angle [o]
c5 -0.25964 Fz : Normal Load [kN]
c6 -1.29724E-03
c7 -0.358348 C = c0 D = c1Fz2 + c2 Fz
c8 3.74476
c9
c10
-15.1566
2.1156E-03

B  c 3 Fz2 + c4 Fz   6 
1 - c g e - c5 Fz
/ CD
E =  c7 Fz2 + c8 Fz + c9   1 - c10 g 
c11 3.4600E-04
c12 9.13952E-03
c13 -0.244556
c14 0.100695 S x = c11 g + c12 Fz + c13
c15 -1.3980
c16
c17
0.44441
0.998344
 
S y = c14 Fz2 + c15 Fz g + c16 Fz + c17
300
250 Fz=2 [kN]
Fz=4 [kN]
200
Fz=6 [kN]
Fz=8 [kN]
150
SAT [Nm]
100
50
0
0 2 4 6 8 10 12 14
Slip Angle [o]
-50
-100
Self Aligning Torque versus slip angle – Magic Formula

Magic Formula
 A major difficulty in using the Magic Formula in
vehicle dynamics arises from the fact that the
longitudinal and lateral tire forces represented
by the relevant equations are applicable for
individual traction/braking or cornering
maneuvers only.
 To obtain longitudinal and lateral tire forces for
combined (simultaneous traction/braking and
cornering) maneuvers, a number of
approaches are suggested in the literature but
as yet there is no universally accepted single
procedure.

Burckhardt Model
 A simpler alternative, Burckhardt Model, to the Magic
Formula (MF) to represent the tire longitudinal force
characteristics may also be useful in some
applications.
μx  s  =
Fx
Fz
 
= c1 1 - e -c 2s - c 3s
 mx : the tire/road longitudinal friction coefficient,

 Fx : the longitudinal force,
 Fz : the normal load on the tire,
 ci : constants to be determined, and
 s : longitudinal slip.

ME 513
VEHICLE
DYNAMICS

CH II - VEHICLE HANDLING
Part 2

Bounce
Yaw
Vehicle
Side Slip
Motions
Pitch
ISO axes
Fore &
Aft
Roll
If considered as a rigid body, the vehicle has six

degrees of freedom, three translations and three
rotations in space.

PLANE MOTIONS
of
VEHICLES

PLANE MOTIONS OF VEHICLES
 The full range of lateral acceleration, aL, to which a
vehicle is subjected under various conditions, can
be divided into two phases.
 In the first phase, the vehicle is subjected up to 0.3g
(~3 m/s2) lateral acceleration* during:
 Soft cornering and normal straight ahead driving;
the side loading being caused by
 cross winds,
 incorrectly banked road, and
 variations of road surface.
* International standart ISO-TR 8726, Road vehicles - Transient open-loop response
test method with pseudo-random steering input, 1988.

Low Lateral Acceleration Range
 In this range of lateral accelerations, vehicle roll is
usually not significant and hence the flexibility of
suspension can be neglected.
 A simplified vehicle model taking into account of plane
motion of the vehicle only is then sufficient to describe
basic handling behavior in terms of main vehicle
parameters such as
 weight distribution, Suspension stiffness
 tire specification, is not effective !
 tire pressures, and Vehicle is practically a
 tire camber angles. rigid body.

Low Lateral Acceleration Range
 With a given vehicle design, changes to

ground stability during the first phase
can only be made by alterations to these
main parameters.
 In the earlier design stages of a new
model when suspension data are not yet
available, such information will be
extremely valuable.

High Lateral Acceleration Range
 The second phase, from 0.3 g

upwards, involves average to hard
cornering and maneuvers usually
producing visible roll.
 In this phase, changes to spring
Suspension
rates, roll axis, and the addition of
stiffness is
auxiliary roll stiffening devices
now effective!
affecting the distribution of roll
stiffness are effective.

High Lateral Acceleration Range
 The handling of the vehicle must, therefore,

be checked using a three degree of freedom
(yaw, side slip, and roll) model, as soon as
the suspension data are available.
 Such an analysis normally requires a
mathematical model of the vehicle
including the suspension system and the
solution can be performed on a computer.

Vehicle Model for Plane Motions
 Note that on a plane road surface,

the motion of a vehicle on a curve
can be represented by the
superposition of three motions.
 translation in the longitudinal
(forward) direction,
 translation in the lateral
direction, and
 rotation about the center of
gravity.

O x
Z X Vehicle
Model
r
l C d  A simple four-
y wheel vehicle
v z u moving on a
x
smooth level
Y h
road surface.
 OXYZ is the inertial reference frame fixed to the road

(SAE axes here).
 The moving reference frame Chzx is attached to the
center of gravity of the vehicle.

O
Z X Vehicle
Model
r
l C
The variables
z d describing the
v
u
x motion of the
vehicle may be
Y y restricted to the
three variables :
 u (longitudinal velocity),
 v (lateral velocity), and
 r (yaw velocity)
of the center of gravity of the vehicle.

Vehicle Model
 The use of a body fixed reference frame located at
the centre of gravity with the axes corresponding to
the principal axes of inertia, allows us to define
 fixed rotational inertia properties,
 zero products of inertia,
together with
 no inertial coupling between translational and
rotational coordinates
resulting in the simplest set of governing equations.

O x
Z X Vehicle
Model
l C d
It is obvious
y that the
trajectory of
the vehicle
Y can be plotted
if the
 location of the center of gravity (x and y) and

 direction of the vehicle longitudinal axis l (heading
angle)
are known at any instant.

O x
Z X Vehicle
Model
r
Location and
l C d heading of the
y
z u vehicle w.r.t.
v x the inertial
Y h reference
frame
λ = r y = vcosλ + usinλ x = ucosλ - vsinλ

Note that it is not possible to obtain the x and y position
of the center of gravity simply by integrating u and v !
You need to integrate x and y. 
O x
Z X
r Vehicle
l C d
y Model
v z u
x
Y h
 x  cosλ -sinλ 0   u  Velocity and

heading of the
     

 y  =  sinλ cosλ 0   v  vehicle w.r.t. the
l   0 0 1  r  inertial reference
     frame

Assumptions
1) The external forces, cornering and traction or

braking, are assumed to act from ground to
the four tires in their contact centers,
disregarding the self-aligning torques.
2) The effects of tractive (or braking) force on
tire side force characteristics are neglected.
3) The wheels of the vehicle are supposed to
remain in a vertical position, i.e., zero
camber.

Assumptions
4) Suspension is assumed to be sufficiently
rigid so that out of plane motions (roll and
pitch) can be disregarded.
5) Load transfer from inside wheels to
outside wheels during cornering and the
corresponding change in tire side force
characteristics are negligible.
6) Steering control with equal steering angles
d of the left and right front wheels is
assumed.

C 
x i
r VC Acceleration
h      
O z j VC  ui  vj ω  rk
 X  
k   dV   dV   
a XYZ       ωxV
 dt  XYZ  dt C
  du  dv     
Y a XYZ   i   
j   rk x u i  v j 
Z  dt dt 
 
i, j, and k The acceleration of the vehicle
are the unit vectors cog in the inertial reference frame:
  
a  u  vr i  v  ur  j
in the body fixed
reference frame.

Longitudinal and
 Note that the longitudinal and
lateral accelerations of the
vehicle body cog are not u and v.

ax  u a y  v
 They are given by:
a x = u - vr a y = v + ur

L
b a
Fx22
Fx12 Two
Fy22
r
u
Fy12 Track
t
d Vehicle
Fx21
v Model
Fx11
Fy11
Fy21
m(u  vr )  (Fx11  Fx12 ) cos d  (Fy11  Fy12 ) sin d  (Fx 21  Fx 22 )

m( v  ur )  (Fx11  Fx12 ) sin d  (Fy11  Fy12 ) cos d  (Fy 21  Fy 22 )
 
Jr  (Fx11  Fx12 ) sin d  ( Fy11  Fy12 ) cos d a  (Fy 21  Fy 22 ) b 
 t

 (Fx11  Fx12 ) cos d  ( Fy11  Fy12 ) sin d  (Fx 21  Fx 22 )
2
t
2

Bicycle Model

Bicycle Model
 Even with the simplifying assumptions

already made in the modeling process,
the resulting set of equations
describing the motion of the full (four
wheel) vehicle model are still too
complex for a general analysis.
 Considerably simpler linear equations
can be obtained, however, by further
observations and assumptions.

Slip Angles - No steer angle is given (for clarity!)
L
b a
22 u+rtr/2 u+rtf/2
12
v-br v+ar
a22 r
dom a12
tr u  dom
i tf
v
u-rtr/2 u-rtf/2
21 v-br 11 v+ar

a21 dom j a11 dom
     t f     v  ar
   
V11 = ui + vj + rk x  ai +
t
j =  u - r f  i +  v + ar  j
tanα11 
 2   2  tf
ur
2

v  br v  ar
tan a 22  tan a12 
ur r
t Slip Angles ur f
t
2 2
L
b a
u+rtr/2 u+rtf/2
v-br v+ar
a22 r
dom a12
u dom
tr tf
u-rtr/2 v u-rtf/2
v-br v+ar
a21 dom a11 dom
v  br No steer angle is v  ar
tan a 21  tan a11 
t t
ur r given (for clarity!) ur f
2 2

Slip Angles
v  ar
tan a11 
tf
ur For small slip angles :
2 tan a  a v  ar
v  ar and af 
tan a12 
t rtf / 2  u u
ur f
2
v  br
tan a 21 
tr For small slip angles :
ur
2 tan a  a v  br
v  br and ar 
tan a 22 
t rtr / 2  u u
ur r
2

Bicycle Model
 The representation of a slip angle for
the front and rear axles, af and ar,
allows the combination of tire side
forces for each axle, i.e.,
 Fyf = 2Fy11 = 2Fy12 front axle

side force
 Fyr = 2Fy21 = 2Fy22 rear axle
side force

L
b a
Fxr u u u
v-br Fxf v+ar d
ar r
v - af dom
dom F
Fyr yf
Bicycle Model
The vehicle model
considered may now Fyf  Cright
s a f  Cs
left
a f  Cfaxle
s af
be reduced to a two-
wheel (single track) Fyr  Cright
s a r  Cs
left
a r  Craxle
s ar
model.

L
b a
u u u
v-br v+ar d
ar r
v - af dom
dom
v  ar v  ar
 af  d   af  d Slip
u u
v  br Angles
ar 
u

L
b a
Fxf
Fxr
u
r d
Fyf
Fyr v
Equations of Motion
mu  vr   Fxf cos d  Fxr  Fyf sin d
mv  ur   Fyf cos d  Fyr  Fxf sin d
Jr  aFyf cos d  bFyr  aFxf sin d
mu  vr   Fxf cos d  Fxr  Fyf sin d
Linearized
mv  ur   Fyf cos d  Fyr  Fxf sin d
Equations
Jr  aFyf cos d  bFyr  aFxf sin d of Motion
 If the steer angle, d, is
assumed to be small
such that
 cos d = 1 and
 sin d = 0 m  u - vr  = Fxf + Fxr
 The three equations of
motion for the vehicle
m  v + ur  = Fyf + Fyr
model reduce to :
Jr = aFyf - bFyr

Reduction of Equations of Motion
 The term vr on the left hand side is a product of two
variables of small magnitude and hence can be
neglected.
mu  vr   Fxf  Fxr
 Then this equation is uncoupled from the other two
equations and can be used to study acceleration
performance of the vehicle in straight motion.
mv  ur   Fyf  Fyr

 The remaining two equations
involve the yawing and side
slip motions of the vehicle
Jr  aFyf  bFyr
and will be examined in detail.

Reduce number of variables
mv  Ur   Fyf  Fyr  The forward speed of

the vehicle will be
Jr  aFyf  bFyr kept as a constant
parameter, U, and
thus the degree of
freedom of the
Variables  v and r vehicle model is now
reduced to two.

Express Cornering Forces in terms of
Slip Angles
 Under normal driving mv  Ur   Fyf  Fyr
conditions, slip angles < 4o
Jr  aFyf  bFyr
 Therefore, assumption of small Fyf  Cf a f

slip angles is justified.
Fyr  Cr a r
 Consequently, tire cornering
force characteristics in the mv  Ur   Cf a f  Cr a r
linear range can be expressed
in terms of tire cornering Jr  aCf a f  bCr a r
stiffnesses.

Insert Slip Angle Expressions
mv  Ur   Cf a f  Cr a r
Jr  aCf a f  bCr a r
v  ar
af  d
u
v  br
ar 
u
. v r
m(v + Ur) = (Cf + Cr ) + (aCf - bCr ) - Cf δ
Equations of U U
Motion v r
Jr = (aCf - bCr ) + (a Cf + b Cr ) - aCf δ
2 2
U U
. v r
m(v + Ur) = (Cf + Cr )
U
+ (aCf - bCr ) - Cf δ
U Equations
v r
Jr = (aCf - bCr ) + (a2Cf + b 2Cr ) - aCf δ
U U of Motion
 M x  =  N x +  P u in matrix form
 (Cf + Cr ) (aCf - bCr )  mU 2 

m 0   v    v  - C 
 0 J  r   
U U      f δ
     (aCf - bCr ) (a 2Cf  b 2Cr )  r  - aCf 
 
 U U 
x   M  N x  M
1
 
 1
P u
x  =  A x +  B u

State Equation Form
  =  A  x(t) +  B  u(t)   E  z ( t )
x(t)
y( t ) =  C x(t) +  D u(t)
 {x(t)} : Vector of state variables
 {u(t)} : Vector of control input
 {z(t)} : Vector of disturbance input
 [A] : System matrix
 [B] : Control input distribution matrix
 [E] : Disturbance input matrix
Set of 1st order differential equations.
 {Y} : Output vector
 [C], [D] : Output selection matrices

State Equation Form
The general matrix equation of motion for multi-degree
of freedom lumped parameter systems is given by
 Mq +  Cq  +  K q = f 

where [M], [C], and [K] are the mass, damping and
stiffness matrices and {q} and {f} are the vectors of
mass displacements and applied external forces,
respectively.
This form can be converted into the state equation form
x  =  A x +  B u

 Mq +  Cq  +  K q = f  Conversion to
Introducing a state vector State Equation
 x1  q 
x =   =    Form
 x2  q 
 x 1   q   x 2 
and using the fact that x  = x  =  =   
 2  q   q 
the state space form
x  =  A x +  B u
can be obtained in terms of the mass, damping, and
stiffness matrices
  0  I    0 
x  =   x +  -1  u
 -  M   K  -  M   C  M  
-1 -1

State Equation Form
 Cf  Cr aCf  bCr   Cf 
 U -
 v   mU mU  
v  m 
 = 2 2    +  aC d
 r   aCf  bCr a Cf  b Cr   r  - f 
 JU JU   J 
State Variables : v and r

Control Input : d
Forward Speed : U
Vehicle Parameters : m, J, a, b, Cf, Cr

Stability Analysis

  =  A  x(t) +  B u(t)
x(t) Stability Analysis
 To study the stability, the roots of the characteristic
equation is examined.
lI  A  0
a1ln  a 2ln 1  ...  an l  an 1  0

 The roots of the characteristic equation are also the
eigenvalues of the system matrix A.
 Eigenvalues of a square matrix A can be obtained using
the Matlab command eig(A).

Roots of the Characteristic Equation on Complex Plane
-λ R
λ = λ R + jλ I = -ζωn + jωn 1- ζ 2 2
ωn = λ R + λ I2 , ζ =
ωn

 C  Cr aCf  bCr   Cf 
.  f  U
v
. = 
mU
2
mU
2
v  - m 
   +  aC d
Stability
 r   aCf  bCr a Cf  b Cr   r  - f 
Analysis
 JU JU   J 
lI  A  0
Cf  Cr aCf  bCr
l  U
mU mU
2 2
0
aCf  bCr a Cf  b Cr
 l
JU JU
 This determinant is going to result in

the quadratic characteristic equation.

C  Cr aCf  bCr
l f
mU

mU
U
0
Stability
2 2

aCf  bCr
JU
l
a Cf  b Cr
JU
Analysis
This gives a second order polynomial characteristic
equation.
2
a1l  a 2l  a 3  0 a1  1
C C a 2
C  b 2
C 
a 2   f r  f r 
 mU JU 
 
 Cf  Cr   a Cf  b Cr   aCf  bCr
2 2
  aCf  bCr 
a3      U 
 mU   
JU   mU  JU 

a1l2  a 2l  a 3  0
a1  1 Stability
C C a 2
C  b 2
C  Analysis
a 2   f r  f r 
 mU JU 
 
 For a system with a quadratic characteristic

equation, the necessary and sufficient condition is
that all three coefficients be positive (or of the
same sign).
 It is clear that a1 and a2 are always positive -
remember that both Cf and Cr are both negative
quantities!
 Thus, the stability of the system will depend on a3.

a1l2  a 2l  a 3  0 Stability Analysis
 Cf  Cr   a Cf  b Cr   aCf  bCr
2 2
  aCf  bCr 
a3       U  
 mU   
JU   mU  JU 
 After some manipulations, the

expression for a3 can be reduced to :
 L2C C  aC  bC mU 2 

a3   f r f r
2

 mJU 

 L2C C  aC  bC mU 2 
a3   f r f r 
2
 mJU  Stability
Analysis
 Therefore for a3 > 0, the condition to

be satisfied is :
L2Cf Cr  aCf  bCr mU 2  0

Stability Analysis
L Cf Cr  aCf  bCr mU  0
2 2
Divide both sides by the

positive terms
L  a b  m U2
LR and CfCr     0
R  Cr Cf  L R
where R is the radius of
curvature!
1 U2 L  a b W
aL      aL  0
g R R  Cr Cf  L

Stability Analysis
L  a b W
    aL  0
R  Cr Cf  L
Define Understeer Coefficient:

Condition for
 a b W Stability !
K us    
 Cr Cf  L L
 K usaL  0
 Wr Wf  R
K us    
 Cr Cf 

L
 K usaL  0
R Stability Analysis
 a b W
K us    
 Cr Cf  L
Three distinct cases can be identified:
aCf = bCr  K us = 0  Neutral Steering

aCf < bCr  K us > 0  Understeering
aCf > bCr  K us < 0  Oversteering

L
 K usaL  0 Stability Analysis
R
Unconditionally stable Cases :
1. Neutral
L
K us = 0  aCf = bCr  >0 Steer
R Vehicle
2. Understeer L
Vehicle K us  0  aCf  bCr   K usaL  0
R

L U2 Stability
 K usaL  0 L  K us  0
R g
Analysis
3. Oversteer Vehicle
L 2 gL
K us  0  aCf  bCr   K usaL  0 if U 
R K us
L 2 gL
unstable   K usaL  0 if U 
R K us
System becomes unstable  gL

at a forward speed:
Ucr 
K us

Handling
Characteristics

Steady State
 If a steady angle of steer, d, is applied  v  0
and held, a steady state where   
 yawing and slip accelerations
 r  0
become zero
U
is reached. r=
R
1  Cf  Cr (aCf  bCr )  mU2   v 1

     Cf  d
U aCf  bCr a2Cf  b2Cr  r  a

Yaw Velocity and Lateral
Acceleration Gains
1  Cf  Cr (aCf  bCr )  mU2   v 1
 2 2     Cf  d
U aCf  bCr a Cf  b Cr  r  a
 Eliminate v to  
get Yaw r ULCf Cr
ss
= 
Velocity Gain d 2 2
 (aCf  bCr ) mU  L Cf Cr 
 Multiply r by U
to get Lateral U2 / R  U 2LCf Cr 
Acceleration ss
= 
d 2 2
Gain  (aCf  bCr ) mU  L Cf Cr 
 Exercise : Obtain the expression for the side slip velocity gain.
Handling
Characteristics
 The neutral steer vehicle has
yaw velocity gain in direct
proportion to its speed.
 The oversteer vehicle has
yaw velocity gain above taht
of neutral steering vehicle
and reaches an infinite value,
indicating instability, at the
critical speed.
 The understeer vehicle
reaches a maximum gain at
the characteristic speed, and
this value is one-half that of
the neutral steer vehicle.

Handling
Characteristics
 A similar plot is obtained for
the lateral acceleration gain.
 The understeer vehicle
tends toward a constant
value at speeds well above
the characteristic speed.
 In other words, the vehicle
would respond to a
particular input of steer
angle with nearly the same
steady state lateral
u acceleration regardless of
vehicle speed at very high
speeds.

K us  0  aCf  bCr Neutral Steering
  neutral
r ULCf Cr r steer U
= 2 2  
δ ss  (aCf - bCr ) mU + L Cf Cr  d L
ss
For steady state motion : U

r
R is the radius of the curve R
2   neutral
U /R U 2LCf Cr 2
U / R steer U2
ss
=  
d 2 2 d L
 (aCf  bCr ) mU  L Cf Cr  ss

K us  0  aCf  bCr Oversteer
 A critical speed is
Set denominator of the yaw defined for an oversteer
velocity gain equal to vehicle.
zero.  It corresponds to the
forward speed at which
(aCf  bCr ) mU 2  L2Cf Cr  0 the yaw velocity and
lateral acceleration gains
both go to infinity.
2
2  L Cf Cr Thus the vehicle is
Ucr  
aCf  bCr m stable for U<Ucr, and is

unstable for U>Ucr.

K us  0  aCf  bCr Understeer
 A characteristic speed
can be defined for an
UchLCf Cr Uch understeer vehicle as
 that forward speed at
2 2 2L
(aCf  bCr )mUch  L Cf Cr which the yaw velocity
and lateral acceleration
d  UchLCf Cr  gains are one-half the
 2 2 =0
dUch  (aCf - bCr )mUch + L Cf Cr  corresponding gains of
a neutral steer vehicle.
 Alternatively,
characteristic speed
2 corresponds to max.
2 L Cf Cr yaw velocity response
Uch 
aCf  bCr m of an understeer vehicle.

Summary

b U

Exercises

Exercises
1. Show that the stability condition
(aCf  bCr ) mU 2  L2Cf Cr  0
can be obtained by setting the determinant of

the coefficient matrix of the state variable
vector in the steady state equations equal to
zero.
1  Cf  Cr (aCf  bCr )  mU2   v 1

 2 2     Cf  d
U aCf  bCr a Cf  b Cr  r  a

Exercises
2. a) Show analytically that the peak of the yaw
velocity gain for an understeer vehicle is
obtained at characteristic speed and
obtain the expression for the peak value.
b) Derive the limiting value of the lateral
acceleration gain of
- understeering , and
- neutral steering
vehicles at very high forward speeds
using vehicle and tire specifications.

Exercises
3. Relevant parameters are specially tuned to give a
car neutral steering handling behavior. If the front
wheels are given a steering angle of 1.8o, while
travelling at a constant speed of 90 kph:
a) What will be the resulting steady state yaw
velocity of the car ?
b) What will be the radius of the turning circle ?
Car data : Mass: 1320 kg, L = 2.47 m,

Cf = 72000 N/rad, Cr = 80000 N/rad
Ans.: a) 18.2 [o/s], b) 78.6 [m]

Exercises
4. A conventional car is following a steady state
circular trajectory of 120 m radius. What is the
maximum speed below which the use of a
model for the plane motions of the vehicle (i.e.,
2 degree of freedom bicycle model) is adequate
so that reasonably representative information
on the motion of the vehicle may be obtained
without the need of suspension data? Briefly
explain your reasoning.
Ans: 67.7 [kph]

Exercises
5. The yaw velocity gain expression (with respect to the
steering angle d in [rad] of the front wheels) of a car is
given as a function of forward speed U in [m/s], below.
a) Determine the value of the critical or characteristic

speed (specify !), whichever is appropriate, of the car.
b) Calculate the maximum yaw velocity gain of the vehicle.
c) Calculate the lateral acceleration of the car at a forward
speed of 90 kph for a steering angle of 2o of the front
wheels.
Ans.: a) 107 [kph], b) 6.0 [s-1], c) 5.15 [m/s2]

Exercises
6. Consider a 22 ton truck
with three axles as
illustrated in the figure (all
dimensions are in [mm]).
It has a single set of tires on the front axle. Tandem wheels are used
on the intermediate and rear axle. The cornering stiffnesses of tires
are given below.
Front axle : Csf = 2.00 [kN/o] per tire
Intermediate axle : Csi = 1.25 [kN/o] per tire
Rear axle : Csr = 1.50 [kN/o] per tire
- Use a suitable model to estimate the handling behavior of the truck
and calculate the corresponding critical or characteristic speed.
- Ans.: uch = 188 [kph]
- Can you comment if this characteristic speed is suitable or not?
Hint: For a truck the maximum speed is likely to be around 100 kph.
Another look at
Handling

Wr a L L
b a Wf a L
dom
r
Front and
ar
u
Rear Slip
dom
v Fyf af
d
Angles
Fyr
R
L
L
d  af
d   af  ar
R R
ar
d  af

W r aL L
b a Wf a L
dom
ar
r Front and
u
dom
d
Rear Slip
af
Fyr
v Fyf
Angles
R
d  af
L
R Wr aL  Fyr cos a r
Wf aL  Fyf cos d  a f 
ar
d  af
 For small steering and slip angles :
Wr aL  Fyr Wf aL  Fyf

Thus, the front and rear slip Front and Rear
angles can be related to the
lateral acceleration. Slip Angles
Fyf
Wf b W b mg
αf = = aL = aL = aL
Cf Cf L Cf L Cf L
d   af  ar
Fyr
Wr a W a mg
R
αr = = aL = aL = aL
Cr Cr L Cr L Cr
L  a b  mg L
d      aL d   K usaL
R  Cr Cf  L R
Note here that the cornering stiffnesses Cf and Cr must
have negative values.
L L
d   a f  a r d   K usaL Front and
R R Rear Slip
Angles
a f  a r  K usaL
 A comparison of these equations reveals that the

handling behavior of a vehicle depends on the
difference between the front and rear slip angles.
 The handling behavior of a vehicle can thus be
interpreted in terms of the relative magnitudes of the
front and rear axle slip angles.

af
VF
Front and Rear
Slip Angles
Vehicle
Side path
force af > a r
 For an understeer
Understeer
Centrifugal
vehicle the resulting
a VR
r force front slip angle is larger
than the real slip angle
and thus the generated
centrifugal force tends
to cancel the effect of
I the side force.

V
a F
f
I
Front and Rear
Centrifugal
Slip Angles
force
Vehicle  For an oversteer
Side force
path vehicle the centrifugal
ar force tends to add to
a < a
f r
VR the side force as now
Oversteer the rear slip angle is
larger than the front
slip angle.

af V
F
Front and Rear

Slip Angles
Vehicle
path
a =ar
 A neutral steer vehicle has
f
Side the rear slip angle equal to
force a Neutral the front slip angle. So the
r V steer
R vehicle simply moves
sideways and forward
without rotation.

Conclusion
 The effect of any changes in the vehicle
specifications on the handling behavior
can be judged by examining the changes
induced in the front and rear slip angles.
 Any change
 increasing the front slip angle will thus
increase the understeer tendency,
whereas
 increasing rear slip angle will add an
oversteering tendency.

Suspensions
and
Anti-roll Bar
(ME 465)

Leaf
Damper Spring
Beam axle Rigid Beam Axle

Suspension
Differential
Leaf Helical
Spring Coil
Spring
Panhard
Rod
A-bracket
Longitudinal
arm
Longitudinal
arm

Strut
McPherson Strut Type
Helical
coil
spring
Independent Suspension
Damper
Track
control
arm

Double Wishbone
Upper
wishbone
Type of Independent
Suspension
Anti-roll bar
Upper
wishbone
Lower
wishbone
Anti-roll bar

Anti-Roll Bar
Helical
coil
Damper spring
Anti-roll bar
 Rigid beam axle suspension with anti-roll bar.

Anti-Roll Bar
 Anti roll bar
connects one side of
a suspension with
the other and acts in
torsion.
 It is attached to the body such that it can rotate in

rubber bushings. It is mounted across the vehicle with
ends attached to the suspension arms on each side.
 Anti-roll bar can be fitted to the front suspension, rear
suspension, or both.

Anti-Roll Bar
 Anti roll bar may have a
significant influence on
the handling behavior of
a vehicle.

Anti-Roll Bar
 The function of the anti-roll bar is to
reduce the roll motion of the vehicle
body.
 During equal vertical wheel travel, the
anti-roll bar simply rotates about the
bushings at the body attachment points
without twisting.
 Thus the anti-roll bar has no effect on
the body bounce motion, front end lift
during acceleration, or dive during
deceleration.

Anti-Roll Bar
 The addition of an anti-roll bar to a front or
rear suspension originally without one, or
replacement of the existing one with a
stiffer one will increase the slip angle of
the tires.
 Thus, an additional anti-roll bar
- at the front will introduce an
understeering tendency,
- at the back will introduce an
oversteering tendency.

Example
 What would be the effect on vehicle handling of
 underinflated front tires
for an originally neutral steering vehicle?
 Ans.:
Underinflated front tires will result in reduced
cornering force characteristics. Thus for the
same cornering maneuver, the required cornering
force can only be generated at a higher slip angle
which in turn requires a larger steering angle.
This means the term af - ar will get larger

introducing a tendency towards understeering.
Exercises
1. What would be the effect on vehicle
handling of
 front wheel drive as opposed to
rear wheel drive,
 overloaded luggage compartment,
 overinflated rear tires,
 introduction of a positive camber
to front (or rear) wheels,
 load transfer from inside to
outside wheels during cornering?

Exercises
2. Consider a passenger car equipped with
front and rear anti-roll bars. The handling
tests showed that the vehicle rolls very
little but understeers rather heavily.
What would you do to improve handling
characteristics of the vehicle, if anti-roll
bars with a range of diameters are
available? Explain briefly and clearly
indicate the best choice if you have more
than one solution.

Vehicle Dynamic
Testing and Simulation

Skid
Pad
Testing
Skid pad testing is an important part of a new vehicle

design. Valuable information on tire characteristics,
traction, and vehicle handling can all be determined
through skid pad testing.
Skid Pad Testing
According to ISO 4138,

the vehicle is driven
through a turn of a known
diameter at increasing
speeds until the test
driver cannot keep the
vehicle in the turn*.
* SkidPadTest_9s.wmv
The vehicle under/over-steer behavior, roll angle at steady state, steering

as a function of lateral acceleration, and side slip as a function of lateral
acceleration can be obtained by this test.
It can be performed using constant radius or steering angle at varying
speed, or at constant velocity with varying steering angle.

The vehicle is first driven at test speed on a straight line. A step
steering wheel input is then applied as rapidly as possible.

The vehicle is driven at the test speed on a straight
line. It must then follow the designed route using
steering wheel angle and constant throttle.
The consists of an entry
and an exit lane each with a length of 12m and a side lane with a
length of 11m. 2m after the start of the entry lane, the throttle is
released so that the entire maneuver is completed with the top gear
and an engine speed of at least 2000 rpm.
At the end of the entry and exit lane the velocity is measured. The
entry velocity is increased gradually. If no cones are overturned, the
test is passed.

Passenger cars - Vehicle dynamic simulation and
validation - Steady-state circular driving behavior
specifies a method for comparing computer
simulation results from a vehicle mathematical
model to measured test data for an existing vehicle
according to steady-state circular driving tests as
specified in ISO 4138 or the Slowly Increasing Steer
Test that is an alternative to ISO 4138.
It is applicable to passenger cars as defined in ISO
3833.

GM Proving Ground

Common Mistakes
to be avoided!

Critical Mistakes
 A mistake you should avoid:
Critical speed is defined only for an oversteer

vehicle!
- Thus, do not talk about critical speed for an
understeer vehicle!
Characteristic speed is defined only for an

understeer vehicle!
- Thus do not talk about characteristic speed
for an oversteer vehicle!

Critical Mistakes
 Another mistake you should not make :
As defined in this course, the cornering

stiffness in vehicle models is always a
negative quantity.
Thus, when entering a value for the cornering

stiffness in a vehicle model equation, place a
‘minus’ sign in front of it, even if it is specified
as a positive value.
Note, however, that this sign convention does
not apply for nonlinear tire models such as
Dugoff and Magic Formula models.
Critical Mistakes
 Basic assumptions for the linear bicycle model:
- Lateral acceleration < 0.3, otherwise roll

degree of freedom must be added!
- Slip angles < 4o, otherwise linear tire

characteristics cannot be used and the use
of a nonlinear tire model is necessary!

Critical Mistakes
 Lateral acceleration is not v

Lateral acceleration is: v +ur
 For all equations either
- Compatible units must be used,

or
- Unit conversions must be made.

Useful Matlab codes

Example
 Write down the equations m  v + Ur  = Fyf + Fyr
of motion for the bicycle
model in the state Jr = aFyf - bFyr
equation form.
Ans.:  m 0   v   0 U   v   1 1   Fyf 
 0 J   r  =  0 0   r    a  b   F 

         yr 
1 1 
 U
 v   0  v   m m   Fyf 
  m      
 r   0 0   r   a  b   Fyr 
   J J 

2
Input Generation 1.5
Steeringangle (deg)
0.5
% input: full sine with delay 0
% for lane change maneuver -0.5
-1
input_amp=1.9 % [deg]; -1.5
delay1=1; [s]; -2
0 1 2 3 4 5 6 7 8 9 10
Time (s)
input_duration=8; % [s]
deltad=input_amp*sin(pi*(t-delay)/(input_duration)*2); % [deg]
deltad(t<delay)=0.0; % [deg]
deltad(t>input_duration+delay)=0.0; % [deg]
deltar= deltad*pi/180; % [rad]

Numerical Solution of LTI
System Response
% lsim Solution for linear bicycle model
sys = ss(A,B,C,D); % state matrices
[x,T] = lsim(sys,deltar,t];
v = x(:,1); r = x(:,2);
sys = system definition,
t = time vector,
deltar = input vector

EXAMPLE: Plot the trajectory of vehicle cog.
% Solve state equations for v and r, integrate r to obtain
% heading angle lamda using t=0:deltat:tmax
% Numerical integration by trapezoidal rule :
lamda=cumtrapz(r)*deltat;
% Calculate xdot and ydot.
xdot=U*cos(lamda)-v.*sin(lamda);
ydot=v.*cos(lamda)+U*sin(lamda);
% Integrate xdot and ydot to find position of vehicle cog
% relative to inertial reference frame.
xt=cumtrapz(xdot)*deltat;
yt=cumtrapz(ydot)*deltat;
% Plot the trajectory of the vehicle cog
plot(xt,yt,'b.')
grid

EXAMPLE - ODE & lsim Solutions for Bicycle Model
% ODE & lsim solutions for the linear bicycle model
% --------------------------------------------------------------------------------
close all; clear; clc
% --------------------- ONLY FOR ODE SOLUTION ----------------------
% Define global variables
global A B stepA
% --------------------------------------------------------------------------------
Ukph = 60; % Forward speed of the vehicle [kph]
U = Ukph/3.6; % [m/s]
% --------------------------------------------------------------------------------
% Vehicle specifications
m = 1400; % [kg]
J = 1750; % [kg.m2]
a = 1.42; % [m]
b = 1.03; % [m]
% --------------------------------------------------------------------------------
% -----------------------------------------------------------------------------------
% Tire cornering stiffness – here defined for each tire
Cft = -37500; % [N/rad/tire]
Crt = -50000; % [N/rad/tire]
Cf = 2*Cft; % [N/rad/front tires] – right and left tires
Cr = 2*Crt; % [N/rad/rear tires]
% -----------------------------------------------------------------------------------
% Step steering input
stepA = 1.8; % [deg] % Amplitude of step
tf = 2; % Final time [s]
n=200; % number of points
delta = zeros(1,n);
delta(1,:) = stepA*pi/180; % steering angle [rad]
% -----------------------------------------------------------------------------------

% -----------------------------------------------------------------------------------
% States: {x}=[v r]'
% System matrix
A = [(Cf+Cr)/m/U (a*Cf-b*Cr)/m/U-U;
(a*Cf-b*Cr)/J/U (a^2*Cf+b^2*Cr)/J/U];
% Input distribution matrix
B = [-Cf/m -a*Cf/J]’;
% -------------------- ONLY FOR LSIM SOLUTION -------------------------
% Output vector
C = eye(2);
D = [0;0];
sys = ss(A,B,C,D);
% -----------------------------------------------------------------------------------

% --------------------------------------------------------------------------------------------
% SOLUTION
% ------------------------- ONLY FOR LSIM SOLUTION -----------------------------
t=linspace(0,tf,n); % linearly distributed time vector
[x,T] = lsim(sys,delta,t);
% ------------------------- ONLY FOR ODE SOLUTION ------------------------------
% Steady State Initial Conditions for the states
x_ss = zeros(2,1);
% Call function for solution
[t, x] = ode15s(‘BicycleModel_Function’,[0 tf], x_ss);
n=size(t); % time array is decided by the nonlinear solver!
delta=zeros(1,n(1,1));
delta(1,:)=stepA*pi/180; % steering angle [rad]
% --------------------------------------------------------------------------------------------

% RESULTS
v = x(:,1);
r = x(:,2);
beta = v/U;
vdot = A(1,1)*v+A(1,2)*r+B(1,1)*delta’;
alfaf =( (v+a*r)/U – delta’)*180/pi; % [deg]
alfar = (v-b*r)/U*180/pi; % [deg]
% -------------------------------------------------------------------------------------
% Plot results
figure(1)
plot(t,delta*180/pi,'m-','LineWidth',2);
xlabel(‘time [s]’)
ylabel(‘\delta [ô]’)
grid on
title('Steering Input')

figure(2)
plot(t,v,’k-','LineWidth',2);
xlabel('time [s]’)
ylabel('Side Slip Velocity [m/s]')
grid on
figure(3)
plot(t,r*180/pi,'r-','LineWidth',2);
xlabel('time [s] ‘)
ylabel('Yaw Velocity [ô] ')
grid on
figure(4)
plot(t,beta*180/pi,'m-','LineWidth',2);
ylabel('Side Slip Angle [ô]')
grid on

figure(5)
g = 9.81 % [m/s2]
plot(t,(vdot+U*r)/g,'m-','LineWidth',2);
ylabel(‘Lateral Acceleration [g]')
grid on
figure(6);
plot(t,alfaf, t, alfar,’--’,'LineWidth',2);
xlabel('time [s]’)
ylabel(‘Tire slip angle [ô]’)
legend(‘Front’,’Rear’,’Location’,’Best’)
grid on
% ----------------------------------------------------------------------------

% ----------------------- ONLY FOR ODE SOLUTION ----------------------

% Function for solution – NOTE that this part must be saved
as a Matlab function, separate from the main m-script file!
% -------------------------------------------------------------------------------------
function xdot=BicycleModel_SolutionFunction(t,x)
% -------------------------------------------------------------------------------------
global A B stepA
% State equation (differential equations to be solved)
xdot = A * x + B * stepA*pi/180;
% -------------------------------------------------------------------------------------
end

RESULTS:

EXAMPLE - ODE Solution for Bicycle Model
RESULTS:

EXAMPLE - ODE Solution for Bicycle Model
RESULTS:

Films &
Animations
- MAN-tga_sicherheit_ger_NEU3816.swf
- MBGüvenlikSistemleri_1.flv
- RussianTunnel.wmv

ME 513
VEHICLE
DYNAMICS

Prof. Dr. Y. Samim
Ünlüsoy
CH III - Articulated
Vehicles

 Rigid Truck + Full Trailer
Articulated
Vehicles
 The articulated vehicle is made up  Tractor + Semitrailer

of two or more units which are
hitched to each other in series.
 The tractor is a conventional two
or three axle vehicle which
supports the front end of a
 Car + Trailer
semitrailer on a platform, or 'fifth (Caravan)
wheel' which has freedom to
rotate in yaw and pitch; the rear
end of the semitrailer is mounted
on an axle.
 Tractor + Semitrailer
Fifth Wheel
The tractor is a conventional two or three axle vehicle

which supports the front end of a semitrailer on a
platform, or 'fifth wheel' which has freedom to rotate
in yaw and pitch; the rear end of the semitrailer is
mounted on an axle.

Fifth Wheel
The tractor is a conventional two or three

axle vehicle which supports the front end
of a semitrailer on a platform, or 'fifth
wheel' which has freedom to rotate in yaw
and pitch; the rear end of the semitrailer
is mounted on an axle.

Articulated Vehicles

Articulated Vehicle Models
 In most of the studies single track models
of articulated vehicles are developed. These
models are essentially the extensions of the
bicycle model applied to vehicles
consisting of two or more units hitched
together.
 In these models, plane motions only are
considered, and roll motion is not included.
it is usual to start with a single virtual tire
representing each axle of the semitrailer
which may be lumped together in case
further simplification is considered.

Articulated Vehicle Models
 First studies involve linearized models

where the tire characteristics are
represented by cornering stiffnesses,
assuming small slip and trailer (hitch)
angles.
 These models can then be extended to
nonlinear models which involve large
trailer angles and/or large slip angles
requiring the use of nonlinear tire
characteristics.

Modeling Articulated Vehicles
 As at least two units are involved, a

minimum of three degrees of freedom
are required to describe the system
when travelling at a steady forward
speed.
 However, the classical example of a
one degree of freedom semitrailer
model provides an interesting starting
point for more detailed studies.

Single degree of
freedom model

d
Single
O
m t , Jt
c Degree of
y U
Freedom
H
Model
 The classical example of the simple one degree of

freedom semitrailer model disregards the effects
of semitrailer motion on the towing vehicle.
 Therefore, it represents the case where a light
weight semitrailer is towed by a very large vehicle,
i.e., of infinite mass

d
Single
O c Degree of
m t , Jt
Freedom
y U
H Model
 The hitch point H has a constant forward velocity U,

which is the velocity of the towing vehicle.
 mt is the semitrailer mass and Jt is the yaw moment of
inertia of the semitrailer about its center of gravity.
 The only degree of freedom is the articulation (trailer)
angle y.
 The single wheel in the model is intended to represent all
the wheels on the semitrailer axle.

d
O
Single Degree
c
mt , Jt of Freedom
y U Model
H
c  d y  c + d  ψ + Usinψ dom

U
y  t
U siny
U cos y Ucos y
 c + d ψ
 + U sin ψ For small c+d
tan  -α t  = trailer & slip -α t =  +ψ
ψ
U cos ψ U
angles

Fc
Single
d
Degree of
c Freedom

J Hy
y
Model
H
The equation of motion can be obtained by taking moments

about the hitch point, H .
  c  d Fc  0 J H  J t  mt c2
JH y
Fc  C t  t   c  d C t  t  0
JH y
 2
c + d  Ct c+d
 -
JH ψ  -  c + d  Ct ψ = 0
ψ -α t =  +ψ
ψ
U U
2
 c + d  Ct
JH  -
ψ -
ψ  c + d  Ct ψ = 0
U Single Degree
 2
c + d  Ct  c + d  Ct of Freedom
 -
ψ -
ψ ψ= 0
JH U JH Model
x  2 n x  n2 x  0
- Thus the motion of the trailer is governed by a second

order differential equation with constant coefficients.
- Such a motion will be stable provided that all the
coefficients are of the same sign.
- Remembering that the cornering stiffness, Ct, is by
definition a negative quantity, this condition is satisfied.

2
 -
ψ
 c + d  Ct
-
ψ
 c + d  Ct ψ= 0 Single Degree
JH U JH
of Freedom
x  2n x  n2 x  0 Model
n  
c  d C t
Undamped Natural
JH Frequency
ζ=
 c + d
-
 c + d  Ct
Damping Ratio
2U JH
- Effects of parameters U, Ct, JH, c, and d on relative

stability of the system ?
- Note that mt does not appear in the equations, but it
is included in JH !


 2
c + d  Ct

 c + d  Ct Single Degree of
ψ- ψ- ψ= 0
JH U JH Freedom Model
x1  ψ x 1  x 2
(c  d) Ct (c  d) 2 Ct

x2  ψ x 2  x1  x2
JH J HU
 0 1 
 x 1     x1 
  =  (c + d) Ct
2
(c + d) Ct   
x 2    x2 
 JH J H U 
Matlab command “damp” can be used to obtain natural frequency
and damping of LTI model poles.
[Wn, Z] = damp(A)
returns vectors Wn and Z containing the natural frequencies and
damping factors of the LTI model with the system matrix A.

Example
Data for 6 semitrailers:
Semitrailer No. 1 2 3 4 5 6
mt [kg] 1823 365 906 2160 23354 19130
Jt [kg.m2] 6100 542 830 7759 159583 132000
c [m] 0.38 0.15 0.213 3.87 5.45 3.24
d [m] 4.27 2.29 2.387 0.43 4.28 2.96
Ct [N/rad] -124500 -89000 -44480 -73400 -427000 -460000
JH [kg.m2] 6363.2 550.2 871.1 40109 853260 332820
fn [Hz] 1.52 3.16 1.84 0.45 0.35 0.47
 [U=100 kph] 0.89 0.87 0.54 0.22 0.39 0.33

1 2 3
4 5 6
Response to 5 degrees initial trailer angle

Conclusions
 The single degree of freedom model cannot predict
any instabilities. The deficiency of this model is that it
assumes a semitrailer mass negligibly small in
comparison to the mass of the tractor. Thus the
interaction between the tractor and the semitrailer is
neglected.
 It can still give some indications as to the stabilizing
or destabilizing effects of the parameters involved.
 In most cases of interest, however, the mass of the
semitrailer is comparable to that of the tractor if not
greater.
 Therefore, a more complicated model including the
interactions of the two parts of the combination is
needed.
Three Degree of Freedom
Tractor-Semitrailer Model

 The tractor unit is modelled
using the bicycle model. Three Degree of
 The semitrailer is represented Freedom Model
by a single equivalent wheel.
e
c b a
d y H m, J u
u d
mt , J t t r
dom r F dom
rt 
O v cf f
F
dom v cr
t t
Fct
 The units are connected at the hinge point H, called the

fifth wheel.
 Angle between the tractor and the semitrailer is called the
articulation angle and is denoted by y.
Three Degree of Freedom Model
 One can write three equations for each of the units
coupled by the forces X and Y at the fifth wheel (hitch
point).
X H X H
tractor tractor
On the On the
y y
Semitrailer Tractor
semitrailer Y semitrailer Y

Three Degree of Freedom Model
e
b a
X H m, J u
y d
r r F
dom 
v cf f dom
Y Fcr
mu  vr    X  Fcf sin d

Equations of
motion for the mv  ur   Fcf cos d  Fcr  Y
tractor
Jr  aFcf cos d  bFcr  eY

c Y
d X y H
mt , J t u
t
Three Degree of
O
dom V
rt Freedom Model
t t
Fct
m t u t  v t rt   X cos y  Y sin y

Equations of
motion for the m t v t  ut rt   Fct  X sin y  Y cos y
semitrailer
J t rt  dFct  cY cos y  cX sin y

v
c rt
Three Degree of
H u Freedom Model
compatibility relations
y
Velocity of point H must be
ut vt the same when considered in
er
the truck or in the semitrailer.
The accelerations can
u t  u cos y  v sin y  er sin y
now be obtained by
v t  crt  u sin y  v cos y  er cos y differentiating these
equations w.r.t. time.
u t  u cos y  uy sin y  v sin y  vy cos y  er sin y  ery cos y
v t  u sin y  uy cos y  v cos y  vy cos y  er cos y  ery sin y  crt

e
c b a
d y H m, J u d
u dom
mt , J t t r
dom r F 
rt v cf f
O Fcr
dom v
t t
Fct
y  y o   r  rt  d
t
0
Three Degree of
dy Freedom Model
y   r  rt
dt Articulation (trailer) angle y
  r  rt
y is an important variable.

ut
d rt ut
Slip Angles O O
t
dom
vt
y y
v  ar vt  d r t
f  d
u
For the tractor
v  br
r 
u
v t  drt
For the semitrailer t 
ut
u sin y  v cos y  e r cos y  (c  d ) rt

t 
u cos y  v sin y  e r sin y

Linearization
For small trailer angles, take
cos y = 1 and
sin y = y, and
similarly for small steering angles, take
cos d = 1 and
sin d = d
Neglect all products of small variables
(such as vr, vtrt, ry, vy, etc., but not uy ! )

Linearization
u t  u cos y  v sin y  er sin y ut  u
v t  crt  u sin y  v cos y  er cos y v t  u ψ + v - er - crt
v t  u sin y  uy cos y  v cos y  vy cos y  er cos y  ery sin y  crt
v t  u y  v  er  crt
ψ  r  rt v t  v - er - crt + u  r - rt 

Linearization
After linearization, the first equations for the

tractor and the semitrailer are uncoupled from
the rest. They do not include lateral motion
and are related mainly to longitudinal motion.
Thus, one can leave the two equations for
longitudinal motion aside and get involved
with the second and third equations for the
tractor and the semitrailer.

Linearization
m v  ur   Fcf  Fcr  Y
For the tractor
J r  a Fcf  b Fcr  e Y
m t v t  u t rt   Fct  Y For the semitrailer
J t rt  d Fct  c Y
Y  m t v t  u t rt   Fct The hitch force is solved from

the first equation for the
semitrailer.
Y  mt ur  mt v  mt er  mt crt  Fct
WHY?
Exercise
m v  ur   Fcf  Fcr  Y
For the tractor
J r  a Fcf  b Fcr  e Y
m t v t  u t rt   Fct  Y For the semitrailer
It is possible to solve for the hitch force from any one

of the above equations to obtain three equations
independent of Y.
Do you think you will get the same vehicle behavior in
each case? If your answer is “No”, what will be your
justification for this answer?

Y  m t Ur  m t v  mt er  mt crt  Fct
Linearization
That leaves three equations for
the combination in terms of the
variables :
m  v + Ur  = Fcf + Fcr - Y
v
J r  a Fcf  b Fcr  e Y r, and
rt
However, articulation (trailer )

ψ  r  rt angle, y, is an important
variable and usually a fourth
equation is also added.

Linearization
Fcf  Cf  f
Fcr  Cr  r Cornering forces
Fct  Ct  t
v - er -(c + d)rt
αt  ψ + Linear form of trailer slip
U angle

Linearized equations of motion
Mx   Ai x Bi d
m  m t  e mt  c mt 0 v   Cf 
  em J  e 2m t 0
r   a C 
M    t e c mt
x    B i    f

  c mt e c mt J t  c 2m t 0 rt   0 
   y   0 
 0 0 0 1
 Cf  Cr  Ct aCf  bCr eC  c  d  Ct 
   m  mt  U  t  Ct 
 U U U U 
 aC  bC
r  eCt a2Cf  b 2Cr e 2 Ct  c  d  eCt 
 f  em t U  eCt 
 A i    U U U U U 
 
   c  d  Ct  c  d  eCt  c  d 2 Ct
 cmt U (c  d)Ct 
 U U U 
 0 1 1 0 
 

 M x  =  A i x + Bi  δ State equations
Pre-multiply all
terms with the x   M 1 A i x M 1Bi d
inverse of [M],
to get the standard A  M 1 A i 
state space form.
B M 1B i 
x  =  A x + B δ

Stability Analysis

Stability Analysis
One can examine the stability of this
system by applying Routh-Hurwitz criteria
to the characteristic equation
sI  M 1A i  sI  A  0
which will be in the general form :

4 3 2
a4s  a3s  a2s  a1s  a0  0
Coefficients of the characteristic polynomial can be

obtained by MATLAB command ai = charpoly(A).

a4s4  a 3s3  a 2s2  a1s  a0  0
Stability Analysis
A necessary but not sufficient condition for
stability (Hurwitz criterion):
ai > 0 i = 0, 2, 3, 4
It can be shown* that a4 = 1 and a3 is always
positive for the tractor-semitrailer model.
However, a0, a1, and a2 may become negative
above some critical forward speed.
* M.F.J. Luijten , Lateral Dynamic Behavior of Articulated

Commercial Vehicles, M.Sc. thesis, Department Mechanical
Engineering, Eindhoven University of Technology, 2010.

a4s4  a 3s3  a 2s2  a1s  a0  0
Stability Analysis – Routh Hurwitz Criteria

s4 a4 a2 a0 Stability criteria:
ai > 0, i = 1,2,3,4
s3 a3 a1 0
(Necessary but not
a2a3  a1a4 sufficient!)
s2 a0 0
a3 a2a3 - a1a 4 > 0
a1a2a3  a0a23  a4a12
s1 0 a1a2a3 - a0a32 - a 4a12 > 0
a 2a 3  a1a4
s0 a0 0

a4s4  a 3s3  a 2s2  a1s  a0  0 Stability Analysis
It should be noted that when the coefficient
of the zeroth power of s changes sign, i.e.
a0<0 ,
a nonoscillatory instability will be indicated.
An oscillatory instability, on the other hand,

may be associated with a1<0 and a2<0, as
well as for the case*:
a1a 2a3  a0a23  a4a12  0
Note: if all the coefficients are positive, then

a 2a 3 - a1a4  0 this condition is contained in the above
condition, and thus need not be checked.
* Review Routh-Hurwitz stability criteria, Hurwitz determinants.

sI  M 1 A i  sI  A  0 Stability Analysis
Another approach for stability analysis
might be the examination of the roots of
the characteristic equation
a4 4  a3 3  a 2 2  a1  a0  0

Matlab which are, in fact, the eigenvalues of the
commands: matrix [A].
eig (A) If [M] is singular, then it is necessary to
or use an algorithm to find the generalized
eig (Ai, M) eigenvalues of the matrices [Ai] and [M].

Stability
Analysis
Reduce to a single
virtual tyre
EXAMPLE
Plot the roots of the characteristic equation, for the

tractor-semitrailer combinations specified below,
- in the forward speed range from 20 to 120 [kph]
- with 10 [kph] increments.
Comment on the stability of the combination.

Stability Analysis
Parameter Unit 1 Unit 2
Tractor m [kg] 1562 1496
J [kgm2] 2685 3004
a [m] 1.271 1.25
b [m] 1.411 1.55
e [m] 2.804 2.83
Cf [kN/rad] 33.8 29.0 axle
Cr [kN/rad] 49.8 64.0 axle
Semitrailer mt [kg] 906 2160
J [kgm2] 830 7759
c [m] 0.213 3.87
d [m] 2.387 0.43
Ct [kN/rad] 44.48 73.4 axle

Stability
4
3
Analysis
2
1
Imaginary Axis
0 Unit 1
-1 Tractor only
-2
-3
-4
-25 -20 -15 -10 -5 0
Real Axis
Thus the tractor is stable in the forward speed range of

interest.
15
Stability
10
Analysis
5
Imaginary Axis
0
Combination
-5 ≈35 kph
-10
Unit 1
-15
-140 -120 -100 -80 -60 -40 -20 0 20 Nonoscillatory
Real Axis
instability
(a0) changes sign at that speed !

5
3 Stability
2
Analysis
Imaginary Axis
-1 Unit 2
-2 Tractor only
-3
-4
-5
-35 -30 -25 -20 -15 -10 -5 0
Real Axis
Thus the tractor is stable in the forward speed range of

interest.
2.5
2
Stability
1.5
Analysis
1
Imaginary Axis
0.5
Combination
0
-0.5
75 kph
-1
-1.5
Unit 2
-2
-2.5
Oscillatory
-20 -15 -10
Real Axis
-5 0 5
instability
(a2 a3a4 - a1a42 -a22a5) changes sign at that speed !

Steady State Behavior

Steady State
For the steady state motion of the combination, all
the derivatives of v, r, and y are set equal to zero.
Thus, for a constant forward speed of U for the
tractor and the semitrailer:
dy
 y  r  rt  0
dt
r  rt
The expression for the trailer slip angle becomes
v  c  d  e r
t  y 
U

Steady State
m Ur  Fcf  Fcr  Y
For the tractor
0  a Fcf  b Fcr  e Y
m t Ur  Fct  Y
For the semitrailer
0  d Fct  c Y
The hitch force is solved

from the first equation for
the trailer.
Y  mt Ur  Fct

an 1sn  ansn 1  an 1sn  2  ...  a1s  a0  0 Steady State
 C f  Cr  C t a Cf  b Cr cde 
  m  m t  U  Ct Ct 
U U U
   v   Cf 
 a Cf  b C r  e C t a C f  b Cr e (c  d  e )
2 2
   
 e mt U  Ct  e C t   r   a C f  d
 U U U U    
 (c  d ) (c  d ) (c  d  e )  y   0 
  Ct c mt U  Ct  c  d  Ct 
 U U 
The forward speed above which possible nonoscillatory

instability of the tractor-semitrailer unit occurs can be
determined by setting the determinant of the above
matrix equal to zero*. This means :
a0  0
* See exercise 1) in Handling 2

Steady State
To simplify the manipulations, one can
divide the third column by U, and
subtract from the first column.
Further simplification results if
the third column is multiplied by (c+d+e)/U, and
added to the second column.
Thus,
Cf  Cr a Cf  b Cr    
m mt U Ct
U U
a Cf  b Cr a2 Cf  b2 Cr   e Ct  0
e mt U
U U
0 c mt U  c  dCt

an  1sn  an sn  1  an  1sn  2  ...  a1s  a0  0 Steady State
Handling
Behavior
Thus the determinant will give the coefficient of the
zeroth power of ‘s’, i.e., a0 in the case articulated
model here, of the characteristic polynomial. Set it
equal to zero to get:

U 2  m(aC f  bC r )  m t
d
(a  e)Cf  e  b Cr    (a  b) 2 Cf Cr  0
 c  d  
Note : semitrailer axle cornering stiffness does not

appear in this equation !

Steady State Handling Behavior
Thus, an expression for the forward speed at

which the coefficient a0 changes sign, i.e. vehicle
becomes unstable, can be obtained.
 ( a  b ) 2
Cf Cr
U cr 
2
d
m(aCf  bCr )  m t (a  e)Cf  e  b Cr 
cd
Note that this critical speed is associated with a
nonoscillatory instability. Oscillatory instability
can be obtained from the previous stability
analyses using the transient model.

 ( a  b ) 2
C f Cr
2
U cr 
d
cd
Nonoscillatory instability
It is obvious that this critical speed is real only if the
denominator of the expression is negative. Thus, the
condition for the combination to have a nonoscillatory
instability is given by:
d
m(aCf - bCr ) + m t (a + e)Cf +  e - b  Cr  < 0
c+d

 ( a  b ) 2
C f Cr
2
U cr 
d
cd
Nonoscillatory instability
Further, it is clear that this critical speed is
imaginary if the denominator of the expression is
positive. Thus, the condition for the combination
not to have a nonoscillatory instability is given by:
d
m(aCf - bCr ) + m t (a + e)Cf +  e - b  Cr  > 0
c+d

Articulated Vehicle
Instabilities

Articulated Vehicle Instabilities
The best known types of yaw instability of articulated
vehicles are the so called jack-knifing and trailer swing.
Analyses of jack-knifing and trailer swing require the
complete set of equations of motion including
- the longitudinal dynamic equilibrium and traction
together with
- nonlinear tire characteristics representing tire
behavior up to the point of skidding.
It is common to come across studies, however, in which
poor relative stability is also accepted to indicate these
rather terminal cases of instability.

Articulated Vehicle
Instabilities
In case of jack-knifing, the

tractor rotates suddenly and
uncontrollably around the fifth
wheel with blocked rear axle
wheels during braking.
Articulated
Vehicle
Instabilities
In trailer swing, the semitrailer

rotates suddenly and uncontrollably
around the fifth wheel with blocked
rear axle wheels during braking.
Articulated
Vehicle
Instabilities


 Rearward Amplification (RA) refers to the fact
that at high speeds, the lateral acceleration of
the trailer unit may be higher than that of the
tractor unit in steady-state and transient
maneuvers.
 RA is defined as the ratio of the peak lateral
acceleration of the semitrailer to that of the
tractor, measured at the CG points during a
single lane change maneuver. A higher RA
ratio increases the risk of trailer rollover,
causing the towed unit to rollover in severe
driving conditions.

Lateral Accelerations for 120 km/h Lane Change

Articulated
Vehicle
Maneuverability
PFOT
Performance
The radial offset
between the trajectory
of tractor’s front axle
and that of the
rearmost trailer’s rear
axle is called Tracking
Error or Path
 PFOT of a tractor-semitrailer Following Off-
combination during a 360o turn Tracking (PFOT).

Optimization of Articulated Vehicle
Maneuverability and Stability Performance
Thus the optimization of articulated vehicle behavior

involves:
- preventing rearward amplification, i.e., reduction
of RA, and
- improvement of the low speed maneuverability,
i.e., reduction of PFOT
by passive and/or active means.

Jason Stetham running away with the Russian Mafia leader’s daughter

ME 513
VEHICLE
DYNAMICS

Prof. Dr. Y. Samim
Ünlüsoy
CH IV - FOUR WHEEL STEERING
(4WS)

Four Wheel Steering
 Fierce competition between vehicle designers
and producers results, from time to time, in the
review of previous ideas which have not been
implemented widely for some reason. Four wheel
steering (4WS) is one of these ideas.
 There is a wealth of recent research studies in the
literature concentrated on
 the possible application of 4WS to production
vehicles and
 the relative merits of various control strategies
for the determination of proper rear wheel
steering angle.

Four Wheel
Steering
1937
MercedesBenz_4WS
The idea is not new!

Four Wheel Steering
 The potential benefits from the application
of 4WS include :
 Faster and more stable response at high
speeds.
 Improved maneuverability at low speed.
 Reduced sensitivity to lateral wind gusts.

Front Wheel Steering
 If one follows the motion of a front

wheel steering vehicle, originally
travelling straight ahead, from the
instant the driver gives a steering
command, i.e., rotates the steering
wheel, till the vehicle turns round,
the following transient response will
be observed.

Front Wheel Steering
1. Steering command is given,
2. Slip angles at front tires are developed,
3. Cornering forces at front tires are developed,
4. The vehicle starts yawing and side slipping,
DELAY !
5. Slip angles at rear tires are developed,
6. Cornering forces at rear tires are developed,
7. Vehicle turns around the center of curvature,
8. Lateral acceleration is developed,
9. Vehicle settles down to a steady state turn.

Four Wheel Steering
- Thus conventional front wheel steering
vehicles must develop a sideslip angle at
the center of gravity to be able to generate
slip angles, and hence cornering forces, at
the rear wheels.
- Four wheel steering vehicles, in contrast,
can generate rear wheel slip angles and
cornering forces by simply steering the
rear wheels. Thus improved transient
response is expected.

Four Wheel Steering
4WS technology could be traced back to the
Japanese Society of Automotive Engineers
and Technology conference in the 1960s. The
major car companies had started the
research of 4WS technology.
The first significant and longest lasting
application of 4WS was Honda Prelude which
stayed in production from 1987 to 2001.
Due to the limitations of electronic
technology, automation technology and
control technology, 4WS technology research
was temporarily shelved at that time.

Honda Prelude. Slow speed maneuvrability 1m08s.mp4

Four Wheel Steering
Since then, the 4WS technology has returned
to the people’s vision.
- In 2008, the new BMW 7 Series equipped
with 4WS system came into market. At high
speed, the steering angles of front wheels
and rear wheels are in the same phase.
While at low speed, they are in the opposite
phase. The maximum steering angle of rear
wheel is up to 3°.
- In the same year, Infiniti also launched
M37S equipped with 4WS system and it had
a great improvement on steering actuator.

Four Wheel
Steering
In April 2015, Cadillac presented CT6 using 4WS

technology in the United States show. It added a
similar independent steering mechanism driven
by a 12 V servo motor, with deluxe five link
suspension system to achieve precise control
for rear wheel steering.
Production discontinued in 2020.

Four Wheel Steering
- In May 2018, VW introduced all-wheel
steering system with new Touareg.
The vehicle is advertised to reduce turn
radius from 12.19 m to 11.19 m by
turning front and rear wheels at
opposite angles at speeds below 37 kph.
- Above 37 kph, front and rear wheels turn
in the same direction for smoother and
more stable driving.

ZF Active Rear Axle_1min27s.mp4

Four Wheel Steering
At present, the research on the

4WS vehicle is mainly focused on
the control strategy and its
application to achieve best
performance, i.e., shorter turning
radius at low speeds and stable
with controllable understeering
tendency at higher speeds.

Rear Wheel Steering
 Rear wheel steering (RWS) has been considered to be
unstable for high speed vehicles and as such is
restricted to low speed vehicles such as forklifts.
 A problem of maneuverability,
related to leaving a parking
position associated with rear
wheel steering has also been
another major setback.

Rear Wheel Steering
 Rear wheel steering (RWS) tends to be unstable
because, in turns the steering geometry changes
decreasing the turn radius (oversteer tendency).
 A rear wheel steered automobile exhibits non-
minimum phase behavior. It turns in the direction
opposite of how it is initially steered.
 A rapid steering input will cause two accelerations,
first in the direction that the wheel is steered, and
then in the opposite direction: a "reverse
response". This makes it harder to steer a rear
wheel steered vehicle at high speed than a front
wheel steered vehicle.

Four Wheel Steering
Vehicle Model

4WS Vehicle Model
 4WS motor vehicle as represented by the

bicycle model with yaw and sideslip
degrees of freedom .

4WS
Vehicle
f
Model
 As usual, assuming small slip angles, the

expressions for the front and rear slip
angles are obtained as :
v  ar v  br
f   f r   r
u u

4WS
f
Vehicle
Model
 A so-called vehicle sideslip angle can be

defined as the angle between the
longitudinal axis of the vehicle and the
velocity vector at the center of gravity.
v

u

Equations of Motion
If the front and rear
m(v&  ur )  Fyf  Fyr

wheel steering angles

are assumed to be
small, then the 4WS Jr&  aFyf  bFyr
vehicle bicycle model
simplifies to :
m(v&  ur )
Jr
b a
Fr Ff

m(v&  ur )  Fyf  Fyr
Equations of
Jr&  aFyf  bFyr Motion
 If slip angles are assumed to be small, linear

relations can be used.
Fyf  Cf  f m(v&  ur )  Cf  f  Cr  r
Fyr  Cr  r Jr&  aCf  f  bCr  r

m(v&  ur )  C f  f  C r  r
Jr&  aC f  f  bC r  r Equations
of Motion
v  ar
f   f
u
v  br
r   r
u
v r
mv& = (Cf + Cr ) + (aCf  bC r  mu 2 )  Cf  f  Cr  r
u u
v r
&
Jr = (aCf - bC r ) + (a Cf  b Cr )  aCf  f  bC r  r
2 2
u u

v r
mv& = (Cf + Cr )
u
+ (aCf  bCr  mu 2 )  Cf  f  Cr  r
u State
v r
Jr& = (aCf - bCr ) + (a 2Cf  b 2Cr )  aCf  f  bCr  r
u u Equations
 Cf + Cr aCf - bCr   Cf -Cr 

- u -
 v&   mu mu  
v  m  δ f 
 =  + m  
 r&   aCf - bCr a 2 Cf + b 2 Cr   r   aCf
-
bCr  δr 
  
 J J 
 Ju Ju 

v r
mv& = (Cf + Cr )
u
+ (aCf  bCr  mu 2 )  Cf  f  Cr  r
u State
v r
Jr& = (aCf - bCr ) + (a 2Cf  b 2Cr )  aCf  f  bCr  r
u u Equations
v
β  State equations can also be
u
written in terms of vehicle
&β  v& side slip angle.
u
 Cf + Cr aCf - bCr   Cf Cr 
.  - 1 - -
mu 2  mu  δf 
β
  
.=
mu    +  mu
β
 
 r   aCf - bCr a 2 Cf + b 2 Cr  r   - aCf bCr  δr 
     J J 
 J Ju 

 C f  Cr aCf  bCr   Cf  Cr 
.  u
 v 
. = 
mu mu  v  - m
 + m   f  Steady
2 2
a Cf  b Cr   r  - aCf bCr   r 
 r   aCf  bCr  State
 Ju Ju   J J 
.
 v  0  U
.=  uU r
 r  0 R
1  C C aC f  bC r  mU 2   v   C f C r   f 
 f r
2 2     
U aC f  bC r a C f  b C r   r  aC f  bC r   r 
1  ( Cf + Cr ) U aCf - bCr - mU  β   Cf
or 2
Cr   δ f 
in terms   =   
U ( aC - bC ) U a 2
C + b 2
Cr   r  aCf -bCr  δr 
of   f r f

Four Wheel Steering
Strategies

Four Wheel Steering Strategies
 The driver obviously cannot control the front and
rear wheel steering angles simultaneously.
 It is logical to let the driver steer the vehicle through
the front wheels as in the conventional FWS vehicles.
 The first and simplest idea is to make the rear wheel
steering angle proportional to the front wheel
steering angle through the relation :
r  k f
so that the driver does not have to worry about
steering the rear wheels.

r  k f Four Wheel Steering Strategies
 The coefficient k can be a constant or a

function of vehicle parameters and speed.
 A number of control strategies have been
proposed and the rear wheel steer angle is
expressed by combinations of
 front wheel steer angle,
 front wheel steer angle rate,
 yaw rate, and
 lateral velocity,
 etc.

r  k f Steady State Equations
1  C C aC f  bC r  mU 2   v   C f C r   f 
 f r
   
U aC f  bC r
2 2
a C f  b C r   r  aC f  bC r   r 
1  C C aC f  bC r  mU 2   v   C f  kCr 
 f r
2 2     f
U aC f  bC r a C f  b C r   r  aC f  k bC r 
Thus the number of control inputs is reduced to one.

r  k f Selection of a value for k
Low
speed
 The commonly used condition to
FWS, determine the expression for the
(+) parameter k, is to set the vehicle
sideslip angle equal to zero.
 This results from favorable
4WS subjective ratings received for
(0) smaller vehicle side slip angles.
 Geometrically, the forward speed
of such a vehicle will always be
tangent to the curve described by
High
speed the center of gravity of the
FWS, vehicle.
(-)
Selection of
a value for k

 U

1  C C aC f  bC r  mU 2   v   C f  kCr 
 f r
2 2   aC  k bC   f
U aC f  bC r a C f  b C r   r   f r
r  k f Selection of a value for k
 Eliminating r from the two equations, one

obtains the vehicle side slip angle. Then
v
β  U  ( Cf + kCr ) ( a Cf + b Cr ) - ( aCf - bCr - mU ) ( aCf - kbCr )
2 2 2
 =
δf δf ( Cf + Cr ) ( a2 Cf + b 2 Cr ) - ( aCf - bCr - mU 2 ) ( aCf - bCr )

v
β  U  ( Cf + kCr ) ( a Cf + b Cr ) - ( aCf - bCr - mU ) ( aCf - kbCr )
2 2 2
 =
δf δf ( Cf + Cr ) ( a 2 Cf + b 2 Cr ) - ( aCf - bCr - mU 2 ) ( aCf - bCr )
Selection of a value for k

 Set vehicle sideslip angle, , equal to
zero and solve for k.
( Cf + kCr ) (a2Cf + b 2Cr ) = (aCf - bCr - mU 2 ) ( aCf - kbCr )
a 2
bL  mU
Cr
k
b 2
 aL  m U
Cf

a 2
bL  m U
Cr  It is observed that k is a function of
k vehicle parameters and forward
b 2
 aL  m U speed.
Cf

8
aL [g] Yaw 2WS
Yaw Velocity [deg/s]

Lateral 2WS Velocity
Acceleration 6
4
0.15
r  k f r  k f
0.1
2
0.05
0 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 0
0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2
Time [s] Time [s]
 This control strategy for the rear wheels

Four improves maneuverability at low speeds.
Wheel  At high speeds, however, it results in
increased stability at the expense of
Steering heavy understeering.
Further, the vehicle will have a drop in its
r  k f

lateral acceleration gain resulting in
reduced controllability.

Generalized Formulation
 A general algorithm to relate front and rear
wheel steering angles would be written as
 r  H x  k  f *
 where [H] is a matrix to relate rear wheel angle to

state variables and k* is a parameter relating the
rear wheel steering angle to front wheel angle.
 Thus, different control strategies can be
obtained by using specific definitions for the
matrix [H] and the coefficient k*.

 C  Cr aC f  bC r   Cf    Cr 
.  f
v
.=
mu mu
 u
v  - m
 +
  m 
  f   bC  r
Generalized
 r   aC f  bC r a 2 C f  b 2 C r   r  - aC f   r 
 Ju Ju   J   J  Formulation
x&  =  A x + Bf  δf + Br  δr

δr =  H x + k*δf
x =  A  + Br  H  x + Bf + k Br δf

& *
 
Note that system matrix [A] is modified !

Note on Handling Behavior
     r     f
&
x =  A + B H  x + B + k *
Br δf  
 Since the system matrix [A] is modified by an
additional term, the stability characteristics of
the 4WS vehicle may be different than the
classical bicycle model!
 It is now possible for the second coefficient of
the characteristic equation to change sign,
resulting in a oscillatory instability.
 Thus the definitions (i.e., oversteer, neutral steer,
understeer) may be inadequate to classify the
handling behavior of the 4WS vehicle.

A Selection of Control
Strategies

Strategy 1 : Front Wheel Steer
(FWS) (2WS)
 H  = 0 0 This case is to
be used as
the reference
*
k 0 for evaluating
various 4WS
strategies.
r  0

 H  0 0 Strategy 2 :
a 2 Zero Steady State

bL  m U
* Cr Vehicle Sideslip
k 
b 2 Angle
 aL  m U
Cf
This control strategy ensures zero sideslip

angle during steady state (only) cornering.

Strategy 3: Rear Wheel Angle as a
function of Front Wheel
Angle only
 H = 0 0
*
k  a specified function of  f
In this case, the front wheel angle determines the
rear wheel angle according to the specified
function.

Strategy 3: Rear Wheel Angle as a function
of Front Wheel Angle only
 H = 0 0 k *  a specified function of  f
 At low speeds involved in city traffic or parking, large

front wheel steering angles are required.
 At high speeds, small front wheel steering angles are
sufficient for lane changes or negotiating gentle curves.
 Therefore, the appropriate function must provide for
steering of the front and rear wheels in the same
direction with small steering inputs and in the opposite
directions with large steering inputs.

Strategy 3 : Rear Wheel Angle as a
function of Front Wheel
Angle only
δr = f ( δ f )

Strategy 4 - version 1:
Proportional plus yaw velocity feedback
δr = CUr - δf  m b a  
 H  = 0 -  +  U 
* 
 L  Cf Cr  
k = -1
 This control strategy ensures (almost) zero
vehicle sideslip angle during steady state
cornering.
 Further, the 4WS vehicle will always be stable
irrespective of the forward speed and the handling
characteristics (neutral, under, or oversteering) of
the original 2WS vehicle.
* See example 3 for proof
Strategy 4 - version 2 :
Proportional plus yaw velocity feedback
*
δr = Cr + k δf  (aCf  bCr )  mU 2 
 H = 0 
* Cf  Cr U 
k =-
Cr
 This control strategy will give zero vehicle sideslip

angle even in the transient part of cornering.
* See Example 2 for proof.

Strategy 5 : Yaw velocity feedback
δr = CUr k *  0  H  = 0 m b a  
 -  U
 L  Cf Cr  
 An interesting property of this control strategy is

that it gives neutral steer characteristics to the 4WS
vehicle when applied, irrespective of the forward
speed and the handling characteristics (neutral-,
under-, or oversteering) of the original 2WS vehicle.
 It is to be noted that, the resulting vehicle will be
stable for the originally neutral- and oversteered
vehicles. For the originally understeered vehicles,
however, the vehicle will be unstable above a certain
critical speed.
Strategy 5 : Yaw velocity feedback
 m b a  
δr = CUr  H = 0 L  C - C  U  *
k 0
  f r  
Can you explain the reason behind the apparent

contradiction between the statements :
 “This control strategy gives neutral steer
characteristics to the 4WS vehicle.”
 “For the originally understeered
vehicles, the vehicle will be unstable
above a certain critical speed.”

Examples

Example 1:
Using a ramped step input steering wheel angle
 sw  0.501 * 1  tanh (25t  3 )

plot the time variation of the vehicle sideslip angle,
yaw velocity, and lateral acceleration for the five
control strategies at a vehicle speed of 90 kph.
J Cf Cr
m [kg] a [m] b [m]
[kg.m2] [N/rad] [N/rad]
1310 1760 1.2 1.4 -34870 -31730
Note that this is an understeering vehicle.

1.4  Ramped
1.2
δsw = 0.501* 1 + tanh ( 25t - 3 )  step input
steering
wheel
Front Wheel Angle [deg]
1
angle
0.8
 It is good
0.6 practice to
plot the
0.4
input first
0.2 and make
sure that it
0
0 0.5 1 1.5 2 is correct!
Time [s]

 Rear wheel angle as a
function of front wheel
steering angle for
Strategy 3 :
𝟐 𝟑 𝟒
𝐫 𝐟 𝐟 𝐟 𝐟 𝐟

Vehicle side slip
angle, 

12
Yaw Strategy 5
10
Yaw Velocity [deg/s]
Velocity
2WS
8
Strategy 3
6
Strategy 4 v2
2
Strategy 2
0
0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2
Time [s]

0.45
0.4 Lateral Strategy 5

Lateral Acceleration [g]
0.35 Acceleration
0.3 2WS
0.25
0.2 Strategy 3
0.15
0.1 Strategy 4 v2
0.05 Strategy 2
0
0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2
Time [s]

Example 2:
Obtain the general form of a
‘Proportional + yaw velocity feedback’
control strategy, i.e. expressions for the C
and k*, that will give the 4WS vehicle zero
steady state sideslip velocity.
δr = Cr - k * δf  H  = 0 C
v  *
δr =  H    + k δ f
r 

 Cf + Cr aCf - bCr   Cf -Cr  In such problems,
- U -
 v&   mU mU v   m m   δf  you are advised to
 =   +  
 r&   aCf - bCr
2 2 bCr  δr 
a Cf + b Cr   r   - aCf start using ai, bi, etc..
 JU JU   J J  At the end you may
insert the their
expressions from the
Write the steady state equations: state equations.
0   a11 a12   v   b11   b12 
 =    +  δf +   δr
0   a21 a 22   r   b 21   b 22 
Set v=0 in the first equation, and solve for r.

a12 b11
0  a12r  b11δf  b12δr δr   r δf
b12 b12
From the state equations:
aCf  bCr Cf Cr
a12  - U, b11  , b12 
mU m m

aCf  bCr Cf  Cr
a12  - U, b11  , b12 
mU m m
Insert coefficients into the equation for r.
a b
δr   12 r  11 δf
b12 b12
 aC  bC  mU 2   m   C   m 
  f r   r f
   f
 
 mU   Cr   m   Cr 
 aC  bC  mU 2   C 
 f r  r    f  f
 Cr U   Cr 
 
 aCf  bCr  mU 2   v   Cf 
r   0    f δr = Cr - k * δf
 Cr U   r   Cr 
aCf  bCr  mU 2 C
k* = - f
C
Cr U Cr

δr = CUr - δf  m b a  
 H  = 0 -  +  U 
*
k = -1  L  Cf Cr  
Example 3:
Show that the 4WS vehicle with
‘Strategy 4 – version 1: Proportional + yaw
velocity feedback’
will always be stable irrespective of the forward
speed and the handling characteristics (neutral,
under, or oversteer) of the original FWS vehicle.
 Cf + Cr aCf - bCr   Cf -Cr 
- U -
 v&   mU mU v   m m  δf 
 =   +  
 r&   aCf - bCr
2 2 bCr  δr 
a Cf + b Cr   r   - aCf
 JU JU   J J 

.
 v   a11 a12   v   b11   b12 
.=   +  δf +   δr
 r   a21 a22   r   b21   b22 
 
x&    A  x + b1 δf + b2 δr

Insert r in the state equation. δr =  H  x + k * δf
x&    A  x + b1 δf + b2 Hx  k * b2 δf

x&   (  A  + b2 H ) x + (b1  k * b2) δf
Note that the system matrix is now modified.

m b a   b12 
with  H  = 0 C C=-  + U b2 =  
L  Cf Cr   b 22 
The 4WS vehicle system matrix becomes:
 a11 a12   0 b12C   a11 a12  b12C 

 A   A  + b2 H  a a   0 b C  a a  b C
 21 22   22   21 22 22 
Obtain the characteristic equation for the 4WS vehicle.

s  a11  ( a12  b12C )
s  I    A = 0  0
a21 s  a22  b 22C
( 1) s2  ( a11  a22  b22C ) s  ( a11a22  a12a21 )  ( a11b22  a21b12 ) C  0

a1s2  a 2s  a3 0

( 1) s2  ( a11  a22  b22C ) s  ( a11a22  a12a21 )  ( a11b22  a21b12 ) C  0
a1s2  a 2s  a3 0
For instability a2 or a3 must change sign. Thus,

 ( a11  a22  b 22C )  0 or a11  a22  b 22C  0
is the first stability condition.

From the state equations:
C  Cr a2Cf  b 2Cr  Cr bCr

a11  f , a 22  , b12  , b 22 
mU JU m J
m b a 
C=-  + U
L  Cf Cr 

Insert in the stability condition. a11  a22  b 22C  0
Cf  Cr a2Cf  b 2Cr  bCr   m  b a  

       U  0
mU JU  J   L  Cf Cr  
Cf  Cr a2Cf  b2Cr  mb2   Cr   mab 

    U    U0
mU JU  
 LJ   Cf   LJ 
( ) () (+) (+)
Since
Cf < 0, Cr < 0
all the terms in the parenthesis are always negative.
Thus a2 is always positive irrespective of the handling
behavior of the original 2WS vehicle.

Now a3>0 must also be checked.
a3  ( a11a22  a12a21 )  ( a11b 22  a21b12 ) C  0
The term in the first parenthesis is the same as a3 of

the 2WS bicycle model. The effect of 4WS is
represented by the term in the second parenthesis
multiplied by C.
Insert aij expressions in the stability condition.
  a Cf  b Cr   aC  bC
2 2
 Cf  Cr r  U   aCf  bCr  
 mU    f  
  JU   mU   JU 

 Cf  Cr   bCr   aCf  bCr   Cr    m  b a  
   J    
  m  L  C    U  0
  mU     JU      f Cr  

a3  ( a11a22  a12a21 )  ( a11b 22  a 21b12 ) C  0
The term in the first parenthesis was already obtained

in simplified form for the 2WS vehicle. Simplify the term
in the second parenthesis to obtain:
 L2C C  ( aC  bC ) mU 2   LC C   m  b a  
 f r f r  f r
     U  0
 mJU 2   mJU   L  Cf Cr  
 
Carry out the

a3  L2 Cf Cr  2mbCr U 2 0
multiplications
to get: (+) (-)
Thus, it has been proved that a3 is also always positive

irrespective of the handling behavior of the original
2WS vehicle.

Exercises

EXERCISES
1. Prove that control strategy ‘Proportional Plus
Yaw Velocity Feedback’ (strategy 4, version 2)
will give zero vehicle sideslip angle even in the
transient part of cornering.
2. Show that, in case of the control strategy ‘Yaw
Velocity Feedback’ (strategy 5), the H matrix in
the generalized form can be expressed as
 m b a  
 H = 0 L  C - C  U 
  f r  
Hint: Remember that this control strategy gives neutral

steer characteristics to the 4WS vehicle when applied.

EXERCISES
3. Prove that, in the case of the control strategy
‘Yaw Velocity Feedback’ (strategy 5), the
resulting 4WS vehicle will be stable for the
originally (2WS) neutral steer or oversteer
vehicles. It will be unstable, however, above a
certain critical speed for the originally
understeer vehicles.
Derive the expression for this critical speed.
4. Prove that the 4WS vehicle with the control
strategy ‘Yaw Velocity Feedback’ (strategy 5),
will exhibit neutral steering behavior.
Hint: Consider the yaw velocity gain for neutral
steering vehicles.

ME 513
VEHICLE
DYNAMICS

Prof. Dr. Y. Samim
Ünlüsoy
CH V - Roll centers & Roll Axis

Sprung & Unsprung Masses
Sprung Mass
In vehicle dynamics
terminology, the mass of
the vehicle which is
supported by the
suspension springs is
called the sprung mass,
and the part below the
suspension springs is
called the unsprung mass.
Unsprung Mass

Roll center Approach
 In the earlier stages of vehicle design,

when the suspension system is yet to be
designed, the suspension data which
provide the basic information on the roll
motion of the vehicle are not available.
 It is, therefore, natural to confine the
handling analysis to a two degree of
freedom (yaw and side slip) vehicle
handling model such as the one
presented in the previous sections.

 The handling of the vehicle must,

however, be checked using a three
degree of freedom (yaw, sideslip, and
roll) model, as soon as the suspension
data are available.
 Such an analysis normally requires a
mathematical model of the vehicle
including the suspension system and
the solution is normally too
complicated for hand calculations.

 An approach to include the roll motion

into the simple vehicle handling model,
still keeping the analysis sufficiently
simple to obtain useful results by hand
calculation is possible through the
concepts of roll centers and roll axis.
 These concepts can be used to extend
the bicycle model or may be used in the
more complex models.

Roll centers
Roll center is defined as the point in the lateral

plane of the front or rear suspension about which
the sprung mass of the vehicle rolls under the
action of a side force, for example the centrifugal
force during cornering.

Roll
Centers
For a two-axle vehicle

there are two roll
centers, one for the front
suspension and, one for
the rear suspension.

Roll Axis and Roll Moment
 These two roll centers define a roll axis which is then
the axis of rotation of the sprung mass under the
action of a side force.
L
a b
hs
an + bm n
m L
 hs : sprung mass cog, m : front roll center, n : rear roll center

L
F
hs
an+bm
n
No Roll ?
m L
 The roll moment is given by the

product of the centrifugal force  an + bm 
M r =  hs - F
and the distance of center of  L 
gravity to the roll axis.
 Thus, it is theoretically possible to eliminate roll by
locating the center of gravity on the roll axis which
requires the roll centers to be high.
 A benefit of reducing roll is that the associated
camber change of the wheels is also reduced.

Roll Center Heights
 High roll centers - Advantages
 result in reduced roll angle, and
 as a further benefit, the associated camber

change of the wheels is also reduced; thus
handling is improved.
 High roll centers – Disadvantages
 result in considerable lateral wheel
displacements, i.e. track variations, during
wheel bump and rebound, and thus ride is
degraded during travel.

Roll Center Heights
 The actual heights of the roll centers

are, therefore, decided to give the best
compromise between the roll angle and
lateral wheel deflections.
 In practice, a certain amount of roll is
desirable to give the driver an
indication of the severity of cornering.

Roll Center Heights
 In practice, the rear roll center is

designed to be slightly higher than the
front.
 This introduces a larger load transfer at
the front, resulting in an understeering
tendency.
 WHY ?

Effect of load
transfer on
Fc kN
Cornering Force
-DW +DW
CONCLUSION
Total cornering force of
Gain an axle at a given slip
+DFc
angle decreases with load
Loss transfer!
-DFc
α [o ]
Exception:
Total cornering force of
an axle at a given slip
angle increases with load
transfer if negative
W kN
camber angle is given to
the outer wheel!

Permissible Roll Angles
 For racing cars, the maximum permissible roll angles
are specified as 3 to 4o at a lateral acceleration of
0.75 g.
 For passenger cars, permissible roll angles are of the
same magnitude but attained at lower values of
lateral acceleration. For example, in the development
of Austin Allegro, the basic overriding design
parameter was that the roll angle at a lateral
acceleration of 0.5 g should be no more than 3 1/2o.
 The current trend in passenger cars is to use fairly
low roll centers together with anti-roll bars to reduce
roll.

Permissible Roll Angles
 For racing cars, the maximum permissable roll angles
are specified as 3 to 4o at a lateral acceleration of 0.75
g.
 For passenger cars, permissible roll angles are of the
same magnitude but attained at lower values of lateral
acceleration. For example, in the development of
Austin Allegro, the basic overriding design parameter
was that the roll angle at a lateral acceleration of 0.5 g
should be no more than 3 1/2o.
 The current trend in passenger cars is to use fairly low
roll centers together with anti-roll bars to reduce roll.
Prof. Dr. Y. Samim Ünlüsoy ME 513Vehicle Dynamics 17

Roll Angle versus
 Roll angle versus

Roll Angle [o]
4
lateral acceleration
range for some
2 sports cars.
0
0 0.2 0.4 0.6 0.8
Lateral Acceleration [g]

Basic Assumption
 The roll centers are determined by the relative

positions of the suspension linkages.
 It must be kept in mind that, as a result of
suspension deflection during roll, the positions
of the roll centers are also variable.
 The basic assumption in the roll axis approach
is that the roll center locations do not change.
 In spite of this rather severe assumption, the
approach still gives an insight into the
directional behavior of vehicles.

Graphical Determination
of Roll Center Height

Rubber
mount
Spring McPherson Strut
type of independent suspension
Damper
Lateral
control
arm

A
1 McPherson
. Strut
2 M P type of
independent
3 m suspension
B
O
 The vehicle is assumed to be symmetrical with

respect to the longitudinal plane and thus the roll
center is always on the vehicle centerline of a
lateral section.

Double Wishbone
type of
independent suspension

A
1
P 2
3 B
O m
M
Double Wishbone
type of independent suspension
Exercise: What if the upper wishbone is sloping in the other

direction?
Swing Axle
type of independent
suspension
M
P
O
Trailing Arm Semi-Trailing Arm
Trailing and Semi-Trailing Rear

Suspensions

Rear view
Semi-
P
b Trailing
2 4
B 5 M Axle
O type of
ax A a independent
az
suspension
3 2
1
P
1
Top view

M
Rigid Beam Axle
with
M n Leaf Springs

Rigid Beam Axle
M with
n
A-bracket

M
Rigid Beam Axle
with
n
Panhard Rod

Calculation
of
Steady State Roll Angle

Position of Center of
W
h Gravity
 The horizontal position and
Wf a b
Wr the height of the center of
L
gravity of a vehicle can be
a
Wr
L b
Wf
L
obtained by weighing the
W W front and rear axles
 at level floor and then
 on a sloping surface.
h  The centers of gravity of the

W
unsprung masses whose
Wf Wr
weights are denoted by Wuf
and Wur are assumed to lie
at the wheel centers.

a* b*
Position of
Wuf
Wur
h h
Sprung Mass
s
rw Ws W Center of
a b Gravity
Wf Wr
L
Ws  W  Wuf  Wur
* Wr  Wur
a  L
Ws hW  rw  Wuf + Wur 
hs 
Wf  Wuf Ws
*
b  L
Ws

Roll angle of the sprung mass
 For steady state cornering at a lateral

acceleration aL, there will be tilting moments
about the roll axis are due to the action of
the inertia forces of the sprung and
unsprung masses of the vehicle.
 Similarly, restoring moments to balance the
tilting moments are produced by the
suspensions and the antiroll bars.
 Balance of these moments will decide the
steady state roll angle produced during
cornering.

Tilting Moments
 For steady state cornering at a lateral

acceleration aL, there will be tilting
moments about the roll axis are due to
the action of the inertia forces of the
 sprung and
 unsprung
masses of the vehicle.

h*sinf
WsaL
Tilting Moment of
Sprung Mass
Ws
h*cosf
h*
Mst =  WsaL  h*cosf + Ws h*sinf
f
Roll axis
For small f Mst = h* Ws  aL + f 
where * a *n + b *m
h  hs 
L

Tilting Moments of Unsprung Masses
Tilting action for
Independent
suspensions :
Front M ft  Wuf a L h uf
Rear M rt  Wur a L h ur
No tilting action
for Beam Axle M ft  M rt  0
suspensions

Mst = h* Ws  aL + f 
M ft  Wuf a L h uf M rt  Wur a L h ur Tilting Moments

 Separating the tilting moment into two parts;
 one proportional to lateral acceleration, and
 the other proportional to roll angle,
the balance of the tilting and restoring moments
about the roll axis gives :
t *
 *
 M i = Wsh + Wuf huf + Wurhur aL + Wsh f
i

t   *
 i  i L
M = K a + Ws f
h
i i 

Restoring Moments
 The restoring moments to balance the tilting
moments are produced by the suspensions and
the antiroll bars.
 The restoring moments can be expressed in the
general form where j refers to restoring moments
due to
 1-front suspension
 2-rear suspension r  
 j
M =   Cj  f
 3-front anti-roll bar, and j 
j  
 4-rear anti-roll bar
respectively, and Cj is the roll stiffness of the
particular suspension system used.

t   * Steady State
 M i    K i  aL + Wsh f
i  i 
Roll Angle
*
 K i = Wsh + Wuf huf + Wur hur
r t
i 
 j  i
M M
r   j i
 M j    C j  f
j  j 
The roll angle of the  


  K i 
sprung mass at a given i
lateral acceleration can f * L
a
now be calculated from   C j  Wsh 
the expression :  j 

CALCULATION
of
ROLL STIFFNESS

t F
M
F EXAMPLE
f D
M  Derive the
kr
expression
k ts/2
r for the roll
stiffness of
ts
a beam axle
suspension.
The force generated by each

suspension spring in response
to the tilting moment Mt, F  krD
resulting in a roll f, is given as :

t F
M
F
f D
EXAMPLE
M  For small roll
kr
angles, D is
k ts/2
r proportional to
the body roll.
ts
ts
D f
2
 The tilting moment is thus
balanced by a resisting  ts 
moment generated by the M    2F   t s F
r
suspension.  2

t F
M
F EXAMPLE
f D
M Mr  tsF F  krD
k ts/2 kr
r ts
D f
2
ts
 t 2 
 Inserting the expression for F : Mr   s k r  f
 2 
 
2
 Therefore the roll stiffness is given ts
by the coefficient : Cj  kr
2
M kr  Beam Axle
Suspension
ts
B
tr
Calculation of Roll  t2 
Cj =  s  kr
Stiffness  2 
 

b
A
kf
P
B
a
O
M
d f1
f2
 Double Wishbone Suspension
2 2
Calculation of Roll  bd   f1 
Cj = 2     kf
Stiffness  a   f2 

b
A Calculation of
P Roll Stiffness
kf
M
 McPherson Strut
B
Suspension
O
d
a 2
 bd 
Cj = 2   kf
 a 

A b Calculation of
Roll Stiffness
kf
M
P  Swing Axle
Suspension
O
a
d 2
 bd 
Cj = 2   kf
 a 

A kr
Calculation of
Roll Stiffness
P2
B M
O
d  Semi-Trailing
f1
f2 Arm Suspension
2
P1 2  f1 
C j = 2d   kr
 f2 

DYNAMIC WHEEL
LOADS

Dynamic Wheel Loads
 During cornering, there is a weight transfer
from the inside wheels to the outside
wheels.
 This weight transfer can be calculated by
summing up the contributions from
 body roll,
 sprung mass inertia force distributed
to the front and rear roll centers,
 action of the antiroll bars, and
 inertia forces of the unsprung
masses.

front
ms
Wsa
L
n
m ur
h* Dynamic
Wur a L
h ur Wheel
m
uf m
Loads
b*
Wuf a h uf
rear a* L
During cornering centrifugal forces act on the

sprung and unsprung masses.

front
ms
n
m ur
t
Ms
Dynamic
Wur a L
h ur Ws a
L
Wheel
m
uf m
Loads
b*
Wuf a h uf
rear a* L
Sprung mass inertia force, WaL, is shifted to the roll axis

by the addition of a tilting moment, Mts.
t

* *

Ms  WsaL h cos f + Ws h sin f

a*
W a
t L s L front
Msr ms
m ur
Dynamic
n b* W a
t
M sf L s L Wheel
h ur
Wur a L
Loads
m
uf m
b*
Wuf a h uf
rear a* L
- Sprung mass inertia force, WaL, is replaced

by two inertia forces acting at the roll
centers.
- Tilting moment, Mts, is replaced by two tilting
moments acting at the roll centers.

RESTORING MOMENTS
This tilting couple is
f t Cf 
resisted by the front and Mr  M s  
C
 f + Cr
rear suspensions.
Cf
  Wsh*  aL cos f + sin f  
  Cf + Cr
Resisting moments may
be distributed to the front
and rear suspensions, Mr  M t  Cr 
r s
according to the relative  Cf + Cr 
front and rear roll Cr
stiffness (including the  Wsh  aL cos f + sin f 
 * 
  Cf + Cr
antiroll bar stiffness).

t
M sf
b*
L
Ws a L Load Transfers:
Front
m
uf
huf
m 1) Due to the sprung mass
inertia force at the front
t
f
roll center :
Ri Ro=Ri+DWf
 b *   1
DWf 1   Ws a L  m 
 L   tf
 
2) Due to the restoring

moment by the front
r 1
DWf 2  M f
tf
 
suspension :

b*
t
M sf L
W a
s L Load Transfers:
W a Front
uf L m
uf
m
huf
3) Due to the action of the
t front antiroll bar :
f
 
Ri Ro=Ri+DWf
r 1 1
DWf 3  M fa   Cfaf 
tf tf
4) Due to the inertia forces

1
of the front unsprung DWf 4  Wuf a L h uf 
mass : tf

Load Transfers: Front
The sum of the four components gives the total
front load transfer.
4
DWf   DWfi
1
The tilting moment of unsprung mass in case of

beam axle type of suspensions is taken as zero
for the calculation of the sprung mass roll
angle. Note, however, that the load transfer due
to the inertia forces of the unsprung masses,
given by DWf4 is not zero.

t
M sr
a*
L
W a
s L
Load Transfers:
Rear
m
ur
hur
n 1) Due to the sprung mass
inertia force at the rear roll
t
r
center :
Ri Ro=Ri+DWf
 a *   1
DWr1   Ws a L  n 
 L   tr
 
2) Due to the restoring

moment by the rear r 1
DWr 2  M r
tr
 
suspension :

t
M sr
a*
L
W a
s L
Load Transfers:
W a Rear
ur L m
ur
hur
n 3) Due to the action of the
rear antiroll bar :
t
 
r
1 1
Ri Ro=Ri+DWf r
DWr3  Mra   Craf 
tr tr
4) Due to the inertia forces 1

DWr 4  Wur a L h ur 
of the rear unsprung mass : tr

Load Transfers: Rear
The sum of the four components gives the total

front load transfer.
4
DWr   DWri
1
The tilting moment of unsprung mass in case of

beam axle type of suspensions is taken as zero for
the calculation of the sprung mass roll angle. Note,
however, that the load transfer due to the inertia
forces of the unsprung masses, given by DWr4 is
not zero.

Dynamic Wheel Loads
Therefore, the inner and outer wheel loads

are given by :
Wf Wf
Wfi   DWf Wfo  + DWf
2 2
Wr Wr
Wri   DWr Wro  + DWr
2 2

Understeer ?
Dynamic Wheel Loads
 The vehicle will be in a critical state
when one the wheel loads becomes
zero.
 For an understeer vehicle the front
inner wheel is going to lift up first.
 On the other hand, the rear inner

wheel will lift up first in the case of
an oversteer vehicle.
 For the vehicle to start rolling over,

Oversteer ? both of the inner wheels must lift up.

Example
 Calculate, for the specified

minibus size armored vehicle
prototype, the
 roll angle and
 individual wheel loads
for steady cornering at a lateral
acceleration of 0.5 g.

2275
1420 2400 1800

4660 1890
Mass : 4040 kg
Unsprung mass (f/r) : 215/410 kg
Load distribution (f/r) : 0.43/0.57 %
Height of centre of gravity : 1200 mm
Tire radius : 366 mm
Suspension system : Beam axle - front and rear
Roll center height (f/r) : 350/450 mm
Suspension stiffness (f/r) : 180/275 kN/m
Distance between suspension springs (f/r) : 770/900 mm

Roll Angle at 0.5 g
 Calculate sprung mass and axle loads.
Ws  W  Wuf  Wur  (4040  215  410) 9.81  33501[N]
b
Wf  W  0.43 4040(9.81)  0.43 39632  17042 [N]
L
a
Wr  W  0.57 39632  22590 [N]
L

Roll Angle at 0.5 g
 Calculate sprung mass center of gravity position.
Wr  Wur 22590  4022

*
a  L 2.4  1.33 [m]
Ws 33501
Wf  Wuf 17042  2109
*
b  L 2.4  1.07 [m]
Ws 33501
hW  r Wuf + Wur  1.2 39632  0.366 2109 + 4022

hs    1.353 [m]
Ws 33501
a *n + b *m 1.33 0.45 + 1.07 0.35
h  hs 
*
 1.353   0.948 [m]
L 2.4

Roll Angle at 0.5 g
 Tilting moment of the sprung mass
proportional to lateral acceleration :
 K i = Wsh = 33501  0.948    31579   Nm

*
 Tilting moments of the front and rear

unsprung masses are zero since both front
and rear suspensions are of the beam axle
type.
M ft  M rt  0

Roll Angle at 0.5 g
 The roll stiffness of the front and rear

suspensions :
 t2   0.77 2 
C f   sf k  
f
 180000  53361 Nm / rad
 2   2 
   
 t2   0.9 2 
C r   sr k  
r
 275000  111375 Nm / rad
 2   2 
   

Roll Angle at 0.5 g
 Now the steady state roll angle can be calculated.
 K a 

f i
i L  

 31579  0.5 
 0.12rad  6.8o
*  
  C j  Wsh   (53361 + 111375)  33501(0.948) 
 j 
 This roll angle at a lateral acceleration of 0.5 g is
somewhat excessive, and some means of limiting
the roll angle without impairing the ride comfort
may be implemented.
 One possible solution would be the use of antiroll
bars.

Front Load Transfer at 0.5 g
 To calculate the individual wheel loads, the front and
rear load transfers must be determined first.
 b *   1  1.07   1
DWf 1    
Ws a L m     335010.50.35  1452 [N]

 L 
  t f  2.4   1.8
DWf 2   
Mrf
1 
 Wsh  aL cos f + sin f  
tf 
*   Cf
  Cf + Cr
1

 tf
 53361  1
  33501 0.948  0.5cos 0.119 + sin 0.119   
   1.8
 53361 + 111375 
 3516[N]

Front Load Transfer at 0.5 g
 Since there are no antiroll bars :
𝐫
𝐟𝟑 𝐟𝐚
𝐟
1 1
DWf 4  Wuf a L h uf   21090.50.366  214[N]
tf 1.8
 The front load transfer is obtained by summing the

individual load transfer components.
4
DWf   DWfi  1452 + 3516 + 0 + 214  5182 [N]
1

Rear Load Transfer at 0.5 g
 a *   1  1.33   1
DWr1   


Ws a L n

      
 33501 0.5 0.45   2942 [N]
 L   t r  2.4   1.42
 
DWr2  M rr
1 
 Wsh  aL cos f + sin f  
tr 
*   Cr
  Cf + Cr
1

 tr
 111375  1
   33501 0.948  0.5cos 0.119 + sin 0.119    
 53361 + 111375  1.42
 9302[N]

Rear Load Transfer at 0.5 g
 Since there are no antiroll bars :
 t1  0
DWr 3  M rra
r
1 1
DWr 4  Wur a L h ur   40220.50.366  518[N]
tr 1.42
 The rear load transfer is obtained by summing the

individual load transfer components.
4
DWr   DWri  2942 + 9302 + 0 + 518  12762 [N]
1

Wheel Loads at 0.5 g
 The individual wheel loads can now be found.
Wf 17042
Wfi   DWf   5182  3339 [N]
2 2
Wr 22590
Wri   DWr   12762  1467[N]
2 2
Wf 17042
Wfo  + DWf  + 5182  13703 [N]
2 2
Wr 22590
Wro  + DWr  + 12762  24057 [N]
2 2

Wheel Loads at 0.5 g
 It is obvious that the rear inner wheel of

this vehicle has already lifted up the
ground at a lateral acceleration below 0.5 g.
 There is still some safety margin, however,
before the front inner wheel also lifts up
resulting in vehicle turn-over.
 The fact that the rear inner wheel lifts up
first indicates an oversteering behavior at
high lateral accelerations.

ME 513
VEHICLE
DYNAMICS

Prof. Dr. Y. Samim
Ünlüsoy
CH 6 - HIGHER ORDER
HANDLING MODELS
 A bicycle model is adequate to for the investigation

of road vehicle handling dynamics in the low lateral
acceleration range, i.e., lower than about 0.3 g.
 However, in the case of higher lateral accelerations,
most of the assumptions made in developing the
bicycle model are not valid.
Thus, models with added degrees of freedom,
nonlinearities, and other characteristics become
necessary. The roll center/axis concept would be
particularly useful in the extension of the bicycle
model.

HIGHER ORDER HANDLING MODELS
 The first extension would be the addition of the roll
degree of freedom, resulting in a three dof model.
h*sin
Thus, the roll motion
of the sprung mass WsaL h*
around the roll axis W h*cos
 (roll velocity) 
will now be included. x
Fco Fci y
t/2 t/2
Wo z Wi

 Next, the load transfer should be taken into

account. The use of a four wheel model instead of
a bicycle model may thus be appropriate.
L
b a
Fx22
Fy22 Fy12 Fx12

Rg Frf
tr x
Frr m, Ra d
Izz, Ixx
Fx21
Fy11
Fy21 y Fx11

 Nonlinear tire behavior will be necessary

for cornering at high lateral
accelerations.
 Dugoff Model,
 Magic Formula
 Simplified Analytical Models

 Further, wheel rotational degrees of freedom may be required

to obtain tire slips and the forces generated with nonlinear tire
models.
 w = TT - rw Fx - rw R r I w ω
I wω  w = -Tb + rw Fx - rw R r
mw : Wheel mass,
TT Tb Iw : Wheel spin moment of
inertia
ww aw
TT : Tractive torque
mwg rw Tb : Braking torque
U
Fx : Tractive force
Rr Fx Rr : Rolling resistance
Fz rw : Tire radius

 Thus models with

- 3 dof (side slip, yaw, roll),
- 5 dof (side slip, yaw, roll, and wheel rotations
for 2 wheels),
- 7 dof (side slip, yaw, roll, and wheel rotations
for 4 wheels), or
- 8 dof (longitudinal, side slip, yaw, roll, and
wheel rotations for 4 wheels),
may be obtained and used according to the
requirements of specific problems.

3 dof Vehicle Handling
Model

3 dof Vehicle Handling Model
The basic bicycle h*sin
model can be extended

to a 3 dof handling WsaL h*
W h*cos
model by including the

roll degree of freedom
x
in addition to the side
slip and yaw velocity Fco Fci y
degree of freedoms. t/2 t/2
Wo z Wi

 In the 3 dof handling model the

vehicle mass is divided into the
sprung and unsprung masses.
β
Thus now the vector

x =r

of the state variables is given by:
p   p
where is the roll velocity.


 The velocity of the center of gravity of the

vehicle will be given in terms of the unit
vectors along the body fixed set of axes.
  
V = ui + vj
 These axes rotate with the yawing vehicle.
 
ω = rk

 The acceleration of the vehicle center of gravity in terms

of the unit vectors along the body fixed set of axes is
given by :

 dV       
a=
dt C
 + vj
+ ωxV = ui 
 + rk x ui + vj    

 dV  
a= =  u - vr  i +  v + ur  j
dt


dV  
=  u - vr  i +  v + ur  j 3 dof Vehicle
dt
Handling Model

 Thus the equations of motion in  dV
translation for the longitudinal and ΣF = m
lateral directions will be given by : dt
ΣFx = m  u - vr 
ΣFy = m  v + ur 

ISO z
The equations of motion in rotation will be
r 
given by :
p
 
x  dH   dω  
y ΣM = + ωxH = I + ωxH
dt dt
y where
p H : angular momentum of the vehicle,
I : inertia matrix, and
x r w : rotation vector.
z SAE

 The angular momentum of the vehicle with roll
and yaw rotational degrees of freedom is
expressed as :  
H = Iω
 Ixx -Ixz    
where I= and ω = pi + rk
 -I zx I zz 
Assuming symmetry with

respect to x-z plane  I xx 0 
I=
(i.e., Ixz = Izx = 0)  0 I zz 

 
H = Iω
 I xx 0 
I= 
 0 I zz 
    Carry out the intermediate calculations.
ω = pi + rk    
H = Iω = Ixx pi + I zzrk
     
dω dω
 + rk
= pi  I  + I zzrk
= I xx pi 
dt dt
     
  
ωxH = pi + rk x Ixx pi + I zzrk 

=  -I zz pr + Ixx pr  j

  
dω
I  + I zzrk
= Ixx pi 
dt

3 dof Vehicle
 
ωxH =  -I zz pr + Ixx pr  j Handling Model
 The equations of motion in rotation will be
given by :
 
 dH   dω  
ΣM = + ωxH = I + ωxH
dt dt
   
 +  -I zz pr + Ixx pr  j + I zzrk
ΣM = Ixx pi 
 Considering the longitudinal

ΣM x = Ixx p
(roll) and vertical (yaw) axes : ΣM z = I zzr

L
b a
3 dof
Fxr
u
Vehicle
Fxf d Handling
r
Fyr v F Model
yf
 Along the longitudinal axis x : m  u - vr  = ΣFx
 
m  u - vr  - ms h*rp = Fxf cosδ + Fxr - Fyf sinδ
 Here, the inertial term “msh*rp” arises due to the rotation of

the roll momentum vector by the yaw motion.

Acceleration
p+dp
dq p x
dp
dp = pdθ
 vr
u-
dp dθ
a'' =  h * =  h * p =  h * pr h*pr
dt dt
 Here, the longitudinal acceleration term h * rp arises due

to the rotation of the roll vector by the yaw motion.

L
3 dof
b a
Fxr Vehicle
u Handling
Fxf d
r Model
Fyr v F
yf
 Along the lateral axis y : m  v + ur  = ΣFy ms h*p
m  v + ur  + ms  h * p  = Fyf cosδ + Fyr + Fxf sinδ h*

p
 Here, the term msh * p stands for the inertial
force of the sprung mass due to roll
acceleration around the roll axis.
Roll axis
L
3 dof
b a
Fxr Vehicle
u
Fxf d
Handling
r
Fyr F Model
v yf
 Around the vertical axis z : I zzr = ΣM z
I zzr =  Fxf sinδ + Fyf cosδ  a - Fyr b

ms h*  v +ur 
3 dof Vehicle
msg
h*cos
Handling Model
h* 
Ixx p = ΣM x
 Around the longitudinal
(kf+kr)+(cf+cr)p axis x :
  
I xxp + msh*  v + ur  = ms gh*sin - k f + k r  - cf + cr p 
 Here, the term msh *  v + ur  stands for the moment (around
the roll axis) of the inertial force due to lateral acceleration of
the sprung mass . The terms on the right hand side are the
tilting moment of the sprung weight and the opposing
moments of the suspension spring and damper forces around
the roll axis. ki and ci represent vehicle roll stiffness and
damping.

 Thus the four equations of motion (longitudinal, lateral, yaw, and
roll) are :
m u - vr  - msh*rp = Fxf cosδ +Fxr - Fyf sinδ
m  v +ur  + msh * p = Fyf cosδ +Fyr  Fxf sinδ
Izzr = Fyf cosδ  Fxf sinδ  a - Fyrb
  
Ixxp + msh*  v + ur  = msgh* sin - k f +k r  - cf + cr p

 As before we may assume small front wheel steering angle

(i.e. sind = 0 and cosd = 1) and small roll angle (sin = ):
m u - vr  - msh*rp = Fxf +Fxr
m  v +ur  + msh * p = Fyf +Fyr
Izzr = Fyf a - Fyrb
  
Ixxp + msh*  v +ur  = msgh* - k f +k r   - cf + cr p
  

 The first equation is basicly related to longitudinal motion
(vr and rp are products of small quantities). Thus it can be
removed and three equations are left. Forward velocity u is
taken as a parameter (U) rather than a variable, leaving four
variables, v, r, p, and .
m  v +Ur  + msh * p = Fyf +Fyr
Izzr = aFyf - bFyr
  
Ixxp + msh*  v +Ur  = msgh* - k f +k r   - cf + cr p
  

 Due to the term involving , it is better to use an augmented
vector of variables;
 Introducing   = p
m  v +Ur  + msh * p = Fyf +Fyr
Izzr = aFyf - bFyr
  
Ixxp + msh*  v +Ur  = msgh* - k f +k r   - cf + cr p
  
 = p

 The augmented equations are given in matrix form below.
v 
r 
 M x  =  A ii x + Ff  1 1 
a -b 
 
x =    m 0 msh * 0  F =  
p  0 0 
 0 I zz 0 0  
   M =    0 0 
 ms h * 0 Ixx 0
  Fyf 
 0 0 0 1
f  =  
0 -mU 0 0  Fyr 
0 0 0 0 
 A ii  =   k roll = k f + k r
 0 -msh * U -croll msgh * -k roll 
  croll = c f + c r
0 0 1 0 
M x  =  Aii x + F f 

1 1 
a -b  Fyf 
F =   f  =  
0 0  Fyr 
 
 0 0 
Linear forms of the v + ar

tire cornering Fyf  Cf a f αf = -δ

U
forces are
Fyr  Cr a r v - br
introduced: αr =
U

3 dof Vehicle
M x  =  A i x + Bi  d
Handling Model
 Cf + Cr aCf - bCr 
 - mU 0 0 
U U
 
 aCf - bCr a 2Cf + b 2Cr 
 i 
A = 0 0 
U U
 
 0 -msh * U -croll msgh * -k roll 
 0 0 1 0 
 
 m 0 msh * 0  -Cf 
 0 I zz 0 0 -aC 
M =   Bi  =  f 
 ms h * 0 Ixx 0  0 
   0 
 0 0 0 1

 Now you can obtain [A] and {B} by
premultiplying [Ai] and {Bi} by the
inverse of [M].
 M x  =  A i x + Bi  d

-1 -1
x =  M   A i x +  M  Bi  d

x  =  A x + B d

Exercises
1) Obtain the steady state equations for the 3 dof
linear bicycle model.
2) Show that:
- The addition of the roll degree of freedom
does not alter the basic steady state
handling behavior of the vehicle, i.e. over-,
neutral, or under-steer characteristic as
well as the critical or characteristics speed.
- It will, however, introduce another
possibility for instability. Check the
characteristic equation.

ME 513
VEHICLE
DYNAMICS

Prof. Dr. Y. Samim
Ünlüsoy
CH 7 - RIDE COMFORT
Part 1

RIDE COMFORT
 Vehicle ride refers to the vibrations

of a vehicle which adversely affect
the passengers and/or the load.
 The vibrational activity of a vehicle
may be excited by a variety of
sources. Road irregularities,
however, is almost exclusively
considered to be the main excitation
in ride comfort studies.

Bounce
Yaw Yaw RIDE COMFORT
Sideslip  If considered as a rigid body, the
Pitch vehicle has six degrees of freedom,
Fore and three linear motions and three
Aft
Roll
rotations, in space.
 The relative emphasis on degrees of freedom in

subsections of vehicle dynamics :
Freedom Bounce Subject Ride
Sideslip Handling
Longitudinal (Fore & aft) Performance
Yaw Handling
Pitch Ride
Roll Handling, Ride
RIDE COMFORT
 In reality, a vehicle is a system of

bodies and thus can have a large
number of degrees of freedom.
 For this reason, analysis of vehicle ride
can be very complex and simplifying
assumptions must be made to be able
to obtain results.
 The frequency range of interest in most
ride comfort studies for motor vehicles
is usually taken as 0.5 to 25 Hz.

RIDE COMFORT
 The main disturbance input causing

vehicle vibrations comes from the
irregular road surfaces.
 Excitation may also arise from the
unbalance of the engine, torsional
fluctuations in the engine output,
unbalanced wheels, tires,
environmental inputs, etc.

RIDE COMFORT
 It is usually convenient to separate the ride

movements from the handling responses and
performance estimates.
 The justification for this separation is simply that
road irregularities mainly cause the car body to
pitch, bounce and roll. Although it can be seen from
the detailed suspension analyses that these body
displacements result in steering inputs, these are
usually second order effects such as roll steer.
 On the other hand, the steering inputs do not cause
large pitch and bounce motions. They do, however,
result in body roll.

Vibrational Characteristics of Road
Vehicles
 Major resonant frequencies of common road vehicles,
i.e. passenger vehicles, light and heavy commercial
vehicles :
 Frequency Range : 0.5 to 2.5 [Hz]
Body Bounce :
Motion of the sprung mass on
suspension system approximated
by bounce, pitch, and roll modes.

Vibrational Characteristics
of Road Vehicles
Frequency Range : 7 to 16 [Hz]

Wheel hop : Motion of the unsprung mass on
combined tire and spring
elasticity. Also includes "tramp"
of beam axles.
Wheel hop Tramp

Vibrational Characteristics of
Road Vehicles
Frequency Range :
7 to 18 [Hz]
Engine Vibrations :
Motion of the engine
mass on engine mounts.

Road Vehicles

Body structural vibrations
Body vibrations in beam and torsion modes.

Road Vehicles

Vibrational
Characteristics
of Road
Vehicles
 Frequency Range :
From 10 [Hz] upwards
Body cavity resonances

Road Vehicles
 Frequency Range : From 25 [Hz] upwards

 Tire flexural vibrations.

Vibrational
Characteristics of
Road Vehicles
 Frequency Range : From 40 [Hz]
upwards
 Higher flexural modes of
structure and panel vibrations.
Vibrations of a rectangular
membrane in some of its modes.
Vibrations of a circular membrane in some of its modes.

Road Vehicles
Forced vibrations excited by the IC engine:
- Frequency Range : 10 to 100 [Hz]
Source : Engine primary unbalanced
forces.

Source : Engine secondary unbalanced
forces, Hooke's joints.

Vibrational
Characteristics of
Road Vehicles
Drive line and engine bending resonances.

Vibrational
Characteristics
of Road
Vehicles
Road Noise.
Source : Road irregularities
High frequency vibrations
Source : Machinery, gears etc. .

BODY BOUNCE

Sprung
mass M Body Bounce
y
(0.5 –2.5 Hz)
ks c
Suspension
is the up and down (vertical)
motion of the sprung mass
Unsprung on the suspension spring.
mass
 The sprung mass is represented as an ideal

point mass and the suspension system
consists of a linear spring and a viscous
damper in parallel.
 The unsprung mass is assumed to be
stationary and tire elasticity is not considered.

Sprung
mass M Body Bounce
y
(0.5 –2.5 Hz)
ks
Suspension c
The bounce frequency of
the body is given by the
Unsprung
mass
expression :
 It has been known from the times of

horse carriage that the lower the ks
frequency of vibration of the car body bb 
on its suspension springs, the greater M
the comfort.
 However, the possibility of sea sickness below 0.5 Hz
puts an absolute lower limit on body bounce
frequency.
ks Example :
bb 
M
Typical mid-sized saloon car
M (laden) = 1280 kg, ks = 75 kN/m
N
75000  
loaded
bb   m   kg.m   7.65  rad 
1280 kg  N.s 2   s 
 
loaded  rad   cyc 
bb  7.65      1.2 [Hz ]
 s   2 rad 

ks
bb  Example
M
Typical mid-sized saloon car
M (unladen) = 1000 kg, ks = 75 kN/m
N
75000  
empty
bb   m   kg.m   8.66  rad 
   s 
1000  kg   N.s 
2  
empty  rad   cyc 
bb  8.66      1.4 [Hz]
 s   2 rad 

Suspension
bb
Damper
y
M
ks
c
c
bb 
2Mbb
 Response of a one-dof mass-spring damper system to

an initial displacement of the mass.
c
bb  Example
2Mbb
 Typical mid-sized saloon car :

M = 1280 [kg],
ks = 15 kN/m/per suspension,
c = 980 N/m/s/per suspension
 Ns 
4  980   
bb  m  0.20
 rad 
2  1280   kg   7.65  
 s 
Question:
Why is the damping ratio low compared with common values
used in other applications such as control and robotic systems?

Body Bounce (0.5 –2.5 Hz)
 Improved ride has always been a priority,

and the suspension design starts with the
selection of a suitable suspension
frequency.
 A usual order of priorities is to choose
the lowest possible body bounce
frequency for ride comfort first, and then
to control roll and tune handling by the
use of anti-roll bars.

 The lower limit for the body bounce
frequency, and hence the suspension spring
stiffness, is dictated by the need to limit the
static deflection (wheel travel) that should be
accommodated when the vehicle load varies
from minimum to maximum.
 The undamped natural frequency and static
deflection are related to each other.
Mg ks ks g g g
 bb    
2
ks M Mg  bb

g
δ= 2 Body Bounce (0.5 –2.5 Hz)
ωbb
1.5
1.45
In modern medium size
Body Bounce Frequency [Hz]
1.4
1.35 to large cars the body
1.3 bounce frequency has
1.25 generally settled at
1.2
about:
1.15
1.1
70 cyc/min [  1.2 Hz].
1.05
1
0.1 0.15 0.2 0.25
Static Deflection [m]

 Before the common use of anti-roll bars

to control roll, the suspension springs
used to be stiffer to be able to limit roll
angles during cornering.
 Thus, frequencies settled at around :
 90 c/min [1.5 Hz] for medium-size
to large cars, and
 120 c/min [2.0 Hz] for small cars.

 The introduction of anti-roll bars has removed
the need to supply this stiffness by the main
springs .
 For this reason :
 small cars seem to become limited at
around 90 c/min [1.5 Hz] unladen.
 In modern medium size to large cars the
body bounce frequency has generally
settled at about 70 c/min [1.2 Hz].
 The off road vehicles conventionally have
frequencies around 120 c/min [2.0 Hz]
unladen.

WHEEL HOP

Sprung
mass
Suspension ks c
Wheel Hop
y
Unsprung m (Wheel Bounce)
mass
(7 –16 Hz)
Tire kt
 In this mode the unsprung mass (wheel-hub assembly)

bounces on the tire and the motion is damped by the
shock absorber.
 The wheel hop frequency (around 10 Hz) is usually well
separated from body bounce, pitch, and roll frequencies
which are considerably lower (just above 1 Hz).
Sprung
mass
Wheel Hop
Suspension ks c
(Wheel Bounce)
Unsprung m
y (7 –16 Hz)
mass
ks  kt
Tire kt  wh 
m
Note that for wehicles with independent suspensions :
M kt  wh ks  kt M
  10   10
m ks bb m ks

Typical mid-sized saloon car – single suspension :
M = 1280/4 = 320 kg,
m = 46 kg,
ks = 15 kN/m, Example
kt = 145 kN/m
N
15000  145000  
ks  kt  m   59  rad 
 wh  
m 46 kg  s 
 
59
 wh   9.4 Hz 
2

QUARTER CAR MODEL

Quarter Car Model
Combined Body Bounce
& Wheel Hop
y
Sprung 2  In many ride studies, a
M
mass model involving both the
body bounce and wheel
k
hop modes is required.
Suspension c
s The so-called 'Quarter car
model' is the best known
y vehicle ride model with
Unsprung 1
m two degrees of freedom.
mass
 It usually represents one
Tire k
t suspension together with
y
0 its share of the sprung
mass.

y
Sprung 2
mass
M Quarter Car
Model
Suspension k c
s
ks
y bb 
Unsprung
m
1 M
mass
ks  kt
Tire k
t
y
 wh 
0 m
 It is possible to find the approximate body bounce

and wheel hop frequencies using the previous
approximate expressions.
 As these frequencies are well separated, the results
are normally quite close to exact values.

Amplitude Ratio
Amplitude Ratio
Single dof models Quarter car model
 It is noted that at 1 Hz the response amplitude of the

wheel hop is negligible. Similarly at about 9 Hz the
response amplitude of the sprung mass is negligible.
 Thus the body bounce and wheel hop vibrations are
effectively uncoupled.
y
Sprung 2
mass M
Damper
Suspension k c
s
Unsprung
y
1
Note that common use
m
mass of “shock absorber” is
Tire k
incorrect.
t
y It is the spring that
0
absorbs the shocks.
 The conventional damper has to damp the

sprung mass on its springs around 1 Hz, and the
unsprung mass on the combined elasticity of
suspension spring and tires at around 10 Hz.

c c
 bb   wh  Damper
2Mbb 2m wh
 However, with independent suspensions, the ratio

(sprung mass / unsprung mass) is similar to the ratio
(wheel bounce frequency / suspension frequency) so that
the percentage of critical damping is similar for both
modes with which the damper is trying to handle.
 Thus, the suspension damper can damp body bounce
and wheel hop vibrations almost equally well.
 bb 2Mbb  M   bb 
       1.0
 wh 2m wh  m    wh 
 bb   wh

c c
 bb   wh  Example
2Mbb 2m wh
 Typical mid-sized saloon car – single suspension :

M = 1280/4 = 320 kg,
m = 46 kg,  Ns 
ks = 15 kN/m, 980  
 bb  m  0.22
c = 980 N/m/s
 rad 
kt = 145 kN/m 2320kg6.9  
 s 
Note that the damping ratio  Ns 
low compared with 980  
 m
common values used in  wh   0.18
other applications such as  rad 
control and robotic systems. 246kg59  
 s 
y
Sprung 2
mass M
Quarter Car Model
k c
Suspension s The model can be used to
y
study the motion of a
Unsprung
mass m
1 vehicle on an irregular
road surface.
k
Tire t
y
0 𝐌𝐲̈ 𝟐
𝐲𝟐
M
Equations of motion for this system are :
𝐤 𝐬 𝐲𝟏 − 𝐲𝟐 𝐜 𝐲̇ 𝟏 − 𝐲̇ 𝟐
𝒌𝒔 𝒚𝟏 − 𝒚𝟐 𝐜 𝐲̇ 𝟏 − 𝐲̇ 𝟐
𝐲𝟏
m
𝐦𝐲̈ 𝟏
𝐤 𝐭 𝐲𝟎 − 𝐲𝟏

2  cy 2  k s y 2  cy 1  k s y1
My
1  cy 1  (k s  k t )y1  cy 2  k s y 2  k t y 0
Quarter Car
my
Model
Equations of motion in matrix form :
 M y + Cy  +  K y = f 

where
y1  k t 
y    f     y 0
y 2  0
m 0   c  c k s  k t  ks 
M    C   K   
 0 M 
   c c 
   ks k s 

y2
M State Equation Form
 The equations of motion for
ks c the quarter car model can be
converted into state space
y1 form by choosing a state
m vector consisting of the
 sprung mass velocity,
 suspension deflection,
kt  tire deflection, and
y0
 unsprung mass velocity
(arbitrary).

y2
M State Equation Form
 Choose* the state variable
ks c vector as :
Sprung mass
y1  y 2  velocity
m y - y 
 1 2 Suspension
x =   travel
 y 0 - y1  Tire
kt  y 1  deflection
y0
Arbitrary
* Note that other choices may also be

made.

M y + Cy  + K y = f 
STATE EQUATION
To convert the conventional matrix equations of
motion to state space form numerically:
a) Introduce a state vector  x1   y 
x =  x  =   
 2  y 
 x 1   y 2 
b) Use the fact that  x  =   =   to show that
 x 2   
y
the state space form is then given by:
  0 I    0 
x  =  A x   Bu =    x +  -1 u

 -  M   K  -  M   C  M  
-1 -1
Note: This numerical approach is not needed for the quarter car model as
the conversion to state form is easy and straight forward.
 y 2  State Equation
y - y 
x =  1 2  Form
 y 0 - y1 
 y 1   Equations of motion in state form :
x  =  A x + W y 0
 c ks c Note that this is

- M M
0
M  now a disturbance
  0 
 -1 0 0 1  input. Further it is
 A =  
0 
  the velocity and not

0 0 0 -1
    
W = the displacement of
 c ks kt c 1  the road surface
 m - -  0 
m m m profile.

Handling versus Ride Comfort
The sprung mass acceleration here, should not be
confused with lateral acceleration during cornering. In
this case, the passengers on vehicles may be
exposed to accelerations as high as 3 to 4g at high
speed on very rough (off-road) surfaces.
However, normal humans can withstand no more than
9g's, and even that for only a few seconds. When
undergoing an acceleration of 9g’s, the body feels
nine times heavier, blood rushes to the feet, and the
heart can't pump hard enough to bring this heavier
blood to the brain. The vision narrows to a tunnel,
then goes black. If the acceleration doesn't decrease,
people will pass out and the final result is fatal.

M y + Cy  + K y = f  Quarter Car Model
To determine the undamped natural frequencies of

the system, consider the equations for undamped free
motion.
M y + K y = 0
If a solution of the form y = Y e iωt
is assumed, then :
   M  2   K   Y  0
 

   M  2   K   Y  0
  Quarter Car Model
For a nontrivial solution to this equation, one must

have :
 K  -  M  ω2 = 0
k s  k t  m  ks
0 2
 ks k s  M
where λ=ω
The characteristic equation is then obtained as :
mMλ 2 - k s m +  k s + k t  M  λ + k s k t = 0

mMλ 2 - k s m +  k s + k t  M  λ + k s k t = 0 Quarter Car Model
Dividing all terms by mM and using the previous
definitions :
2   ks + k t  ks  k sk t
λ - + λ + =0
 m M k sm
2
λ -  2
ωbb 2
+ ω wh  2  2
λ + ωbb  ω wh

-
ks
m

=0

The exact undamped natural frequencies for the
quarter car are obtained as :
 2 
2 1  k s + k t k s   ks + k t ks  4k sk t 
ω1,2 =  + ±  +  
2  m M  m M mM 
 

 2 
1  k + kt ks   k s + k t k s  4k sk t 
2
=  s
ω1,2
2 

+ ±
m M

 m
+  -
M mM 

Example
 kN   kN 
M  320  kg  , m  46  kg  , k s  15  ,
 t k  145 m
 m   
The exact undamped natural frequencies are

obtained as :
 
 15000 + 145000 15000  4 15000 145000  
2
2 1  15000 + 145000 15000 
ω1,2 =  + ±  +  -
2  46 320   46 320   46  320  
 
𝐤𝐬 𝟏𝟓𝟎𝟎𝟎 𝐫𝐚𝐝
= = 1.09 𝐇𝐳  rad 
𝛚𝐛𝐛 =
𝐌 𝟑𝟐𝟎
= 𝟔. 𝟖𝟓
𝐬 ω1 = 6.51   = 1.04 Hz 
 s 
 rad 
𝐤𝐬 𝐤𝐭 𝟏𝟓𝟎𝟎𝟎 𝟏𝟒𝟓𝟎𝟎𝟎 𝐫𝐚𝐝 ω2 = 59.01  = 9.4  Hz 
 s 
𝛚𝐰𝐡 = = = 𝟓𝟖. 𝟗𝟖 = 9.4 𝐇𝐳
𝐦 𝟒𝟔 𝐬

2
λ -  2
ωbb + ω2wh  2  2
λ + ωbb  ω wh -

ks 

m
=0 Quarter Car Model
It is clear that if k s << mω2wh or 2

ωbb << ω2wh
2

λ2 - ωbb + ω2wh λ + ωbb
2

ω2wh = λ - ωbb
2
  
λ - ω2wh = 0
and the roots of the characteristic equation will be

the same as the approximate body bounce and
wheel hop frequencies.
Thus, the body bounce and wheel hop motions are
decoupled from each other, and the approximate
equations for body bounce and wheel hop
frequencies may give reasonably accurate results.
Prof. Dr. Y. Samim Ünlüsoy ME 436 Automotive Engineering II 53

Example
M= m=
𝐍 𝐍
𝐬 𝐦 𝐭 𝐦
𝐚 𝐚
𝐛𝐛 𝐰𝐡 – approximate
𝐞 𝐞
𝐛𝐛 𝐰𝐡 – exact
Exact and approximate undamped natural frequencies

are close to each other since
𝟐
𝟐 𝟐
𝐰𝐡

TRANSFER FUNCTIONS

M y + Cy  + K y = f  y o Transfer
Functions
Taking Laplace Transform of each term :
  Y  s  + C s Y  s  + K Y  s  = f  Yo  s 
M s 2
The transfer function matrix is obtained.

 Y1 (s) 
 H  s   =
 Y  s 
=
f Y(s) = 
Y
 2 (s)

Yo  s  M  s2 + C s + K  Y (s)

H1 (s) = 1
Y0 (s)
Note that since there are two outputs, the TF Y2 (s)
H 2 (s) =
matrix will have two transfer functions, H1 and H2. Y0 (s)

Transfer Functions
 Transfer functions for sprung mass acceleration,
suspension travel, and tire deflection for road profile
velocity input can be established as
 Sprung mass acceleration Suspension travel (or deflection)
s2Y2 (s) sk t A Y (s)  Y1 (s)

Hst (s)  2
k (A  B)
 t
Hsma (s) = =
sY0 (s) BC - A 2 sY0 (s)

s BC  A 2 
Tire deflection
where
Y1 (s)  Y0 (s) (k t -C)B+A 2
s = j H td (s)  
sY0 (s)

s BC  A 2 
A  cs s  k s B  Ms2  css  k s C  ms2  css  k s  k t
2
s Y2 (s) sk t A Transfer
Hsma (s) = =
sY0 (s) BC - A 2 Functions
30
 Body acceleration
20
measure of ride
10 comfort
Typical mid-sized
| Hsma| (dB)
0 
saloon car
-10
M = 1280/4 = 320 kg,

-20
m = 46 kg,
-30 ks = 15 kN/m,
-40
kt = 145 kN/m,
10 -2 10 -1 10 0 10 1 10 2 c = 980 N/m/s
Frequency (Hz)

Transfer
Functions
 Body
acceleration
 Related to
ride comfort
Prof. Dr. Y. Samim Ünlüsoy Automotive Engineering 59

Transfer
Functions
 Body
acceleration
 Related to
ride comfort

Transfer
Functions
 Body
acceleration
 Related to
ride comfort

Transfer
Functions
 Body
acceleration
 Related to
ride comfort

Transfer
Functions
 Body
acceleration
 Related to
ride comfort

Y2 (s)  Y1 (s) k t (A  B)
Hst (s)   Transfer
sY0 (s) s BC  A 2   Functions
0
 Suspension Travel
-10
(Rattle space)
-20
important for
-30
packaging
-40
|Hst| (dB)
-50  Typical mid-sized

-60 saloon car
-70
M = 1280/4 = 320 kg,
-80
m = 46 kg,
-90
ks = 15 kN/m,
-100
10 -2 10 -1 10 0 10 1 10 2 kt = 145 kN/m,
Frequency (Hz)
c = 980 N/m/s
Y1 (s)  Y0 (s) (k t -C)B+A 2 Transfer
H td (s)  
sY0 (s) s BC  A 2   Functions
-20
Tire Deflection
-30
measure of road
holding
-40
 Typical mid-sized
|Htd| (dB)
-50 saloon car

-60 M = 1280/4 = 320 kg,
m = 46 kg,
-70
ks = 15 kN/m,
-80 -2 -1 0 1 2
kt = 145 kN/m,
10 10 10 10 10
Frequency (Hz) c = 980 N/m/s

Transfer Function in Matlab
 Obtain the expressions for the [A], [W], [C], and [D]
matrices of the state equations for the vehicle model.
 Specify the frequency range. For example, logarithmically
distributed 250 data points in the frequency range from
0.1 to 100 [Hz]:
f = logspace( -1, 2, 250); w = 2*pi*f;
 Obtain transfer function magnitudes and phases between
the disturbance input and the state variables:
[ T, p ] = bode( A, W, C, D, 1, w )
You can also make use of the definitions:
j=sqrt(-1), and s=j*w

INVARIANT POINTS OF
TRANSFER FUNCTIONS
Reading assignment:
- Hedrick, J.K., Butsuen, T., (1990) ‘Invariant properties of automotive
suspensions’, Proc. Instn. Mech. Engrs., v. 204, pp. 21-27.
- Karnopp, D., (2009) ‘How significant are transfer function relations and
invariant points for a quarter car suspension model?’, Vehicle System
Dynamics, v. 47, n. 4, pp.457-464.

Invariant Points
 Add the two equations of motion:
 2 + c  y 2 - y 1  + k s  y 2 - y1  = 0
My
1 + c  y 1 - y 2  + k s  y1 - y 2  + k t  y1 - y 0  = 0
my
1 + k t  y1 - y 0  = 0
 2 + my
My
 The resulting equation is independent of
suspension forces (i.e., does not involve ks or c)!
 This equation called the basic invariant equation
for automotive suspensions.

My 1 + k t  y1 - y 0  = 0
 2 + my Invariant Points
 Take the Laplace Transform of each term:
Ms2 Y2 (s) + ms2Y1 (s) + k t  Y1 (s) - Y0 (s) = 0
 Reorder:
 
Ms2Y2 (s) + ms2 + k t Y1 (s) = k t Y0 (s)

 
Ms2Y2 (s) + ms2 + k t Y1 (s) = k t Y0 (s) Invariant
Points
 Express in terms of transfer functions
s2Y2 (s) H td (s) 

Y1 (s)  Y0 (s)
Hsma (s) =
sY0 (s) sY0 (s)

MHsma  s  + ms2 + k t H td  s  = -ms 
 and insert s=j:

MHsma  jω  + k t - mω2 H td  jω  = -jmω 
 
MHsma  jω  + k t - mω2 H td  jω  = -jmω Invariant
Points
 It is observed that when
kt
k t - mω12 = 0 ω1 =
m
 Sprung mass acceleration transfer function Hsma
will have an invariant point
mk t
Hsma  ω1  = -j
M
 This invariant point depends only on specified
vehicle parameters and cannot be changed by the
application of passive or active suspension forces.

kt mk t Invariant
ω1 = Hsma  ω1  = j
m M Points
SPRUNG MASS ACCELERATION
30
Invariant
20 point
10
Magnitude (dB)
-10
Passive
-20 Full State Feedback
Velocity Feedback
-30
-40
-2 -1 0 1 2
10 10 10 10 10
Frequency (Hz)

 It can be
shown that
suspension
travel
(deflection)
function
has an
invariant
point at :
kt
ω2 =
M+m
M+m M+m
Hst  ω2  = j Invariant Points
M kt

Invariant Points
 The tire TIRE DEFLECTION
-20
deflection
transfer function,
H td  ω = 0  = 0
-30
on the other
Magnitude (dB)
hand, does not -40
have any
-50
invariant point
except at  = 0, -60
Passive
Full State Feedback
Velocity Feedback
i.e.
-70
H td  ω = 0  = 0 -80 -2 -1 0 1 2
10 10 10 10 10
Frequency (Hz)

BODY BOUNCE AND PITCH
(HALF CAR) MODEL

Body Bounce & Pitch
With combined bounce
and pitch freedom, it is
apparent that these two
modes are close, but the
L wheel hop frequencies
y f are well separated.
a b
So it is reasonable to
consider pitch and
m, Ip
ksf cf cr ksr bounce separately from
qf qr
unsprung mass
movements.
 A two degree of freedom model is thus necessary to
examine the pitch motion of the body in addition to the
body bounce.
L
y f
a b
Body Bounce
m, Ip
& Pitch
ksf cf cr ksr
qf qr
 The vehicle has a mass m, and ksf and ksr are the effective
vertical stiffnesses of the front and rear suspensions.
Denoting by cf and cr the effective damping constants of
the front and rear dampers, the model represents the
bounce, y, and pitch, f, modes of the vehicle. Let Ip be the
moment of inertia of the vehicle in pitch mode. The road
displacement inputs are shown as qf and qr.

L
f y
a b
Body
Bounce &
m, Ip
ksf cf cr ksr Pitch
qf qr
 Equations of motion in matrix form :

 
y  y  y  q f  q f 
M  φ  + C φ  +  K  φ  = F1  q  + F2  q 
       r  r
 Note that for undeformable road surfaces qf
and qr are identical, but qr is applied with a
delay after qf.
y y  y  q f  q f 
M   C    K    F1    F2   Body
f  f  f  q r  q r 
Bounce &
with Pitch
m 0   cf + cr ac f - bcr 
 M  =  0 I   C =  2 2 
 p acf - bcr a cf + b cr 
 k sf + k sr ak sf - bk sr   k1 k 3 
K  =  2 2 = 
ak sf - bk sr a k sf + b k sr   3k k 2
 cf cr   k sf k sr 
 F1  =    F2  =  
ac
 f -bc r ak
 sf -bk sr 

 -  M  ω 2 +  K   Y = 0
  Body Bounce & Pitch
 The characteristic equation is given by :
k sf + k sr - mλ ak sf - bk sr
2 2 =0
ak sf - bk sr a k sf - b k sr - Ip λ
 
mI p λ 2 -  Ip  ak sf - bk sr  + m a2k sf - b 2k sr  λ + L2k sf k sr = 0
 
where λ = ω2

mI p λ 2 -  mk 2 + Ip k1  λ + L2k sf k sr = 0
Body Bounce &
Pitch
 The undamped natural frequencies are given by :

 2 
2 1  mk 2  I p k 1  mk 2  I p k 1  k k  k 2

1,2  1,2      4 1 2 3

2 mIp  mI  mIp
 p  
 
where
k1  k sf  k sr
k 2  a2k sf  b 2k sr
k 3  ak sf  bk sr

k1  k sf  k sr  2 
1 mk + I k  mk + I k  2
k k -k 
k 2  a2k sf  b 2k sr
2
ω1,2 = 
2
2
mIp
p 1
 


2
mI p
p 1
 -4 1 2 3 

 mI p Example
 
k 3  ak sf  bk sr
 m = 1240 kg,
Ip = 1820 kg.m2 k1  38000  42000  80000[N / m]
2
a = 1.37 m, b = 1.14 m k 2   1.37   38000    1.14  2  42000   125905.4[Nm]
ksf = 38 kN/m, ksr = 42 kN/m k 3   1.37  38000    1.14  42000   4180[N]
 2 2
2
1,2
1
 
 1240  125905.4    1820  80000     1240  125905.4    1820  80000    80000  125905.4    4180  
  4 
2  1240  1820    1240 1820    1240   1820  
 
 Note that at this stage, the 1  63.22 [rad / s]  8.0[rad / s]  1.27 [Hz]
relation of these two frequencies 2  70.48 [rad / s]  8.4[rad / s]  1.34 [Hz]
with the body bounce and pitch
modes is not known!

 
mI p λ 2 -  Ip  k sf + k sr  + m a 2k sf + b 2k sr  λ + L2k sf k sr = 0
 
Uncoupling Bounce & Pitch

 The characteristic equation can be written in the form :
2 2  ak sf - bk sr 
2 2
 2

λ - ω y + ω φ λ + ω y ωφ -
mIp
=0
where y and f are the pure* bounce and pure* pitch

mode frequencies, given below.
k sf  k sr a 2 k sf  b 2 k sr
y  f 
m Ip
* Uncoupled

2 2  ak sf - bk sr 
2 2
 2

λ - ω y + ωf λ + ω y ω f -
mIp
=0 Uncoupling
Bounce & Pitch
λ 2 -  ω 2y + ω f2  λ + ω y2 ωf2 = 0
λ - ω λ - ω  = 0
y f
 It is noted that when ak sf  bk sr then the last term in

the characteristic equation is zero and thus the bounce
and pitch modes of the body are uncoupled (a moment
applied at the center of gravity will result in a rotation
about the center of gravity with no bounce and a force
applied at the center of gravity in the vertical direction
will result in only bounce with no rotation).
 The natural frequencies are then given by y and f for
the pure bounce and pitch modes.

L
a b
Pure
D m Bounce
k sf k sr Frequency
 In the pure (uncoupled) bounce

mode, the undamped natural
m frequency is given by:
y
k sf + k sr
ksf + ksr ωy =
m

L
a b mr
mf f
D2 Pure Pitch
D1
Ip
k sr
Frequency
k sf
 Moment balance M = aFf + bFr = a  Δ1  k sf + b  Δ 2  k sr

about the cog:
= a  af  k sf + b  bf  k sr =  a k sf + b 2k sr  f
2
f +  a 2k sf + b 2k sr  f = 0
I p  a 2
k + b 2
k sr

f + ωff = 0  ω f 
2 sf

f+
 sf k sr 
a 2
k + b 2
f=0
Ip
Ip Uncoupled (pure) pitch frequency
 M = 1240 kg, ak sf = 1.37  35000  = 47950  Nm
 Ip = 1820 kg.m2 bk sr = 1.14  42000  = 47880  Nm
 a = 1.37 m, b = 1.14 m ak sf  bk sr
 ksf = 35 kN/m, ksr = 42 kN/m
Example
N
35000  42000  
k sf  k sr  m   7.9  rad   1.25 Hz 
y  
m 1240 kg  s 
 
N
a 2k sf + b 2k sr
1.37 2
 35000  + 1.14 2
 42000     kg.m 
ωf = = m 
 2 
Ip 
1820 kg.m 2  Ns 
 
 rad 
ωf  8.1    1.3  Hz 
 s 
L
a b
Uncoupled
mr
Front & Rear
mf
m
Suspension
k sf k sr
Lateral axis Frequencies
Front and Rear Suspension Frequencies

 Let us now replace the mass of the vehicle by two
unequal parts
b a
mf  m mr  m
m  mf  mr L L
the former concentrated at the front axle vertical plane
and latter concentrated at the rear axle vertical plane.
L
a b Front & Rear
mf mr Suspension
m
k sf k sr Frequencies
Lateral axis
 Then one degree of freedom models at the

front and rear suspensions will be obtained.
The front and rear suspension bounce
frequencies in such a case will be given by :
k sf k sr
 yf   yr 
mf mr

L
a b Equal Pitch &
mf mr Bounce
k sf
m k sr Frequencies
Lateral axis
 When aksf = bksr

a(ak sf )  b(bk sr ) (a  b)ak sf Lk sf k sf
f      yf
mab mab bm mf
(a  b)bk sr Lk sr k sr
    yr
mab am mr
Thus pitch, front, and rear bounce

 f  f  r

frequencies become equal.

L
a b
mf mr Body Bounce
k sf m k sr & Pitch
Lateral axis
 The pitch moment of inertia of the

vehicle given by the expression : I p  a 2m f  b 2m r
b a
mf  m mr  m Ip  m a b
L L
 It is usual to write the moment of
inertia in the form (where k is the I p  mk 2
radius of gyration) :
Thus if aksf = bksr for ω f = ωf = ωr , 2

k ab
one has:
2
Ip Body Bounce & Pitch
k =
m
 The ratio k2/ab is an indication of how far the vehicle
departs from equal pitch and bounce frequencies.
k2/ab < 1  pitch frequency > bounce frequency
k2/ab = 1  pitch frequency = bounce frequency
k2/ab > 1  pitch frequency < bounce frequency
 A typical value of k2/ab is 0.9.
 The smaller k2/ab is, the smaller is Ip, and thus the
higher will be the pitch frequency..
 Pitching introduces fore-and-aft motion to which
people are more sensitive. So, If the almost pure
pitch mode has a higher frequency than the almost
pure bounce mode, this is disturbing to passengers.
k2 k2
< 1.0
ab
 1.0 Body
ab
Bounce &
ω f > ωbb ω f  ωbb Pitch

Body
Bounce
& Pitch
 To avoid pitching as the car goes over a bump on

the road
r  f
Then after a few cycles the two suspensions will
move in phase resulting in bounce of the body only.
Example
 Consider the vehicle with specifications listed below.
a) Can you tell which of the body bounce and pitch
frequencies will be higher ?
b) Is the use of the uncoupled (pure) frequencies
justified ?
c) Determine the pure bounce and pitch frequencies.
M = 1240 kg, Ip = 1820 kg.m2,

a = 1.37 m, b = 1.14 m,
ksf = 35 kN/m, ksr = 42 kN/m.

M = 1240 kg, Ip = 1820 kg.m2,
a = 1.37 m, b = 1.14 m, Example
a) Check the value of k2/ab.
Ip 1820
2
k 1240
= M = = 0.94
ab ab 1.37 1.14 
Since k2/ab<1, the pitch frequency will be

higher than the body bounce frequency.

M = 1240 kg, Ip = 1820 kg.m2,
a = 1.37 m, b = 1.14 m, Example
b) Check if aksf = bksr.
ak sf  1.37  35000   47950 N

bk sr  1.14  42000   47880 N
ak sf  bk sr
Since aksf  bksr, the pitch and bounce

frequencies will be very close to the exact
frequencies (1.25 and 1.29 Hz, respectively).

M = 1240 kg, Ip = 1820 kg.m2,
a = 1.37 m, b = 1.14 m, Example
b) Pure (uncoupled) bounce and pitch frequencies:
35000  42000   N 
k sf  k sr  m   7.9  rad   1.25 Hz
ωy    s 
m 1240 kg  
1.372 35000   1.142 42000   N 

a2k sf  b2k sr  m   kg.m 
 
ωφ 
Ip

1820 kg.m 2   

2 
 N.s 
 rad 
ωφ  8.1 
s   1.3 Hz
 

Body Bounce and Pitch
Oscillation Centers

k1  k sf  k sr
2 2
Body Bounce & Pitch
k 2  a k sf  b k sr
k 3  ak sf  bk sr Oscillation Centers
m 0   k sf + k sr ak sf - bk sr   k1 k3     M  2   K   Y  0
  0 I 
M = K  =  2 2 =
k 2   
 p ak sf - bk sr a k sf + b k sr  k 3
 Rewrite the eigenvalue  k  m2  Y

k3
equation in terms of k1,  1    0
k2, and k3:  k3
 k 2  I pm2   q 
to get two simultaneous

algebraic equations with  
k1  m2 Y  k 3q  0
the unknowns Y and q,
and the parameter 2. 
k 3 Y  k 2  I p 2 q  0 

 
k 1  m 2 Y  k 3 q  0 Body Bounce & Pitch

k 3 Y  k 2  Ip 2 q  0  Oscillation Centers
 Solve for the amplitude ratios r1

and r2 for =1 and then for =2
from each of the two equations :
 Remember
that Y and q r1 
Y 
1

k 3


 k 2  Ip12 
are the q 
1
 k1  m12  k3
amplitudes of
body bounce
r2 
Y 
2

k 3


 k 2  Ip22 
and pitch
motions.
q  2
 k1  m22  k3

Y  Y  Body Bounce & Pitch
1 2
r1  r2 
q  q 
1 2
Oscillation Centers
 The natural modes of vibration are then given by :

 Here m (1) and m (2) are called
the modal vectors.  Y 1 
(1)   (1) r1 
m  q  
 Natural vibration modes of  q 1  1
the system, Y (1) and q (1) and
 Y 2  
Y (2) and q (2), determine the (2)   (2)  r2 
m  q  
shape that the system  q 2   1
executes synchoronous
harmonic motion with frequencies 1 and 2.

Body Bounce & Pitch - Oscillation Centers
 There are two nodes of vibration. These nodes are
called the oscillation centers. With Y positive upwards
and q positive clockwise:
Y c.g.
 First one with positive Y
q Node
and negative q , and l1
q c.g.
Node
Y  the second with
l2 Y and q of the
same sign.

The distance, l1, of the Y c.g.

oscillation center from the
q Node
center of gravity will be l1
smaller than the wheelbase.
c.g.
The distance, l2, of the
q
Y
Node oscillation center
l2 from the center of
gravity will be larger
than the wheelbase.

 When a node is within the wheelbase, then
the ride is called pitchy.
Node
Y c.g. q
l1
Pitchy ride

 When a node goes outside the wheelbase, then the

ride is called bouncy.
l2
q c.g.
Y
Node
Bouncy ride

 As the two modes become uncoupled, i.e.,
aksf = bksr, one of the nodes goes to infinity, while
the other coincides with the center of gravity.
Thus the distances from the center of gravity
become :
l1  0 l2  ∞
l2
Y c.g.
q c.g.
q Y
Node l1 Node
Pitchy ride Bouncy ride

 On the other hand, if
Ip=mab
then the two nodes will be on the suspension to
body connection points.
l1  a l2  b
Y l2
c.g.
q c.g.
q Y
Node l1 Node

k 1  80000[N / m]
 m = 1240 kg,
Ip = 1820 kg.m2
k 2  125905[Nm] Example
k 3  4180[N]
a = 1.37 m, b = 1.14 m
ksf = 38 kN/m, ksr = 42 kN/m 1  63.22 [rad / s]  8.0[rad / s]  1.27 [Hz]
2  70.48 [rad / s]  8.4[rad / s]  1.34 [Hz]
 The natural modes of vibration are then

given by :
Y  1 Y  1
k 3 4180
r1   Bouncy Mode
q 
1
 k1  m12  r1 
q 
1

 80000  1240  63.22  
 2.601[m / rad]
Y  Y 
2 2
k 3 4180
r2   r2    0.565[m / rad]
q 
2
 k1  m22  q 
2  80000  1240  70.48   Pitchy Mode
r1 
Y 
1


 k 2  Ip12  r1 
Y  1

 125905  1820  63.22  
 2.601[m / rad]
Same results
are obtained if
q  q 
1 k3 1 4180
the other two
r2 
Y 
2


 k 2  I p22  r2 
Y 
2

 125905  1820  70.48  
 0.565[m / rad] equations
q  q 
2 k3 2 4180
are used !

L
y f
a b
Body
Bounce &
m, Ip
ksf cf cr ksr Pitch with
yf yr
muf mur Unsprung
ktf qf qr ktr Masses
 Now the model has four degrees of freedom;

namely the bounce and pitch of the body and
one bounce freedom each for the front and
rear unsprung masses.

L
y f
a b
Body
Bounce &
ms, Ip
ksf cf cr ksr Pitch with
yf yr
muf mur Unsprung
ktf qf qr ktr Masses
T
x = y φ yf y r 
 M x +  Cx  +  K  x =  F 
 ms 0 0 0 
 0 I
M =  p 0 0  Body Bounce & Pitch
0 

0 0 muf
 with Unsprung Masses
 0 0 0 mur 
 c f + cr ac f - bcr -c f -cr   0 0 
   0
a2 c f + b2 cr 0  qf 
 C =  acf - bcr -ac f bcr 
F = Hq = k  
-c f -ac f cf 0  tf 0   qr 
   
 -cr bcr 0 cr   0 k tr 
 k sf + k sr ak sf - bk sr -k sf -k sr 
 
 ak sf - bk sr a2k sf + b2k sr -ak sf bk sr 
  
K =

-k sf -ak sf k sf + k tf 0
 
 -k sr bk sr 0 k sr + k tr 

BODY BOUNCE and ROLL

t
y
t/2 g t/2
Body
Bounce &
ms, Ir
kl cl cr kr Roll
ql qr
Left Right
 Note that by a suitable interpretation of the half

car model parameters, a model for the body
bounce and roll motion can be obtained.
 Question: What can you say about the coupling
between body bounce and roll motions?

Engine Vibrations
The main function of engine mounts is
to carry the engine or the
engine+transmission. However,
ENGINE
- Forces and moments generated by
the engine are transmitted to the
vehicle through engine TAKOZ
- Vibrations of the engine and VEHICLE

vehicle due to road inputs interact .
They also locate the engine to ensure
that its motion is restricted so that no
interference occurs with the
surrounding components.

Engine Vibrations
Thus, the optimization of the ENGINE
engine mounts to minimize

TAKOZ
engine –vehicle interaction is
one of the basic subjects of VEHICLE
motor vehicle design.

me
ye Engine Vibrations
ke ce
 In general, the engine has six
degrees of freedom as a rigid
M body in space.
yb
 A third degree of freedom,
ks cs representing the engine bounce
y can be added to the quarter car
w
m model, to study the interaction
kt yo between wheel hop and engine
bounce modes.

Engine Model
Engine is modelled
as a 6 dof (3
translation and 3
rotations) rigid body.

Mount Model
Mounts are modelled by springs with different stiffnesses
in each of the three principal directions.
Z
Springs representing mount

Y
elasticity in three directions
are assumed to be linear X
within a certain deflection r (lateral)
kq q (lateral)
range.
Mount damping is usually kp
p (axiall)
kr
neglected.

Vehicle Models
Models of varying complexity
are used to represent the
dynamic behaviour of the
vehicle.
- Rigid chassis,
- Lumped parameter with 6 dof,
- Modal superposition,
- Detailed vehicle model
(i.e. MBD),
- Interior noise model.

Optimization
Optimization of engine mounts can be defined formally
as the determination of the values of a design
parameters vector, {X}, to minimize a constrained
nonlinear multivariable function, f(X). Thus
Find the values of a design parameters vector, {X}, to
minimize the value of the objective function f(X) subject
to :
gi(X) = 0, i = 1, 2, , p equity constraints,
hj(X) ≤ 0, j = 1, 2, , q inequality constraints
and
Xl ≤ X ≤ Xu upper and lower limits on design
parameters

Design Parameters
Design parameters may be taken as
- The number (usually 3 or 4) of and
locations of the mounts, and/or
- Stifnesses of the mounts in three
directions, and/or
- Orientation of each mount (Euler angles).
Since the mount locations are usually
selected using some other criteria, it is
common to optimize mount stiffnesses and
orientations.

Objective Functions
- Uncoupling roll mode from other modes.
- Uncoupling all modes from each other and keeping natural
frequencies of engine mount system within a frequency
range.
- Minimizing forces transmitted to the vehicle through the
engine mounts.
- Minimizing forces transmitted to the vehicle through the
engine mounts.
- Minimizing displacement or acceleration magnitudes at
prespecified body points.
- Minimizing sound pressure levels at prespecified points in
the passenger cabin.

Engine
Vibrations
 Data set for

a classical
example of
an engine
mounting.
 Ünlüsoy, Y. S., Bilal, H., Çalışkan, K., “Motor Destek Takozlarının

Optimizasyonu”, OTEKON’10, 5. Otomotiv Teknolojileri Kongresi,
Kongre Kitabı, s. 145-156, Bursa, Haziran, 2010.

Engine
Vibrations
 A typical set
of
eigenvalues
and
eigenvectors
of the
example
engine
mounting

Engine
Vibrations
 Eigenvalues
and
eigenvectors
after
optimization
(Uncoupling
all modes
from each
other )

Engine
Vibrations
 Initial
modal
kinetic
energy
distribution

Engine
Vibrations
 Modal
kinetic
energy
distribution
after
optimization

ME 513
VEHICLE
DYNAMICS

Prof. Dr. Y. Samim
Ünlüsoy
CH 7 – Ride Comfort Part 2
ROAD INPUTS TO VEHICLE

Road Inputs to Vehicle
M  An analytical approach to the

y prediction and analysis of ride
motions requires a mathematical
k c model of the vehicle that is
s
z excited by the irregular road
profile.

Road Inputs to Vehicle
 Profile characteristics of actual road

surfaces to be used as road input to
vehicle models may be described as either
 deterministic or
 random
in both the
 time and
 frequency
domains.

Road Inputs  The excitation to a vibratory system
can be classified as either
to Vehicle  deterministic or
 random
In both time and frequency domains.
 Deterministic functions are

identified by the fact that their
values at any future time can be
determined.
 Functions whose value at a future

instant of time cannot be predicted
are called random.

Deterministic Inputs
 As a meaningful deterministic input
for vehicle analysis, one can use the
recorded profile of a representative
road surface.
 However, any given profile may not
be representative of the wide variety
of road surfaces that a vehicle may
encounter.

Deterministic Inputs
 The automotive industry have developed

proving ground surfaces that span the
range of various roads and road conditions
encountered by a given type of vehicle.
 The profile for a set of proving ground
surfaces may then be recorded and these
records may be used as deterministic
“virtual” proving grounds for simulations.

GM Proving Ground

Proving Grounds – IDIADA test track

Proving Grounds

Proving Grounds

RANDOM ROAD INPUTS
Prof. Dr. Y. Samim Ünlüsoy Industrial Fluid Power 12

Random Inputs
 The response of a system

 subjected to random inputs, and/or
 with randomly changing parameters
will also be a random (or stochastic) process.
 In such a case, system response can be
evaluated only in terms of probabilities and a
statistical approach to vibration analysis is
thus necessary.

Random
Inputs
Consider the
surface profile of
a road.
 Measurements of road surface profile along the road

are taken. Each individual set of measurements, xk(t),
started at a different point across the road, is called a
random sample function.
 The ensemble (collection) of all sample functions,
xk(t), constitutes a random (or stochastic) process
and will be denoted by { xk(t) }.

 Time (or temporal) average:  Ensemble average:
of a single random variable xk(t), is the average of a collection of n
belonging to a random process sample functions xk(t), at time t1.
denoted by { xk(t) } is given by 1 n
T
E  xk  t1   = Lim n  xk  t 1 
n k =1
1
E  x k  t   = Lim T   x k  t  dt
T 0  Thus, to obtain the mean value of
the random process at a given
where E stands for "statistical time t = t1, simply sum up the
expectation of ...", E [ xk(t) ] is the values of the n sample functions
mean value of xk(t) and T is the at t = t1, and divide the sum by n.
sample time duration.
Mean Square
 The mean square value of xk(t) is obtained
from the expression
T
1 2
E  xk2  t   = Lim T  xk  t  dt
T0
 The mean square value is a measure of the

energy content of xk (t). The root mean
squared (rms) value is obtained from
1T 2
x rms  Lim T    x k t dt
T0

Variance
 The variance 2 of xk(t) is defined as

σ = E xk  t  - E  xk  t   


2
2
 
where
E [ xk(t) ] and { xk(t) – E [xk(t)] }
represent the static and dynamic component
of the random variable xk(t). Thus, the variance
can be defined as the mean square value of
the dynamic component.
 Variance equation can be simplified to:
 
2
σ = E  x  t  - E  xk  t  
2 2
k

 
2
σ = E  x  t   - E  xk  t  
2 2
k Variance
Note that for a zero mean function the
variance is equal to the expected square
value of the function.
σ 2 = E  xk2  t  
For discrete data, variance is calculated by:
n
x - x
2
i
σ2 = i=1
n -1
x
where is the mean of data.
The standard deviation  is simply the

positive square root of the variance.

Stationary
Process
 When the ensemble average E [xk(t1)] depends on t1, the

random process is called nonstationary.
 When the ensemble average E [xk(t1)], for all possible
averages over the samples, is independent of t1 then the
random process is said to be stationary.
 For a stationary process, the statistical properties are
invariant under a shift of the origin of the x-axis.

Ergodic Process
 If a random process is stationary and if the

temporal mean value E [ xk(t) ] is the same for all k,
then the process is ergodic.
In such a case one single sample function is
representative of the whole random process.
 Obviously an ergodic process is by necessity a
stationary process, while the opposite is not
necessarily true.
 The ergodicity assumption allows the use of a single sample
function to calculate averages describing a given random
process.

Random Road Inputs
 It is necessary to define surface roughness in a
concise and meaningful fashion with a minimum
number of parameters and in a form suitable form
usable for vehicle analysis.
 There has been almost complete agreement on
the choice of Power Spectral Density (PSD) as the
parameter to represent road surface roughness.
 It should be noted, however, that the power
spectral density includes all pertinent information
as to the frequency content of the profile, but
neglects phase relationships.

Auto-
Correlation
Function
 The autocorrelation function provides information about the

dependence of the value of a random variable at one time on
its value at another time separated by  in time domain.
 It is equal to the mean (expected) value of the product x(t)x(t
+ ) over the ensemble of sample functions.
Rx    = E  x  t  x  t + τ   1T
R x k t     Lim T    x k t  x k t   dt
T0

Autocorrelation Function
 Since
Rx() = Rx(-)
the autocorrelation function is an even function of .
 It is also obvious that Rx() must have its maximum
value at  = 0. In this case,
Rx(0) = E [x2]
 It must be noted that if x(t) is periodic, then so will
be Rx().
In this case, the maximum value of Rx() will be
obtained not only at  = 0, but also for  values that
are integer multiples of the period.

POWER SPECTRAL
DENSITY

Power Spectral Density
 The power spectral density function Sx() for a
representative sample function x(t), from an
ergodic random process {x(t)}, is the Fourier
transform of Rx(), divided by 2.
1  -iωτ
Sx ω  =  R x  τ  e dτ
2π 
In this equation  [rad/s] is the temporal frequency.
Then the autocorrelation function is given by the
inverse Fourier transform of the Sx().

R x  τ  =  S x  ω  eiωτ dω



R x  τ  =  S x  ω  eiωτ dω Power Spectral Density

Sx () is always nonnegative and an even

* 
function of .
If, therefore, one lets  = 0 in the above equation,

R x  0  = E  x  =  S x  ω  dω
2

which means that the area under the Sx ()
* curve gives the mean square value of the

random variable x, evaluated over the entire range
of frequencies, -  <  < .


R x  0  = E  x 2  =  S x  ω  dω Power Spectral Density
-
Sx  For a zero mean process, the

area under the Sx() curve gives
* the mean square value of the
random variable x, evaluated
over the entire range of
f1 f2 frequencies, i.e.,
f [Hz]
n [cyc/m] n1 n2 <<
E  x 2  m2  = Sx ? Δω Hz   PSD can be written in a form
 m2  that does not contain negative
 Sx   frequencies. Then it is called a
 Hz 
one sided PSD.
 cyc 
E  x 2  m2  = Sx ? Δn   
 m 
R x  0  = E  x  = 2  Sx  ω  dω
*
2
 m3 
 Sx  
0
 cyc 

Units of Power Spectral Density
Note that there are various definitions for frequency.
In the spatial domain
spatial angular frequency:  [rad/m]
wavenumber: n [cyc/m]
In the temporal (time) domain
temporal angular frequency: w [rad/s]
temporal frequency: f [cyc/s]  [Hz]
Then the displacement power spectral density may
have the units:
In the spatial domain: Sd [m2/rad/m]  [m3]
Sd [m2/cyc/m]
In the temporal (time) domain: Sd [m2/rad/s]
Sd [m2/cyc/s]  Sd [m2/Hz]

PSD of Derived Processes
 If the spectral density Sx() of a stationary

random process is known then the spectral
densities of its derivatives can be found as
follows :
S x  Ω  = Ω 2S x  Ω 
Sx  Ω  = Ω2S x  Ω  = Ω4S x  Ω 
S x  ω  = ω 2S x  ω 
In the time domain
Sx  ω  = ω 2S x  ω  = ω 4S x  ω 

Units of Velocity Power Spectral Density
The velocity power spectral density, i.e., the slope
of the displacement spectral density:

2
  2  
  rad    m    rad   cyc  
2
m 
2
S v m  =  Ω  S  2π n S
  d  rad    cyc   m   d  rad 
  m          
 m   m 
In the time domain
  m 2     
   2
 m    rad   cyc  
2
m 
s     rad  
2 2
Sv   = ω S  2π f  S
 cyc    s   d  cyc    cyc   s   d  cyc 
 s   s   s 
 
Note: disregard any remaining radians.

Sx() White Noise
So  If the PSD is constant for the
whole frequency range of
interest, then the spectrum is
called white, and the random
variable with white spectrum is
 called white noise.
 What is the mean square value of white noise? Is it

possible to generate white noise in the laboratory?
Sx()
So
 Band Limited White Noise
Realistic? Feasible?


Random Road Inputs
 Profile records must be of finite length and therefore it

should be assumed that road profiles are approximate
realizations of stationary random processes (i.e.
statistical characteristics such as mean, variance, etc.
of the road profile do not change with time).
 It is not easy to define S(n) for a given road with great

precision.
 High accuracy is rarely essential, however, since the

aim is to establish a spectral description which can be
taken as typical of a given class of road, rather than to
describe one particular road with precision.

 cyc 
 cyc 
f
 s  Random Road Input
n  
 m  V m 
s A single road surface profile may be used
  
as the road input to a road vehicle,

assuming that the road surface roughness
is ergodic (i.e. the left and right wheels are
subjected to identical inputs). In other
words, any variation of road surface profile
across the road is ignored.
 Note that the
 The road surface profile input is usually
wavenumber
described by a spectral density Sd(n)
contains both
[m2/cyc/m], as a function of the
the velocity and
wavenumber n [cyc/m]. In this form the road
frequency
input does not change with vehicle speed.
information.

 cyc 
Random Road Inputs  cyc 
f
 s 
n  =
 m  V m
 s 
 PSD calculated from profile measurements
must cover the relevant wavenumber range.
 Thus if response spectra are to be
established over the frequency range
 cyc 
0.5 
 0.5 < f < 50 [Hz] and fmin  s   cyc 
nmin = = = 0.01 
Vmax m  m 
50  
vehicle velocities are in the range s
 5 < V < 50 [m/s],  cyc 

50 
fmax  s   cyc 
then wave-number range will be nmax =
Vmin
=
m
= 10 
m 
5  
 0.01 < n < 10 [cyc/m] s

MEASURED
ROAD PSD
 Where profile measurements are
available, power spectral
densities S(n) for a particular
road can be computed directly.
 Note that the basic numerical

calculation of the PSD is
performed by using the FFT of the
sample function of duration T.
2
2 FFT(f )
S(f) =
T

Parameters from Measurements
Where profile measurements are not available and a
spectrum representative of a class of roads is
required, or where analytical convenience is
important, it may be useful to employ a mathematical
model.
It was found, that an expression from the power
spectral densities obtained from the measured road
displacement profiles, provides a good fit to them.
-w
n 
S  n  = S  no   
 no 
In this equation S(n0) is a basic roughness coefficient
and w is a parameter called the waviness.

S(no)
Parameters Surface Quality 10-6 m3/cyc w
from Cement-concrete Very good 0.6 2.29

Good 4.5 1.97
Measurements Average 8.7 1.97
Poor 56.3 1.72
Asphalt-concrete Very good 1.3 2.20
Good 6.0 2.18
-w
n  Average 22.3 2.18
S  n  = S  no    Macadam Good 8.9 2.26
 no  Average 20.8 2.26
Poor 42.9 2.15
Very poor 158.0 2.15
Paved road Good 13.7 1.75
(Cobble-stone) Average 22.8 1.75
Poor 36.4 1.81
Very poor 323.0 1.81
Uncovered (track) Good 31.8 2.25
(foot path) Average 155.0 2.25
Poor 602.0 2.14
Very poor 16300.0 2.14
MATHEMATICAL MODELS for
RANDOM ROAD INPUTS

Mathematical Model for Road Inputs
 ISO 8608* was prepared to facilitate the compilation and

comparison of measured vertical road profile data from
various sources. It therefore specifies a uniform method
of reporting data from one-track and multiple-track
measurements.
 ISO 8608* suggests a straight fit in log-log scale to the
power spectrum generated from the measured road
displacement profile.
* ISO 8608 Mechanical Vibration – Road Surface Profiles –

Reporting of Measured Data, International Standard, 2nd edition,
2016.

 Two alternative versions of the
analytical expression for the road profile
ISO 8608
displacement PSD by ISO 8606 are:
  n 
-w
Gd  no    0.01< n  10  cyc 
Gd  n  =  n
 o no = 0.1 
  m 
 0 otherwise
or  -w
 Ω 
0.1< Ω  100  = 1
Gd  Ωo   rad  We will not
Gd  Ω  =  
 Ωo  o  m  use this!
  
 0 otherwise
Gd(n) , Gd() : displacement power spectral density [m2/rad/m],

Gd(no) , Gd(o) : “basic roughness exponent” [m2/rad/m],
w : “waviness” with a single value 2,
no ,  o : reference “wavenumber” [cyc/m], “angular
spatial frequency” [rad/m]
Gd(n0) [10-6 m2/cyc/m]
Road
class Minimum
value
Mean
value
Maximum
value
w ISO 8608
A (Very using
- 16 32
good) wavenumber
B (Good) 32 64 128
C (Average) 128 256 512
D (Poor) 512 1024 2048
2
E (Very
2048 4094 8192
poor)
F 8192 16384 32768
 cyc 
G 32768 65536 131072 no = 0.1 
H 131072 262144 -
 m 
Note: In ISO 8608 Table C.2, the unit of Gd(n0) on this table
is given as n as [m2/rad/m], possibly a misprint.

Gd(0) [10-6 m2/rad/m]
Road
ISO
class Minimum Mean Maximum
w
value value value
A (Very good) - 1 2
8608
B (Good) 2 4 8
using spatial
C (Average) 8 16 32 frequency
D (Poor) 32 64 128
2
E (Very poor) 128 256 512
F 512 1024 2048
G 2048 4096 8192  rad 

o = 1 
H 8192 16384 -  m 
Note: Gd(n0) = 16 Gd(0), see ISO 8606 Equation (C.8).

ISO 8608
The selection of
the table to be
used will depend
on the units and
the type of the
reference
frequency of the
available data.
Note: The unit of

Gd(n0) on this table
should be
[m2/cyc/m].

ISO 8608
Paved road profiles generally belong to ISO 8608 classes

from A to D. Thus, it is common to use class B road
profiles in simulations of vehicles travelling on roads.
Reference: “Use of generated artificial road profiles in road roughness

Evaluation”, Giuseppe Loprencipe, Pablo Zoccali, J. Mod. Transport v. 25,
n. 1, pp. 24–33, 2017.

 The units and the plot in Figure C.2 of the
standard can be confusing! So, be careful with Example
the units.
----------------------------------------------------------------------
Calculate the displacement PSD in [m2/cyc/m] at
n = 0.01 [cyc/m]
for a worst A-class (or best B-class) road.
 Gd(no) = 32x10-6 [m2/cyc/m],

no = 0.1 [cyc/m].
-2    
2 -2 2
 n  -6  m   0.01  -3  m 
G d = G d  no    = 32x10  = 3.2x10 
n
 o cyc   0.1  cyc 
   
 m   m 

Calculate the displacement PSD in [m2/cyc/s] at
f = 1 [cyc/s] (f = 1 [Hz])
Example
for a mean B-class road. Vehicle speed: V = 60 [kph].
 Gd(no) = 64x10-6 [m2/rad/m],  cyc   m   cyc 
f   V   n 
no = 0.1 [cyc/m].  s  s m 
2
     m 
2  2  -2 2 2 V  
m  m  n   1  s  1  s  m   cyc   s
Gd 
cyc 
= G 
d on         = G  n 
 m  d o  cyc  no2  
  
 cyc 
  n o   V   m  V      m   f  cyc  
 s   m    s 
  
   
m  2
m   m  2  cyc  2
2
1
Gd   = V   G d  no    no  
cyc s cyc  m  2  cyc  2
   
 s   m  f  
 s 
 
 60 
2
-6  m  2 1   2
-5 m
Gd = 
3.6  64x10  cyc   0.1  2  = 1.07x10  Hz 
    1   
 m 

ISO 8608 Road Class Displacement PSD
100 H Note the

units.
-2 G
10
F
E
D
-4
10 C
B
10-6
A
10-8
-2 -1 0 1
10 10 10 10
Wavenumber [cyc/m]

ISO 8608 Road Class Displacement PSD
100
Note the
H
units.
10-2
G
F
E
10-4
D
C
B
-6
10
10-8
10-10
10-1 100 101
Spatial Frequency [rad/m]

Road Surface Velocity Input
Let us use w = 2 and check the PSD of the road

profile velocity excitation.
-2 2
 n    n 
G v  n  =  Ω  G d  n  =  2πn  G d  no   = G d  no   2πn  o  
2 2

 no    n 
G v  n  =  2πn o  G d  n o  = Constant
2
Thus, the road profile velocity input PSD turns

out be constant, i.e., white noise!

 n 
S 0  n   S  no   
w
Alternative Road Input
 no  Models
It is observed that the
equation used in ISO 8608
gives:
lim n0 S0(n)  
which means an infinite
variance will result.
This is not realistic and thus
Breakaway
frequency other models for road surface
roughness with finite values
as n  0 are suggested.
 These models will have a lower slope at low wavenumbers up

to a certain breakaway frequency.

Road Vehicle
PSD model
Output
Shaping
Original Configuration Filters
White First-order Road Vehicle  Transfer function
Output
Noise filter PSD model of a first order
filter. Note that
Equivalent Configuration
V is the corner
frequency*.
Vehicle
White model + Output
Noise Filter 1
H  s 
Final Configuration s  V
* Remember Frequency Response/Bode Plot/ First Order Factor

which is a low pass filter - ME 304 Control Systems.

“First-order Shape Filter” Road Model
Si(n) : Single-sided spatial power
2ασ 2 1 spectral density [m2/rad/m],
Sd  Ω  =  : spatial angular frequency
π α 2 + Ω2 [rad/m],
0Ω<  : spatial breakaway frequency,
2 : variance of road irregularity.
 rad   cyc 
Ω  = 2πn m 
 m    Typical parameter values:
Road Type  rad/m 2 m2
Very good asphalt 0.150 9x10-6
Good asphalt 0.225 44x10-6
Filter Average asphalt 0.300 105x10-6
White Road
Noise PSD Poor asphalt 0.375 190x10-6
Paved road 0.450 300x10-6
Dirt road (terrain) 0.750 750x10-6
“First-order Shape Filter” Road Model
-α ξ
Correlation Function: R ξ  = σ e 2
2 
Single Sided Power
Sd  Ω  =  R  ξ  cosΩξdξ
Spectral Density: π 0
2ασ 2 1
Sd  Ω  =
π α2 + Ω2
0Ω<

Random Road Profile Input
in Time Domain
For ride comfort analysis in time domain, the
random road profile displacement with the selected
model should be obtained as a function of time for
a given vehicle speed.
In ride comfort analysis, it is common to select one
of the road classes from ISO 8608 and obtain the
corresponding road profile displacement.
Note, however, that this time function will not be
unique as there may be infinite number of
functions that have the same PSD.

“First-order Shape Filter”
Road Input Model – time domain
For this model, the variation of the road profile height
with time or with the length of the sample road section,
can be obtained by the solution of the first-order
differential equation
y 0 (t) = -αVy 0 (t) + w(t)
where w(t) is zero mean white noise with covariance of
E [ y0(t) y0() ] = 2V2 d(t- ).

The solution of this differential equation will provide a
zero mean random time function with a variance of 2
representing the road surface profile.
Note that the unit of d(t) is [1/s].
Thus, the unit of covariance is (m/s)^2.

“First-order Shape Filter”
Road Input Model – time domain
State equation form of the shape filter
y 0 (t) = -αVy 0 (t) + Vw(t)

is given by:
 2αV 2 
y 0 = Ay 0 + Bw = (-αV)y 0 + V  σ  w
 π 
Here w is the zero mean white noise of variance
one.
Note: Matlab command “randn” may be used to generate
normally distributed random numbers with zero mean and
variance one.

First-order Shaping Filter Model
 
2 m 2
2αVσ 1
S1  f  =  
π  αV 2 +  2πf 2  rad 
 s 

ISO 8608 Road Input Model – time domain
A similar approach, to the generation of artificial
road profiles in time domain for ISO 8608 road
classes, is the solution of the first order
differential equation
y 0  t  =  -2πnb V  y 0  t  + 2πn o  G d  no  V w(t) 

for y0. αV = 2πfb = Ωb V = 2πnb V
Here nb=  / 2 is 2αVσ 2 = 2πno G d  no  V =  2πno  G d  no  V

2
 
the break-away
 rad   m   m 
2 2
 rad   m  2 m 2rad
wavenumber.  m   s  m  =  m    rad    s   s
       
  m  

The road profile height h, as a function of

road distance x, can be described through
summing simple harmonic functions
according to:
h i  x  = A i cos  2πn i x + φ 
where Ai is the amplitude, n is the

wavenumber, and j is the phase angle.

If the PSD of function of road profile height is known,
then it is possible to generate an artificial road
profile described by:
N
h  x  =  2ΔnG d  ni cos  2πni x + φi 
i=0
0≤x≤L : distance along the road [m]

L : road length [m],
N : number of samples
sI = L / (N-1) : sampling interval [m]
Dn = (nmax - nmin) / (N-1) : spatial frequency increment [cyc/m]

(0 ≤ j ≤ 2) : random phase angle from a uniform
distribution. Matlab function rand can
be used to generate j.
Artificial road profiles for ISO 8608 road

classes can be generated by:
N
 no 
h x =  k
Δn 2 10  -3
 cos  2πni x + φi 
i=0  ni 
k (3 ≤ k ≤ 9 ) : constant corresponding to the
ISO 8608 road classes from A to G.
(0 ≤ j ≤ 2) : random phase angle from a uniform
distribution. Matlab function rand can
be used to generate j.


On the Properties of
PowerSpectral Density
1. The area under the PSD vs frequency plot is
equal to mean square value for a zero mean
process.
Remember (slide 16):
The variance 2 of xk(t) is defined as
 2 = E [ (xk(t) - E [xk(t)])2 ]
which can be simplified to:
2 = E [ xk2(t) ] - ( E [xk(t)] )2
So, the variance is equal to:
2 = mean square value- square of mean

 -w
 n 
Gd  n  =  d  o   no 
G n 0.01< n  10 On the Properties of

 0 otherwise Power Spectral Density
2. Note that the ISO 8608 road profile

expression, in the general case, will have
infinite variance. However, since it is
confined within the wavenumber range
0.01 < n ≤ 10
it will have finite variance, given by
10
σ =
2
d G d  n  dn
0.01

On the Properties of
PowerSpectral Density
3. To find the variance of a random process
such as road profile height, area under the
PSD vs frequency should be calculated.
This is valid for all frequency definitions.
νH
σ =  S x  n  dn
2
x nL
ΩH
σ =
2
x S x  Ω  dΩ
ΩL
fH
σ =  S x  f  df
2
x fL

 -w
 n 
Gd  n  =  d  o   no 
G n 0.01< n  10 On the Properties of

 0 otherwise Power Spectral Density

4. σ =
2
d G d  n  dn For the ISO8608 model:
0.01
-w
10 n  10 dn
σ =  G d  no    dn = G d  no  no 
2
d
2
0.01 n 0.01 n 2
 o
10
1
σ = G d  no  n -
2
d
2
o= 99.9G d  no  no2
n 0.01
Thus for the mean E class road of ISO 8608:
σ = 99.9  4094x10   0.1 2
2
d
-6
= 4.09x10-3 m 2 
rms = 4.09x10-3 = 0.064  m 

 2ασ 2 
Sd  Ω  =   2
1
2
0<Ω< On the Properties of
 π α +Ω

4. σ =  S d    d
2 For the shape filter model let:
d 0
  2ασ 2   1   2ασ 2
  dΩ
σ d2 =    2 2 
dΩ =   0 2 2
0
 π  α + Ω   π  α + Ω

2  2ασ  1 -1 
2
Ω   2σ 2
 π  2
σd =   tan   =    = σ
 π  α α
 0  π  2 
Thus for the good asphalt road of shape function model:
σ d2 = 44x10-6 m 2 
rms = 44x10-6 = 0.0066 m 

Response of Linear Systems to
Stationary Random Inputs

H = hreal +ihimaginary
H* = hreal - ihimaginary Response of Linear Systems to
The response y, of a linear system to a single

input x is described by the relation:
S y  ω  = H iω  S x  ω 
2
=H iω  H * iω  S x  ω 
where H() is the complex frequency response
function of the system and H*() is its complex
conjugate. Si() is the power spectral density
(PSD).

Response of Linear Systems to
The mean square response is then obtained from :
 
E  y =  S y  ω  dω = 2  S y  ω  dω
2
 
- 0
 2
E y = 2  H  iω  Sx  ω  dω
 2
 
0

Example 7.1
Determine the mean-square value of

k c the mass displacement x, if the
x power spectral density of the
exciting force f(t) is white, i.e.,
m constant at So for all frequencies.
f(t)

Example
The equation of motion is given by :
 + cx + kx = f(t)
mx
k c
x Take the Laplace transform of all the
terms.
m
f(t)
ms 2 X  s  + csX  s  + kX  s  = F  s 

ms2 X  s  + csX  s  +kX  s  = F  s 
Example
Solve for the Laplace transform of the mass displacement.
 1 
X  s =  2 F  s
 ms + cs +k 
Then the transfer function is given by:
X  s 1
H(s) = =
F  s  ms2 + cs +k
The frequency response function (FRF)
is obtained by inserting
1
s = i H(iω) =
k - mω 2 + icω
in the transfer function.

1
H(iω) =
k - mω 2 + icω Example
To obtain the frequency response expression in
the general form a + ib, multiply and divide by
the complex conjugate of its denominator.

H(iω) = 
1   k - mω 2
 - icω 


 k - mω + icω   k - mω  - icω 
2 2

H(iω) = 
 k - mω 2   
 i cω


    
 k - mω 2 +  cω    k - mω 2 +  cω  
2 2 2 2


 k - mω  - icω
2
 k - mω  +icω
2
H(iω) =
k - mω  + cω 
2
2 2
H * (iω) =
k - mω  + cω 
2
2 2
Example
The square of the frequency response function is obtained as :
 k - mω2  - icω   k - mω2  +icω 
H iω  = H iω  H * iω 
2
H iω  =  
2
    
 k - mω2 +  cω    k - mω2 +  cω  
2 2 2 2

1
H iω  =
2
k - mω 
2
+  cω 
2 2
Then the expression for S f  ω  = So

the spectral density of
So
Sx  ω  = H  iω  Sf  ω  =
2
the mass displacement
k - mω 
2
+  cω 
2 2
can be found.

1
H(iω) =
k - mω 2 + icω Example
+∞ +∞
 S  ω  dω =  H iω  Sf  ω  dω
2
E  x  =
2
x
-∞ -∞
+∞ 2
1
E  x 2  = So  2
dω
-∞
k +icω - mω
To find the mean square value the

integral should be evaluated.

Example
 2
The general integrals of the form I n   H n  d

B o  iB 1  i2 B 2  ...  in  1 B n  1

H n  
with A o  iA 1  i2 A 2  ...  in A n
are available up to n = 7*.
*D.E. Newland, An Introduction to Random Vibrations, Spectral &

Wavelet Analysis, Third Edition, Dover, 2005.

Example
*Mansour Ataei, Ehsan Asadi, Avesta

Goodarzi, Amir Khajepour, and Mir Behrad
Khamesee, “Multi-objective optimization of a
hybrid electromagnetic suspension system for
ride comfort, road holding and regenerated
power”, Journal of Vibration and Control, 1–12.

 2
I n   H n  d

Example
The comparison of the Bo + iω  B1
H2 iω  =
general form for n = 2 A o + iω  A1 - ω2 A2
with the FRF for the 1

H(iω) =
system in the example k  iω  c  ω2m
will give B0 = 1, B1 = 0,
A0 = k, A1 = c, A2 = m

I2 

 A o B 12  A 2 B o2  Example
A o A1 A 2
With B0 = 1, B1 = 0,
A0 = k, A1 = c, A2 = m
π k  0  c +m 1  π
 2 2
I2 =  =
kcm kc
and since E  x 2  = S f I2
Finally the mean square

π
of the mass displacement  2

E  x  =   So
is obtained as:  kc 

Example
The exciting force spectral density function, square

of the transfer function magnitude, and the spectral
density of the displacement of the mass are shown,
graphically.

CH 7 – Ride Comfort - Part 2b
Vehicle Response to Random Road Input

Vehicle Response to Random
Road Input
The response of a linear vehicle model to a single
road input is described by the relation:
S y  ω  = H iω  S x  ω 
2
=H iω  H * iω  S x  ω 
where H() is the complex frequency response

function of the vehicle model between the input
and the output of interest. H*() is its complex
conjugate. Sx() and Sy() are the PSD of the road
profile input and the output of interest.

Road Input
Let us take the linear quarter car model as the

reference.
If the subject of interest is the ride comfort, then
the parameter of interest will be the sprung
mass acceleration.
Thus, in order to calculate the response of a
vehicle to random road input, the transfer
function between the road profile displacement
or road profile velocity, and the acceleration of
the sprung mass should be calculated first.

Road Input
The transfer function of a linear system

can be calculated in several different
ways.
One approach to the calculation of the
transfer function of a quarter car model
was given on slide 56, and another was
explained on slide 66 of the Lecture
slides CH7 - Ride Comfort – part 1.

Road Input
For ease of calculation and

generality, the Matlab function
“bode” is commonly used.
The basic procedure to obtain the
response of a vehicle model traveling
on a random road surface is
illustrated in the following example.

V
y
m=250 kg
k=16 kN/m c=900 N/m/s Example 7.2

y0
- Obtain and plot the transfer functions between the sprung mass
acceleration and the suspension travel, and the road velocity
input,
- Obtain and plot the response of the single dof vehicle model, at
60 [kph] on an ”average asphalt” road, using the shape filter road
input PSD model, in terms of :
 acceleration of the sprung mass, and
 suspension travel.
Consider a frequency range of 0.1 to 20 [Hz].

V
m=250 kg
y
Example
Equation of motion:
k=16 kN/m c=900 N/m/s
  cy  ky  cy 0  ky 0
my
y0
State equation form:
c k c
x1  y x 1  
y   x1  x 2  y 0
m m m
x2  y 0  y x 2  y 0  y  y 0  x1
 x1   y   c
 x 1   
k
 x 
c
 
      m  1
m    m  y 0

x y  y
 2  0   x 2    x2   
 1 0 1

V
y
m=250 kg
Example
k=16 kN/m c=900 N/m/s
y0
Listing of the code:
% Single dof vehicle model with random road Input

% ----------------------------------------------------------------------
% m = sprung mass (body mass) [kg]
% c = suspension damping coefficient [N/m/s]
% k = suspension spring stiffness [N/m]
% V = vehicle speed [kph]
m = 250; k = 16000; c = 900;

  c k c
 1  
x
  m   1  
x
m    m w w  y 0 Example
 x 2   x
0   2   1 
 1
y  Cx  D u
 y1   1 0   x1  y 
   x  y  y
y
 2  0 1  2 0 
% Single dof vehicle model
% W is the disturbance distribution matrix in
% {dx/dt} = [A]{x} +[W]{w}
V
y
A = [-c/M k/M; -1 0]; m=250 kg
W = [c/M; 1];
% {y}=[C]{x}+[D]{u} k=16 kN/m c=900 N/m/s

C = eye(2);
D = zeros(2,1); y0

 c k c
 x 1     x   
  m
 x 2 
m      m  y 0
 x
1
Example
 1 0   2   1 
% Undamped natural frequencies and damping ratios

wn=sqrt(k/M);
fn = wn / (2*pi);
zeta = c /2 / M / wn;
fprintf('Body bounce frequency is %4.2f [Hz]\n’, fn)
fprintf('Body bounce damping ratio is %4.3f \n', zeta)
% Generate frequency vectors
N=200;
f = logspace(log10(0.1),log10(20),N);
w = 2*pi*f;

% Obtain the transfer functions Example
j = sqrt(-1);
s = j*w;
% Transfer functions: T(:,1)ydot/y0dot and T(:,2)(y0-y)/y0dot
[T,p] = bode(A,W,C,D,1,w);
% Transfer function  ydoubledot/y0dot is needed
for i = 1:N
Tsma(i) =s(i)* T(i,1);
Tst(i) = T(i,2);
end
Tsma=abs(Tsma);
% In decibels
Tsma_db = 20*log10(Tsma); Tst_db = 20*log10(Tst);

% Plot transfer function  ydoubledot/y0dot Example
figure(1)
semilogx(f,Tsma_db, ‘r’,’LineWidth’,2)
hold on
% Plot transfer function  ydoubledot/y0dot
semilogx(f,Tst_db,’b’,’LineWidth’,2)
xlabel('Frequency [Hz]’)
ylabel(‘Magnitude [dB]’)
legend('[ydoubledot/y0dot]','[(y0-y)/y0dot]’,
'Location','Best')
title('Transfer functions')
grid

Example
% Vehicle speed
V = 60/3.6; % [m/s]
% Road Displacement Power Spectral Density for “Average” asphalt
road is given by the expression :
% Syo(n) = S(no)*(n/no)^(-w) [m^2/rad/m]
Sno = 105e-6; % [m^2/rad/m]

Sno = 2*pi*Sno; % convert to [m^2/cyc/m]
Sfo = Sno / V; % convert to [m^2/cyc/s]
no = 0.1 ; % [cyc/m]
fo = no*V; % convert to [cyc/s]
ww = 2;

Example
%---------------------------------------------------------------
for i = 1:N
% PSD of road displacement
PSDd(i) = Sfo*(f(i) / fo)^(-ww); % [m^2/Hz]

% PSD of road velocity
PSDv(i) = (2*pi*f(i))^2 * PSDd(i); % [(m/s)^2/Hz]
end

Example
%-----------------------------------------------------------------
% Plot Road PSD - Displacement and velocity
figure(2)
loglog(f,PSDd,'b:',f,PSDv,’r--’)
ylabel('PSD [(m^2)/Hz, (m/s)^2/Hz]’)
legend(‘Displacement’,‘Velocity’)
title(‘Road Input’)
grid

Example
%-----------------------------------------------------------------
for i = 1:N
% PSD of sprung mass acc. in [(m/s^2)^2/Hz]
PSD_sma(i) = Tsma(i)^2*PSDv(i);
% PSD of suspension travel in [m^2/Hz]
PSD_st(i) = Tst(i)^2*PSDv(i);
end

Example
%-----------------------------------------------------------------
figure(3)
% Plot sprung mass acceleration and suspension travel PSD
loglog(f,PSD_sma, f,PSD_st,’LineWidth’,2)
ylabel('PSD [(m/s^2)^2/Hz], [m^2/Hz]’)
legend(‘sprung mass acceleration’,
’Suspension travel’)
title('Vehicle Response')
grid

Example
Transfer functions
30
20
10
0 [ydoubledot/y0dot]
[(y0-y)/y0dot]
-10
-20
-30
-40
-50
-1 0 1
10 10 10
Frequency [Hz]

Example
-1
Road Input
10
Displacement
Velocity
Note that the 10-2
velocity of the
-3
road profile is 10
white noise.
10-4
See slides 24
and 40 of 10-5
CH7_Ride -6
10
Comfort_Part
2a. 10-7
10-1 100 101
Frequency [Hz]

Example
Vehicle Response
100
Take a look at 10-2 Sprung mass acceleration

Suspension travel
the slide 60 of
CH7_Ride 10-4
Comfort_part
2a.
10-6
10-8
10-1 100 101
Frequency [Hz]

EVALUATION OF RIDE
QUALITY

Evaluation of Ride Comfort
 The ride comfort of a vehicle is related mainly
to the acceleration to which the driver and the
passengers are subjected to.
 Thus, in using basic (simplified) mathematical
models of vehicles, ride comfort is based on
the acceleration levels of the sprung mass.
 Further investigations may be carried out by
adding an additional mass-spring-damper
system representing the driver/passenger
seat.

 Although a large number of criteria to evaluate ride
quality for vehicles travelling on roads have been
proposed, the use of ISO 2631-1 “Evaluation of
human exposure to whole-body vibration” by the
International Standards Organization has been
commonly accepted.*
 It is based mainly on existing data regarding
human exposure to sinusoidal vibrations and thus
is directed towards the evaluation of sinusoidal
and periodic vibrations.
* The original 1985 version of the standard has been withdrawn.
However, in the new version, it is stated that “… the exposure boundaries
recommended in ISO 2631-1:1985 were safe and preventive of undesired
effects. Therefore, it is still commonly used.

 ISO 2631 specifies comfort boundaries for
acceleration described by the rms values
calculated in specified 1/3 octave bands for a
number of central frequency values.
 The so-called reduced comfort boundaries
represent the limiting rms values of acceleration in
the frequency range from 1 to 80 Hz. Here reduced
comfort boundaries for for 1, 16, 25 minutes and 1,
2.5, 4, and 8 hours exposure only are given.
1/3 octave band
For a selected fcenter value, the band is between:
flower = fcenter 21/6 , fupper = fcenter * 21/6 = 21/3 flower

ISO-boundaries for Vertical Vibrations
Frequency a [m/s2]
[Hz] 1 16 25 1 hr 2.5 4 hrs 8 hrs
min mins mins hrs
1 1.778 1.349 1.127 0.749 0.444 0.337 0.200
4 0.889 0.673 0.571 0.375 0.225 0.168 0.100
8 0.889 0.673 0.571 0.375 0.225 0.168 0.100
80 8.89 6.73 5.71 3.75 2.25 1.68 1.00
 Reduced comfort boundaries represent the limiting rms

values of acceleration in the frequency range from 1 to 80 Hz.

10
Vertical Vibrations
Acceleration [m/s2] rms
1 min
1
16 mins
25 mins
1 hr
2.5 hrs
4 hrs
8 hrs
0.1
1 10 100
Frequency [Hz]

ISO-boundaries for Lateral Vibrations

Frequency a [m/s2]
[Hz]
1 16 25 1 hr 2.5 4hrs 8 hrs
min mins mins hrs
1 0.635 0.476 0.397 0.270 0.159 0.113 0.071
2 0.635 0.476 0.397 0.270 0.159 0.113 0.071
80 25.4 19.1 15.9 10.6 6.4 4.4 2.9
 Reduced comfort boundaries represent the limiting rms

values of acceleration in the frequency range from 1 to 80 Hz.

100
Lateral Vibrations
10
Acceleration [m/s2] rms
1
1 min
16 mins
25 mins
0.1 1 hr
2.5 hrs
4 hrs
8 hrs
0.01
1 10 100
Frequency [Hz]

 The equivalent reduced comfort boundaries

represent the limiting values of acceleration
for various hours of exposure to sinusoidal
vibrations, in the frequency range from 1 to
80 Hz.
 The vibrations in motor vehicles are
random, however, and thus the standard has
also been converted to a form usable for
data presented in the form of power spectral
density of the random vehicle vibrations.

ISO-boundaries for Vertical Vibrations

Frequency PSD [(m/s2)2/Hz] PSD [(g2)/Hz]
[Hz]
1 hr 4 hrs 8 hrs 1 hr 4 hrs 8 hrs
1 2.42 0.489 0.173 2.5E-02 5.1E-03 1.8E-03
4 0.152 0.031 0.011 1.6E-03 3.2E-04 1.1E-04
8 0.076 0.015 0.005 7.9E-04 1.6E-04 5.2E-05
80 0.757 0.153 0.054 7.9E-03 1.6E-03 5.82E-04
 Note that 8 hour exposure values can be obtained by

multiplying 1 hour exposure values by 0.071.

Vertical Vibrations
10
1 hr
4 hrs
1
8 hrs
PSD [(m/s2)2/Hz]
0.1
0.01
0.001
4 8 80
1 10 100
f [Hz]

ISO-boundaries for Lateral Vibrations
Frequency PSD [(m/s2)2/Hz] PSD [(g2)/Hz]
[Hz] 1 hr 4 hrs 8 hrs 1 hr 4 hrs 8 hrs
1 0.314 0.055 0.022 3.3E-03 5.7E-04 2.3E-04
2 0.157 0.027 0.011 1.63E-03 2.9E-04 1.1E-04
80 6.10 1.07 0.441 6.3E-02 1.1E-02 4.6E-03
 Note that 8 hour exposure values can be obtained by

multiplying 1 hour exposure values by 0.071.

Lateral Vibrations
10
1.0
PSD [(m/s2)2/Hz]
0.1
0.01 1 hr
4 hrs
8 hrs
0.001
2 80
1 10 f [Hz] 100

 The ISO reduced comfort boundaries may be

superimposed on a plot of the acceleration
power spectral density of the vehicle body.
 Ride quality will be acceptable if the
acceleration PSD falls below the selected
reduced comfort boundary throughout the
range of interest.
 If the acceleration PSD rises above the selected
reduced comfort boundary anywhere within the
range of interest, this will indicate inadequate
ride quality.

ME 513
VEHICLE
DYNAMICS

Prof. Dr. Y. Samim
Ünlüsoy
CH 7 - Ride Comfort Part 3a
Passive Suspension
Parameter Optimization

Passive Suspension Parameter Optimization
 The performance of a suspension system

will be measured by the
 ride comfort,
 suspension travel, and
 road holding
characteristics.

 The acceleration of vehicle body,

either as
 the rms values, or
 the power spectral density of
acceleration response
over a frequency range, is usually
taken as an indication of ride comfort.

Passive Suspension Parameter
Optimization
 Suspension travel (deflection or

stroke) is a design constant and
its extreme values are dictated
by the suspension working
space requirements.

Road holding characteristics depend on the

tire-road contact force characteristics and
relates to
 braking,
 acceleration, and
 handling
performance of the vehicle.
For good road holding, tire-road contact force
variations from the static value should be
minimized.

 The response of dynamic systems

to random input can be investigated
using
 power spectral density analysis
in the frequency domain*, or
 covariance analysis in time
domain.
* See “Response of Linear Systems to Random Inputs”

Approach
(frequency domain)

Quarter Car Model PSD approach in
frequency domain
 y 2 
y - y 
x =  y1 - y2 
 0 1
 y 1  x  =  A x + W y 0
 c ks c
- M 0 

M M
  Note that the road input to
 -1 0 0 1  the vehicle model is the
  
A = 
0 0 0 -1
  road surface profile slope,
 c
-
ks kt c
-  y 0 , and not the road
 m m m m
surface displacement y0 .
0 
0 
W =  
1 
0 

 The transfer functions for sprung mass acceleration,
suspension travel, and tire deflection for the quarter
car model were established previously* as:
s 2 Y2 (s) sk t A s = jω
H sma (i)  
sY0 (s) BC  A 2
A  cs s  k s
Y (s) - Y1 (s) k (A - B)
Hst (iω) = 2 = t
sY0 (s) 
s BC - A 2  B  Ms 2  c s s  k s
Y1 (s) - Y0 (s) (k t - C)B + A 2 C  ms 2  c s s  k s  k t

H td (iω) = =
sY0 (s) 
s BC - A 2 
* See CH7 part 1

Optimum Parameters for Ride Comfort
In the optimization of suspension parameters
two aspects are of importance.
- Ride comfort evaluated according to
sprung mass acceleration,
- Ride safety evaluated according to tire
tire deflection.
The design parameters of interest are the
- suspension spring stiffness,
- damping coefficient of suspension damper,
- tire stiffness, and
- unsprung mass.

Optimum Parameters for Ride Comfort
 For most practical purposes, a white noise
road surface velocity input to the vehicle can
be assumed. Then the constant power
spectral density of the input can be taken out
of the integral.
 The mean square sprung mass acceleration
is then given by :
   2
E  2
y 2 = S0 - Hsma (iω) dω
 

E  y 22  = S0Isma
  Optimum Parameters
 2
Isma = -  Hsma (iω) dω for Ride Comfort
 The result of the remaining integral

is given as :
k k 2
1 
I sma    t c s  s m  M  
 M 2 M2 c s 

k k 2
1 Optimum Parameters
I sma    t c 
2 s
s m  M  
 M M 2 c s  for Ride Comfort
 Obviously in this case the solution is trivial

and a suspension spring of zero stiffness
will minimize the mean square sprung mass
acceleration.
 The damping constant of the damper to
minimize the above expression can be
obtained as :
2
k s (M  m)
c s2 
kt

Optimum Parameters for Suspension Packaging
Suspension Travel
 In a similar manner, the mean square

suspension travel is given by:
 2  2
E (y 2 - y1 ) = S0 -  Hst (iω) dω
 

E  (y 2 - y1 ) 2  = S0Ist Optimum Parameters for
 
Suspension Packaging
Suspension Travel
 The remaining integral is (M+m)
given as :
Ist  
cs
 Thus the mean square suspension travel

 does not depend on the suspension
stiffness, and
 is inversely proportional to the
damping constant of the suspension
damper.

Optimum Parameters for Ride Safety
Tire Deflection
 Finally, the mean square tire deflection is given
by :
 2  2
E (y1 - y 0 ) = S0 - H td (iω) dω
 
and the remaining integral is:
 M  m  k 2 2k k m k 2
m k t cs 
I td    s  s t  t  
 M cs2 Mcs M  m  cs M  m 2
M 2

Optimum Parameters for Ride Safety
Tire Deflection
 The suspension stiffness and damping
constant to minimize the remaining
integral are then given as :
k t Mm
ks 
(M  m) 2
2
cs 
k 2
s M  m  
2k s Mm

k t M 2
m
kt (M  m) (M  m) 2

 The results from the procedure followed

so far shows that the optimum values for
the suspension parameters vary
according to the requirements from the
suspension. The results, however, are in
conflict, e.g., a suspension optimized with
respect to road holding will degrade ride
comfort and vice versa.
 Thus a suspension system consisting of
fixed parameter passive components can
only offer a compromise in performance.

 One option to determine the best possible
compromise is based on the multiobjective
optimization of the suspension parameters.
 In the simple case of the quarter car model,
suspension stiffness and damping are the
variables to optimize.
 In the optimization formulation, an objective
function to be minimized will be defined.
Weights are assigned to each term of the
objective function representing the different
requirements, i.e., ride comfort, suspension
travel, and road holding.

CH 7 - Ride Comfort Part 3b
Active Suspensions
Reading assignment:
- Hedrick, J.K., Butsuen, T., (1990) ‘Invariant properties of automotive
suspensions’, Proc. Instn. Mech. Engrs., v. 204, pp. 21-27.

Active Suspensions
 The requirements from vehicle
suspensions with respect to vehicle
handling and ride comfort are in
conflict.
 Good handling requires a stiff
suspension, whereas a soft suspension
provides better ride comfort.
 Therefore, the design of vehicle
suspensions represents a compromise
solution.

Classical Passive Suspensions
 Classical passive suspensions, i.e.,

suspension systems containing only
passive components have the
disadvantages that they can only :
 store or dissipate energy
through springs and dampers
and thus do not have external
power supplies,
 generate forces in response to
local relative motion only.

Active Suspensions
 Active suspensions do not suffer these
limitations as they can :
 continuously supply and modulate the
flow of energy, and produce forces not
restricted by the energy previously
stored in the springs,
 generate forces in response to many
variables, some of which may be
measured anywhere in the vehicle
system or generated.

Active Suspensions
 With the use of active suspensions
 low natural frequencies (good for ride
comfort) together with small static
deflections,
 low dynamic deflections,
 high speed of response,
 suspension characteristics unaffected
by loading
are the expected benefits.

Active Suspensions
 Active suspensions have the
disadvantages that active
controllers add to the cost and
complexity (installation and
Bose Electromagnetic maintenance) of the vehicle, at the
Active Suspension
same time reducing the reliability
Watch1_1min55s of the system.
Watch2_31s
Watch3_25s  However, the advantages of active
suspension systems overweigh
the disadvantages.

y Active Vibration
M
Isolator
u
 A very simple and general
application of active vibration
u z isolator can also be interpreted
as the simplest possible example
of an active suspension.
 The system may also represent the case of a machine

tool to be isolated from ground vibration in a
workshop for better performance.

Sprung
M
y Active Vibration Isolator
mass
u
 Let us consider a rigid mass
(say, the sprung mass of the
Suspension
vehicle) which is to be
u z isolated from a disturbance
(the road surface input to the
vehicle).
 u is the control force applied to the sprung mass

by the suspension system.
 The displacement of the sprung mass and the
road profile are represented by y and z,
respectively.

Sprung y Active Vibration
mass M
Isolator
u
Suspension  The equation of motion is
simply :
u z
My  u
 The state variables are selected as the velocity of the

sprung mass
x1  y
and the suspension travel, sometimes called the "rattle
space".
x2  y  z

Sprung y Active Vibration
M
mass Isolator
u
 Taking the derivatives of
Suspension
the state variables, and
u z denoting the road surface
velocity as w,
z  w
1
x 1  y  u
M
x 2  y  z  x1  w

1
x 1  y  u
M Active Vibration
x 2  y  z  x1  w
Isolator
 Taking u as the control input and w as the disturbance,

the state equations can be written in the standart form :
x  =  A x +  Bu +  W ξ

 x 1  0 0  x1   1  0
      M u    w
x 2  1 0 x 2   0   1
 Note that the disturbance input is road surface velocity !

Quadratic
Performance Index
 The object here is to determine the required

control input such that a certain performance
index is minimized. Here a quadratic performance
index is to be minimized.

J = 0 x
T
 Qx + u  R u T
 dt
x2 terms u2 terms

J = 0  x T
 Q x + u T
 R u  dt Quadratic
Performance
Index
 [Q] is a positive semi-definite real symmetric matrix
and [R] is a positive-definite real symmetric matrix.
 [Q] and [R] are weighting matrices which determine the
relative importance of errors in the system states and
energy expenditure by the control input.
 Since [Q] is allowed to be positive semi-definite, some
of the state variables may not appear in the quadratic
performance index J.

Quadratic Performance Index
 For the vibration isolation problem being analyzed, let
us try to minimize the acceleration of the sprung mass
and the suspension travel. Then by using
 0 0
Q   
0 1
 the first term in the quadratic index becomes
 0 0  x 1 
x1 x 2  
 x   x 2
2
 0 1  2 
 Thus, the first term in the quadratic performance index
results in the minimization of x2, which is the
suspension deflection.

 If now, [R] is taken to be equal to a parameter, i.e.,
R  = ρ
then, the second term in the quadratic index becomes
u  R  u = u ρ  u = ρ u2
so that the control input u is minimized as controlled

by the weighting factor r. This will result in the
minimization of the acceleration of the sprung mass,
since it is directly proportional to the control input.

 The quadratic performance index, to

minimize the sprung mass acceleration and
suspension travel, takes the form :
  
J = 0 x 22 + ρu 2 dt = 0 x 22 + ρM 2
y 2 dt
 The value chosen for the parameter r is

going to specify the relative importance of
errors in minimizing the two variables.

Control Law
 From the quadratic optimal control

theory, the control law minimizing the
quadratic performance index J is
given as:
u = -  K x
where the feedback matrix [K] is given
by:
-1 T
[K]  [R] [B] [P]

{u} = - [K] {x} [K]  [R]-1 [B]T [P] Control
Law
 [P] is a real, symmetric, and positive definite
matrix which satisfies the algebraic matrix
Riccati equation (ARE) :
[A]T [P]  [P] [A] - [P] [B] [R]-1 [B]T [P]  [Q]  0
Note: Matlab command lqr calculates the
optimal gain matrix K in the state-
feedback law u = -Kx which minimizes
the cost function
J = Integral {x'Qx + u'Ru + 2*x'Nu}dt.
[K,S,E] = lqr(SYS,Q,R,N)
[A]T [P]  [P] [A] - [P] [B] [R]-1 [B]T [P]  [Q]  0
Matrix Riccati Equation
 p1 p2 
P  
p 3 
 If the matrix [P] is taken to be
p 2
 the matrix Riccati equation is obtained as :
 0 1   p1 p 2   p1 p 2   0 0  1  p1 p 2   1   1   p1 p 2  0 0  0
 0 0 p     0    
  2 p 3  p 2 p 3  1 0 r p 2 p 3   M M  p 2 p 3  0 1 0
 0 

p 2 p 3  p 2 0 1  p2 p 1p 2   0 0 
   1 
2  0  1
0
 0  p 3 0 r M 2 p 1p 2 p 2    Solution of
Matrix Riccati
Equation
1
2p 2  p 12  0
rM 2
 Three independent equations
in the three unknowns p1, p2, 1
p3  p p 0
2 1 2
and p3 are : rM
1
 p 22  1
rM 2

 From the last equation, p2 can be Solution of
solved in terms of M and r.
Matrix Riccati
p 2  Mr 1 / 2 Equation
 Then, the first equation will
1
provide an expression for p1. 2p 2  p 12  0
rM 2
3/4
p 1  M 2M r 1
p3  p p 0
2 1 2
 Finally, p3 can be obtained using rM
the second equation. 1
 p 22  1
rM 2
p 3  2M r 1 / 4

M 2M r 3 / 4 Mr 1 / 2 
P    1/ 4 Optimum
 Mr1 / 2 2M r 
Feedback
[K]  [R]-1 [B]T [P]
Matrix
 The optimum feedback matrix can

now be calculated.
1 1   M 2M ρ 3/4
Mρ1/2  
K  = ρ  M 0   = 2M ρ1/4 ρ-1/2 
 
   Mρ1/ 2 2M ρ1/4 

{u} = - [K] {x}
 K    2Mr 1/4 r 1/2  Control Input

 The control input is found as :
-1/2   x1 
u = -  2Mρ -1/4
ρ  
 x 2 
= - 2Mρ -1/4
x1 - ρ -1/2
x2

u   2M r  1 / 4 x 1  r  1 / 2 x 2 Closed Loop System
 x 1   x1   0  Matrix & State
  =  A  - B K    +   w
 
x 2  x 2  -1 Equation
 The closed loop system matrix and the state equation
of the closed loop system are then given by :
 1   2 -1/4 ρ-1/ 2 

 1 0    - ρ -
 A  - B K  =    -  M   2M ρ-1/ 4 ρ-1/ 2   =  M M 

 0 0   
  0    1 0 
 2 -1/4 ρ-1/ 2 

 x1   - ρ -   x1  +  0  w
 = M M    
x 2   1
x
0   2   
-1


 2 -1/4 ρ -1/ 2  Closed Loop
x
 1  - ρ -   x1  +  0  w
   M= M    
x 2   1 0   x 2  -1 Transfer

Functions
 The transfer functions for the closed loop system may
be obtained by taking the Laplace transform of each
term in the state equations.
2 1 / 4 r 1 / 2
sX1 s    r X1 s   X 2 s 
M M
sX 2 s   X1 s   Ws 

2 1 / 4 r 1 / 2 Closed Loop
sX1 s    r X1 s   X 2 s 
M M Transfer
sX 2 s   X1 s   Ws  Functions
 Eliminating X2(s) from the two equations, the transfer
function between the sprung mass velocity and the
road surface velocity is found.
r 1 / 2
X1 s  M

Ws  2 1 / 4 r 1 / 2
s2  r s
M M

r 1 / 2 Closed Loop
X1 s 
 M Transfer
Ws  2 1 / 4 r 1 / 2
s2 
M
r s
M Functions
 If one defines an undamped natural frequency for the
system as : -1/4
ρ
ωn = 1/2
M
 then the transfer function can be rewritten as :
X1 s  2
n

Ws  s 2  2n s  n2

X1 s  2
n Closed Loop

Ws  s 2  2n s  n2 Transfer
Functions
 Following the same procedure, the transfer function
between the suspension travel and the road surface
velocity can be obtained.
X2  s 
= 2
- s + 2ωn  
W  s  s + 2ωns + ωn2

X1 s 

n2 X 2 s 


 s  2 n  ρ-1/4
ωn = 1/2
Ws  s 2  2n s  n2 Ws  s 2  2n s  n2 M
Closed Loop Transfer Functions

 These two transfer functions indicate a damping ratio of
0.707 and a natural frequency depending on the choice
of r.
 So if one emphasizes handling, i.e., a stiff suspension
minimizing rattlespace a small value of r should be
chosen, resulting in a high n.
 If one emphasizes ride comfort minimizing u and thus
sprung mass acceleration, a large value of r is to be
chosen with the benefit of a soft suspension and a low
value of n.
y y
Integrator Integrator
y
Accelerometer
y
Amplifier a M
ay u
+ y-z b(y-z) + ay +b(y-z)
b Actuator
- + Summer
Amplifier u z
z Displacement Transducer
 The implementation of the active vibration isolator,

requires sensors (accelerometer and displacement
transducer), a force actuator and integrators, summers
and amplifiers.
y y
Integrator Integrator
y
Accelerometer
y
Amplifier a M
ay u
+ y-z b(y-z) + ay +b(y-z) Actuator
b
- + Summer
Amplifier z
z u
Displacement Transducer
 The most critical difficulty with the implementation is

the measurement of the road surface displacement.
 Note that the control signal can be expressed as
follows :
u   2M r 1 / 4 x1  r 1 / 2 x 2  ax1  bx 2  ay  by  z 

 Write the equations of motion. Basic Procedure
 Select the state variables.
 Obtain the state equation, i.e., the matrices [A], [B],
and [W].
 Decide on the contents of the weighting matrices [Q], [
R], and [N].
 Solve the Algebraic Matrix Riccati Equation for the
matrix [P].
 Find the optimal feedback matrix [K].
 Find the optimal control input {u}.
 Determine closed loop state equation [A]-[B][K].
 Solve the state equation to see the optimal behavior.
Velocity Feedback
Active Damper
 In the previous section a full set of states are
fed back to obtain the optimum system.
 However, in practice, the measurement of
some of the states may be difficult or
impossible.
 Then, a control law involving only some of
the states, such as sprung mass velocity,
which are easily and accurately measured
can be used.

y Integrator
y
Accelerometer Velocity
y
-g Amplifier
M Feedback
u
Active
- gy
k Damper
Actuator u
 It is observed that the suspension system consists of a

linear spring in parallel with an active component - an
actuator.

y Integrator
y
Accelerometer Velocity
y
M Feedback
-g Amplifier
Active
- gy u
k Damper
Actuator u
 The equation of
z motion for the
system is :
My  u  k ( z  y ) My  ky  u  kz

My  ky  u  kz
Velocity Feedback
Active Damper
 State variables are chosen as : x1  y x2  y  z
 Then the state equation is given by :
 x 1  0  k   x1   1  0
  M     M u    w
x 2  1 0  x 2   0   1
 Note that w is the road surface velocity input !

 x 1  0  k   x1   1  0 Velocity
   M     M u    w
x 2  1 0  x 2   0   1 Feedback
Active Damper
 Now, instead of feeding back both states, let us
feedback the velocity of the sprung mass only.
u   gy   gx1
 The closed loop state equation is then given as :
 x 1    g k  x   0 
  1
  M M x     w
x 2   1 0   2   1

 x 1    g 
k  x   0 
1
Velocity
  M 
M x    w

x 2   1 0   2   1
Feedback
Active Damper
 Defining the undamped natural frequency as
k
n2 
M
 the transfer function between the sprung mass
velocity and the road surface velocity input can
be obtained as :
X1 s  n2

Ws  g
s  s  n2
2
M

X1 s  n2 Velocity Feedback

Ws  g
Active Damper
s  s  n2
2
M
 It is observed that this system is

effectively a system with an active
damper.
 The damping ratio is directly
proportional to the velocity feedback
coefficient g.
 Thus the desired damping ratio can
be obtained by varying the value of
this coefficient.

y
2 Active Quarter
M
Car Model
u  With the use of a quarter

ks car model, it is possible to
c u y extend the analysis of the
1 previous section and obtain
m an active suspension
optimized with respect to
kt y0
ride comfort, suspension
travel, and road holding.

y
2 Active Quarter
M
Car Model
u Considering the quarter car
ks 
model together with an active

c u y
1 element in the suspension, the
m equations of motion can be
written as :
kt y0
My 2  cy 2  k s y 2  cy 1  k s y 1  u
my1  cy 1  (k s  k t )y 1  cy 2  k s y 2  k t y 0  u

y
2
M Active Quarter Car
Model
u
ks
 Ride comfort is usually rated
c u y
1 according to sprung mass
m acceleration.
kt y0  Road holding is related to
normal tire force which is
expressed as :
Ft  k t (y 1  y 0 )

y
2 Active Quarter Car
M Model
 It is usual to select the state
u variables as sprung mass
ks velocity, suspension travel, and
c u y tire deflection.
1
m  The fourth state variable is
chosen somewhat arbitrarily as
kt y0 the unprung mass velocity.
x1 = y 2 x 3  y1  y 0
x 2  y 2  y1 x4 = y 1

w  y 0 Active Quarter Car
Model
 The equations of motion can now be written in
the state space form :
 c ks c   1 
 x 1     0 x  0
x   M M M  1  M  0
 2 1 0 0  x
1  2  0    
        u   w
x 3   0 0 0 1  x 3   0   1
x 4   c ks k
 t
c 1
  x 4     0 
 m m m m  m
 The input to the vehicle is the road surface

velocity input !

Active Quarter Car
Model
 A first order approximation for typical road surfaces is
to take this input as a white noise Gaussian process
with zero mean. In such a case, a steady state
performance index might be :
J  lim T  
1
T
   
E 0T r1y 22  r 2 ( y 2  y1 )2  r 3 ( y 1  y 0 )2  r 4u 2 dt
 where r1, r2, r3, and r4 are weighting factors to be

chosen to emphasize ride comfort, suspension travel,
road holding, or energy expenditure.

 c k c   1 
 x 1   M  s 0  x  0
M M  1  M 
x 
 2
 
 1

0 0  1  x 2   0 
  
0
 
 u   w
Active Quarter
x 3   0 1  x 3   0   1
x 4   c
0
ks
0
k
 t
c 1
   x 4     0  Car Model
 m m m m  m
If one solves for y2 from the state

2

equations:
2
 c k c 1 
y 22  x 12    x1  s x 2  x 4  u 
 M M M M 
It is clear that the acceleration of the
sprung mass is a function of the
products of u and state variables, xi.

 T
   T
J  0 x Q x  u R u dt       Active Quarter
Car Model
x2 terms u2 terms
 It is noted that the acceleration of the sprung

mass is a function of the products of u and the
state variables, xi.
 Since the general form of the quadratic
performance index does not include such
terms, it can be extended as shown below.
1  T  T T  Q  N x 
J  lim T   E 0 x u  T    dt 
T   N  R  u 

1  T  T T   Q  N  x  
J  limT E 0 x u      dt
T  
  
  N  T
 R  u  
Active Quarter Car Model
 The extended quadratic performance index can be

simplified by matrix operations.
1  T
T 
T

J  limT  E  0 x  Q x  u  R  u  2 x  N  u dt 
T T
 
x2 terms u2 terms xu terms

J = limT
1  T
T   T T

E  0 x  Q x + u  R  u + 2 x  N  u dt 
T


 The optimal gain matrix is then given by :
-1
 K  = -  R   N  T
+  B
T
 P 
 where [P] is the unique, real, symmetric, and
positive definite solution of the matrix Riccati
equation :
 -1
 P  A  -  B R   N T
 +  A  -  B  R   N  -
-1 -1 T T
-  P B R   B   P  +   Q -  N  R   N   =  0
-1 T -1 T

K  = -R -1
 T
 N +  B T
 P 
 The matrix Riccati equation takes

the simplified form after some
manipulations:
 P A  +  A  T
 P -  P B +  N  R  -1
  B T
 P +  N  T
 +  Q =  0


 K  = -  R -1  NT +  BT  P 
PA  AT P   P B N R 1 BT P  NT  Q  0
MatLab Command lqr
 MatLab command lqr can be used to
obtain the optimal gain matrix K.
[K, P, E] = lqr (A, B, Q, R, N)
where
E = eig (A-BTK)

Example
 A quarter car model with an active suspension

component is considered. Numerical values of the
relevant parameters for a typical passenger car are given
below.
M : sprung mass = 250 kg
m : unsprung mass = 35 kg
ks : suspension spring stiffness = 15 kN/m
cs : damper damping coefficient = 800 N/m/s
kt : tyre stiffness = 120 kN/m

Example
 The magnitudes of the frequency

response functions between the
a) vehicle body acceleration,
b) suspension travel,
c) tire deflection,
and the road surface roughness velocity
are first plotted for the passive (i.e. u = 0)
system.

Example
 Then the control law, u, which minimizes the quadratic
performance index below
1
   
J  lim T   E 0T r1y 22  r 2 ( y 2  y 1 )2  r 3 ( y1  y 0 )2  r4u 2 dt
T
is used to emphasize ride comfort, suspension travel,

and road holding in turn by choosing appropriate
weighting matrices (r4 = 1 for all cases).
 The frequency response functions of the closed loop
(active) system are then superimposed on the
corresponding functions for the passive system.

Example
 Finally, the control law is selected as the absolute
velocity feedback,
u = - kx1
and the plots of the frequency response functions
of the closed loop (active) system are again
superimposed on the previous plots for the passive
and the full state feedback systems.
The value of the coefficient of x1 obtained in the
optimal feedback matrix of the full state feedback is
used as the value of k.

 Note : Example
It is well known that
1. the sprung mass acceleration transfer function
has an invariant point near the wheel hop
frequency,
kt 120000
ω1 = = = 58.6 [rad / s] = 9.3 [Hz]
m 35
and similarly,
2. the suspension deflection transfer function
has an invariant point at the so called
rattlespace frequency.
kt 120000
ω2 = = = 20.5 [rad / s] = 3.3 [Hz]
M+m 250 + 35

Invariant
r1 =1x107 point
Magnitude (dB)

Invariant
Magnitude (dB) point
r1 =1x107

r1 =1x107
Magnitude (dB)

30
Invariant
20 r2 =1x109 point
10
Magnitude (dB)
-10
Passive
Velocity Feedback
-30
-40
-2 -1 0 1 2
10 10 10 10 10
Frequency (Hz)

SUSPENSION TRAVEL
0
Invariant
-10
point
-20
-30
Magnitude (dB)
-40
-50 r2 =1x109
-60
Passive
-70
Full State Feedback
Velocity Feedback
-80
-90
-100
-2 -1 0 1 2
10 10 10 10 10
Frequency (Hz)

TIRE DEFLECTION
-20
r2 =1x109
-30
-40
Magnitude (dB)
-50
Passive
Velocity Feedback
-70
-80
-2 -1 0 1 2
10 10 10 10 10
Frequency (Hz)

30
Invariant
point
20
10
Magnitude (dB)
r3 =1x1012
-10
-20 Passive
Full State Feedback
Velocity Feedback
-30
-40
-2 -1 0 1 2
10 10 10 10 10
Frequency (Hz)

SUSPENSION TRAVEL
0
Invariant
-10 point
-20
-30
Magnitude (dB)
-40
-50 r3 =1x1012
-60
Passive
Velocity Feedback
-80
-90
-100
-2 -1 0 1 2
10 10 10 10 10
Frequency (Hz)

TIRE DEFLECTION
-20
r3 =1x1012
-30
-40
Magnitude (dB)
-50
-60
Passive
Full State Feedback
-70 Velocity Feedback
-80
-2 -1 0 1 2
10 10 10 10 10
Frequency (Hz)

30
Invariant
20 [Q]11=1x108 point
10
Magnitude (dB)
-10
-20 Passive
Full State Feedback
Velocity Feedback
-30
-40
-2 -1 0 1 2
10 10 10 10 10
Frequency (Hz)

SUSPENSION TRAVEL
0
Invariant
-10
point
[Q]11=1x108
-20
-30
Magnitude (dB)
-40
-50
-60
-70 Passive
Full State Feedback
-80 Velocity Feedback
-90
-100
-2 -1 0 1 2
10 10 10 10 10
Frequency (Hz)

TIRE DEFLECTION
-20
[Q]11=1x108
-30
-40
Magnitude (dB)
-50
-60
Passive
Full State Feedback
Velocity Feedback
-70
-80
-2 -1 0 1 2
10 10 10 10 10
Frequency (Hz)

Despite the ever-increasing studies on State of the Art
active suspensions, with illustrations
of impressive improvements over the conventional suspension
systems, there have been only a few recent real-life applications in
modern production vehicles*.
The drawbacks of active suspensions are high cost, added weight,
added power requirement, and difficulty of maintenance which
require specialized tools and technicians. Thus the applications
have been made exclusively on expensive high-end vehicles.
* https://www.thedrive.com/news/20996/the-crazy-bose-magic-carpet-car-
suspension-system-is-finally-headed-for-production.
https://media.daimler.com/marsMediaSite/en/instance/ko/The-new-Mercedes-Benz-
GLE-under-the-microscope-E-ACTIVE-BODY-CONTROL-suspension-
system.xhtml?oid=41843952.
https://www.audi-technology-portal.de/en/chassis/suspension-control-systems/audi-
s8-predictive-active-suspension

y
m = 250 kg Example 7.3
 Obtain and plot the transfer
k = 16 kN/m u c = 900 N/m/s
functions of the passive
yo system for the
 Acceleration of the
sprung mass, and
 Suspension travel.
Repeat the same for the active system using the

weighing matrices
0 0  With Q 22 =0.5x10 5,
Q  = 0 Q  and R = 1*10-3 1.0x105, and 1.0x106 to
 22 
minimize the quadratic

J =  ( x   Q  x +  u  R u)dt
T T performance index:
0
y
m = 250 kg
Example 7.3
k = 16 kN/m u c = 900 N/m/s
yo
  c  y  y 0   k  y  y 0   u
my
c k c 1
x1  y x 1  
y   x1  x 2  y 0  u
m m m m
x2  y  y 0 x 2  y  y 0  x1  y 0
  c k 1 c
 1  
x   1  
x  
  m m     m  u   m  y 0
 x 2    x2   
 1 0  0  1

y
m = 250 kg
Example 7.3 k = 16 kN/m c = 900 N/m/s

u
yo
% m = sprung mass (body mass)

% c = suspension damping coefficient
% k = suspension spring stiffness
% u = control force
% V = vehicle speed [kph]
m=250;k=16000;c=900;

  c k 1 c
 1  
x   1  
x  
  m m     m  u   m  y 0 Example
 x 2    x2   
 1 0  0  1
% Quarter car vehicle model

% W is the disturbance distribution
w  y 0
% matrix in {dx/dt}=[A]{x}+[B]{u}+[W]{w}
A=[-c/M -k/M; 1 0]; y
B=[1/M 0]’; m = 250 kg
W=[c/M -1]';
% {y}=[C]{x}+[D]{u}
k = 16 kN/m u c = 900 N/m/s
C=eye(2);
yo
D=zeros(2,1);

 c k 1 c
 x 1       x     
  m
1
m   m   u   m  y 0 Example
 x 2    x2   
 1 0  0  1
y   Cx  D u
Output equation:  y1   1 0   x1  y   0 
       u
 y 2  0 1  x2  y  y 0  0
y
% Passive system m = 250 kg
message='passive system’
[wn,z]=damp(A); k = 16 kN/m c = 900 N/m/s
u
wn=wn/(2*pi) yo
Z

Example
% --------------------------------------------------
f=logspace(-1,2,100);
w=2*pi*f;
% Obtain the transfer functions
% for the passive system
%---------------------------------------------------
% This gives the transfer functions
% ydot/y0dot and (y-y0)/y0dot
[T1,p1]=bode(A,W,C,D,1,w);

Example
%------------------------------------------
% Select weighting matrices
%------------------------------------------
Q=[0 0;0 100000];
R=0.001;
%--------------------------------------------------
% Obtain optimum feedback matrix
%--------------------------------------------------
[K,P,E]=lqr(A,B,Q,R);

Example
Af=A-B*K;
message='optimal system’
[wn,z]=damp(Af);
wn=wn/(2*pi)
Z
%-------------------------------------------------------------------------
% Obtain the transfer functions for the active system
%-------------------------------------------------------------------------
% This gives the transfer functions
% ydot/y0dot and (y-y0)/y0dot
[Ta1,pa1]=bode(Af,W,C,D,1,w);

Example
j=sqrt(-1);
s=j*w;
% First find the transfer function ydoubledot/y0dot
% for the passive and active systems
for i=1:100
T2(i)=T1(i,1)*s(i);
Ta2(i)=Ta1(i,1)*s(i);
end
Tdb=20*log10(T2);
Tadb=20*log10(Ta2);

Example
figure(1)
% Plot transfer function ydoubledot/y0dot for
% the passive and active systems
semilogx(f,Tdb,'b:',f,Tadb,'r-’)
title(blue/dashed:Passive - red/full:Active’)
ylabel('Transfer function [ydoubledot/y0dot]')

Example
% Now find the transfer function (y-y0)/y0dot

for i=1:100
T3(i)=T1(i,2);
Ta3(i)=Ta1(i,2);
end
Tdb=20*log10(T3);
Tadb=20*log10(Ta3);

Example
figure(2)
% Plot transfer function (y-y0)/y0dot
semilogx(f,Tdb,'b:',f,Tadb,'r-’)
title('blue/dashed:Passive - red/full:Active’)
ylabel('Transfer function [(y-y0)/y0dot]')

%---------------------------------------------- Example
% What we want for the rest here
% is ydoubledot/y0 !
for i=1:100
T4(i)=T1(i,1)*s(i)^2; % Passive
Ta4(i)=Ta1(i,1)*s(i)^2; % Active
end
%-----------------------------------------------------------------
% ISO reduced comfort boundary for 8 hours
fp=[1,4,8,80];
iso=[0.178,0.011,0.006,0.056];
isodb=log10(iso);

Example
% Road Displacement Power Spectral Density

% [m^2/cyc/m] is given by the expression :
% Syo(n)=S(no)*(n/no)^(-w)
Sno=1.4e-6; % [m^2/cyc/m]
Sfo=Sno/V; % convert to [m^2/cyc/s]
no=0.5/pi; % [cyc/m]
fo=no*V; % convert to [cyc/s]
ww=2.1;

Example
%--------------------------------------------------------
for i=1:100
% PSD of road displacement in [m^2/Hz]
PSDd(i)=Sfo*(f(i)/fo)^(-ww);
% PSD of road velocity in [(m/s)^2/Hz]
PSDv(i)=(2*pi*f(i))^2*PSDd(i);
end

Example
%--------------------------------------------------------------------------
figure(3)
% Plot Road Power Spectral Density
% Displacement and velocity
loglog(f,PSDd,'b:',f,PSDv,'r-’)
title('blue/dashed:Displacement - red/full:Velocity’)
ylabel('PSD [(m^2)/Hz, (m/s^2)^2/Hz]')

Example
%--------------------------------------------------------------
for i=1:100
% Passive System
PSDa(i)=log10(abs(T4(i))^2*PSDd(i));
% Active System
nPSDa(i)=log10(abs(Ta4(i))^2*PSDd(i));
end

Example
figure(4)
% Plot sprung mass acceleration PSD for
% passive and active systems
% Superimpose ISO 2631 8h limits
semilogx(f,PSDa,'b:',f,nPSDa,'r-',fp,isodb,'g--’)
title('blue/dashed:Passive -
red/full:Active')xlabel('Frequency [Hz]’)
ylabel('Acc. PSD [(m/s^2)^2/Hz]’)
text(2.75,-1.3,'Green line : ISO 2631 8h')

Example
The expected results are given below.
You can follow the initial improvement and the final
degradation as the value of Q22 is increased from 0.5x105
to 1.0x106.
Q22=0.5x105
blue/dashed:Passive - red/full:Active blue/dashed:Passive - red/full:Active
30 -10
-15
Transfer function [ydoubledot/y0dot]
25
Transfer function [(y-y0)/y0dot]

-20
20 -25
-30
15
-35
10 -40
-45
5
-50
0
-55
-5 -60
10-1 100 101 102 10 -1 10 0 10 1 10 2
Frequency [Hz] Frequency [Hz]

blue/dashed:Passive - red/full:Active
Example
30
25
20
15
10
Q22=1.0x105
5
-5
10 -1 10 0 10 1 10 2
Frequency [Hz] blue/dashed:Passive - red/full:Active
-10
-15

-20
-25
-30
-35
-40
-45
-50
-55
-60
10 -1 10 0 10 1 102
Frequency [Hz]

blue/dashed:Passive - red/full:Active

30
25
Example
20
15
10 Q22=1.0x106
5
-5 blue/dashed:Passive - red/full:Active
-10
10-1 100 101 102
Frequency [Hz] -15

-20
-25
-30
-35
-40
-45
-50
-55
-60
10-1 100 101 102
Frequency [Hz]

ME513VehicleDynamics LectureSlides 2021

Uploaded by

Copyright:

Available Formats

ME513VehicleDynamics LectureSlides 2021

Uploaded by

Document Information

Original Title

Copyright

Available Formats

Share this document

Share or Embed Document

Sharing Options

Did you find this document useful?

Is this content inappropriate?

Copyright:

Available Formats

ME513VehicleDynamics LectureSlides 2021

Uploaded by

Copyright:

Available Formats

What are the three main areas vehicle dynamics is subdivided into?

What are the three main areas vehicle dynamics is subdivided into?

What does vehicle handling refer to?

What does vehicle handling refer to?

ME 513

Prof. Dr. Y. Samim ÜNLÜSOY

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 4

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 2

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 3

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 4

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 5

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 6

Subject of  accelerate, decelerate, or travel at a

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 7

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 8

 Ride (or Ride Comfort) is related to

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 9

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 10

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 11

Prof. Dr. Y. Samim ÜNLÜSOY

Prof. Dr. Y. Samim Ünlüsoy ME 465 Automotive Engineering 2

 Vehicle handling is an area of vehicle

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 3

Steering wheel Position

Accelerator pedal Velocity

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 4

Steering wheel Position

Accelerator pedal Velocity

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 5

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 6

2  U : forward speed of the vehicle

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 7

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 8

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 9

 During cornering, centrifugal force acts at the center of

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 10

Lateral acceleration at A static Stability Safety Factor (SSF)

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 12

Estimate lateral acceleration at which the vehicle

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 13

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 14

 Realistic values considering suspension and tire

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 15

Prof. Dr. Y. Samim Ünlüsoy ME 513 Vehicle Dynamics 16

Prof. Dr. Y. Samim Ünlüsoy ME 465 Automotive Engineering 17

 Thus, NSH can be calculated if

235 / 55 R 17 80 W Speed Rating [kph]

Prof. Dr. Y. Samim Ünlüsoy ME 465 Automotive Engineering 19

Prof. Dr. Y. Samim Ünlüsoy ME 465 Automotive Engineering 20

Prof. Dr. Y. Samim Ünlüsoy ME 465 Automotive Engineering 21

 When the direction of the cords

Prof. Dr. Y. Samim Ünlüsoy ME 465 Automotive Engineering 22

 If, on the other hand, the cords

Prof. Dr. Y. Samim Ünlüsoy ME 465 Automotive Engineering 23

Prof. Dr. Y. Samim Ünlüsoy ME 465 Automotive Engineering 24

Prof. Dr. Y. Samim Ünlüsoy ME 465 Automotive Engineering 25

Prof. Dr. Y. Samim Ünlüsoy ME 465 Automotive Engineering 26

 Directional stability and

Prof. Dr. Y. Samim Ünlüsoy ME 465 Automotive Engineering 27