DYNAMIC MODELING OF FOUR-WHEEL FULL-SUSPENSION VEHICLE

DYNAMIC MODELING OF FOUR-WHEEL FULL-SUSPENSION VEHICLE S. Murat Ye¸silo˘ glu Hakan Temelta¸s [email protected] [email protected]

Istanbul Technical University, Electrical Engineering Department 34469, Istanbul TURKEY

Abstract: This paper presents an analytical dynamic model of a four-wheel fullsuspension vehicle using 24 independent dof. The methodology presented herein is capable of simulating vehicle dynamics on any terrain so long as it is smooth enough that tires can remain in contact. Given the steering torque, driving torque and terrain structure, the forward dynamics problem is solved to obtain velocities, accelerations, forces and torques of all joints. Among these, contact forces between tires and the road play a crucial role in designing an early sensory feedback system c of imminent vehicle roll-over. Copyright 2006 IFAC Keywords: Order of N dynamics, cooperating manipulators, vehicle dynamics

1. INTRODUCTION Mathematical modeling of motor vehicles has attracted a special interest of researchers due to increasing demand from the automotive industry. Recent improvements in computing power has led to a new phase in which many automotive simulation programs and algorithms were developed. Sayers and Han (1995) drove an 18-dof model using a symbolic equation generator; however, the terrain is limited to flat surfaces. A simplified model has also been introduced by Albagul and Wahyudi (2004). Works in Holmberg and Oussama (2000) and Bhatt et al. (2005) are examples of cooperating manipulator approach in mobile manipulators. The main difference between our work and theirs is in the detail of vehicle modeling. In many applications it is critical to know the tire-road interaction forces which determine how a vehicle turns, brakes and accelerates. It is desired that the methodology providing this information gives us an in-depth insight while generating closed form

equations of dynamics. The algorithm presented herein achieves these goals using Spatial Operator Algebra for multibody system dynamics originally developed by Rodriguez et al. (1992).

2. THE MODEL The proposed model of vehicle will be explained for separate parts first of which is the tire model. Tire characteristics are known to be highly nonlinear and very complicated. Since the motivation here is to demonstrate the use of the algorithm, the point contact tire model Captain et al. (1979) is employed. The algorithm, on the other hand, can be extended to include more complicated tire models. The joint at the point contact has 3 rotational dof. Pneumatic characteristics of tire are represented by spring and damper pairs in both vertical and horizontal directions. The end effector shown in Fig.1 corresponds to the center of the wheel.

end effector

3. ALGORITHM

representation of pneumatic characteristic ball joint =

terrain

3.1 Background Consider a rigid N-link manipulator number i as illustrated in Fig. 4 and 5. i

ENVIRONMENT (LINK N+1)

11111 00000 00000 11111

hN

i

i

O

hk

i

i

N+1

h k-1 i

ON

LINK N

h2

i

Ok

Fig. 1. Tire model

i

i

O k-1

LINK k-1

i

h1

O2 i

O1

i

O0

LINK 1

BASE (LINK 0)

Fig. 4. N-link serial manipulator i

h

k+1

i

θk+1

(a)

(b)

iτ = k+1 i

i

O

Fig. 2. a)Suspension b)Trapezoidal geometry

k+1

f

k+1

i

h

k i

θk

i

As shown in Fig.2a, the suspension model has 5dof. All of the rotational joints are represented with a cylinder whose axis is aligned with the axis of rotation. A small circle in the cylinder indicates that there is a torsional spring and a damper attached to that joint. Figure 2b shows the trapezoidal geometry of the steering mechanism. This design closely approximates the Ackerman condition. Finally, in Fig.3, the full model of a vehicle is given. The total number of dof is 51; however, as shown in section 3.5 only 24 of them are independent.

k,k+1 iτ

i i

i

CM

k

Ok

k

i

f

k

k,c

LINK k

Fig. 5. Vectors related with link k Let i ωk and i vk denote angular and linear velocities at i Ok respectively. Similarly, let i tk and i fk denote torque and force respectively. Spatial velocity and spatial force are then defined as follows: i  i  ωk t i i Vk = i Fk = i k vk fk

arm 7

arm 6

e c

d arm 5

arm 3

arm 8

arm 4

arm 9

11111 00000 00000 11111

11111 00000 00000 11111

f common object (axle) arm 2

arm 1

11111111111111111111 00000000000000000000 00000000000000000000 11111111111111111111

Fig. 3. Full model

Spatial velocity and spatial acceleration propagation from link k to link k+1 of manipulator i is done from base through tip by   I 0 i φk+1,k = i b ℓk,k+1 I where ℓb ≡ ℓ× is a skew symmetric linear operator. Spatial force propagation, on the other hand, is done from tip through base by i φTk+1,k . 3.2 Kinematic Constraints Let Vb+ and Vb− represent the base spatial velocities of all the arms on the arm side, and on the base side respectively. Nonslipage condition

requires both velocities to be the same. Let them be equal to Vb .     V b− I = Vb (1) V b+ I Let V be the stacked up spatial velocities of contact points where the bases and the tips of arms meet. Let V s consist of the time derivatives of each of the four wheels’ vertical displacements due to the surface geometry of the terrain. Since these displacements are known to be vertical, each elements of V s , i vs , is scalar.

Vt = σ t V = J θ˙ + φt,b Vb

(6)

From (3):

c

 T V c = VcT VdT VeT VfT  T V s = 1 vs 2 vs 3 vs 4 vs

−

= J θ˙ + φt,b Abs0 V s + φt,b A0bc V c

V c = A†tc Vt

A†tc = (ATtc Atc )−1 ATtc

(7)

Substituting (7) in (6) Vt = Lt J θ˙ + Ls V s

(8)

 −1 where Lt = I − φt,b A0bc A†tc , Ls = Lt φt,b Abs0 . Let us take time derivative of (8) αt = Lt J θ¨ + Lt J˙ θ˙ + L˙ t J θ˙ + Ls αs + L˙ s V s (9)

Base velocity part of the kinematic constraints can be written as:     Vs Vb− = Abs0 A0bc (2) Vc     I 0 0 0 i 000 I 0 0 0  0 i 0 0    Abc =  0 I 0 0  Abs =    0 0 i 0  0 0 0 φ1,f  0 00 i 0 0 0 φ2,f  T i= 0000 01  T  T Abs0 = ATbs 0 A0bc = 0 ATbc

Using the known equality, J˙ θ˙ = σt φa + at , finally we get αt = Lt J θ¨ + L˙ t J θ˙ + Ls αs + L˙ s V s + Lt (σt φa + at )

(10)

Now let us focus on the base accelerations. Taking the time derivative of (2) αb− = Abs0 αs + A0bc αc + A˙ 0bc V c

(11)

Substituting (8) in (7) V c = A†tc (Lt J θ˙ + Ls V s )

Tip velocity part of the kinematic constraints are:

(12)

Taking time derivative of (7) and substituting (8) 

Vt = Atc V

0 0  0 0 Atc =   0 0 φT1,f φT2,f

I 0 0 0

0 I 0 0

c

0 I 0 0

0 0 I 0

0 0 I 0

Putting (2) and (3) together,     s  V b− Ab V = Vt At Vc

0 0 I 0

T

(3)

0 0  I 0

(13)

Substituting (12) and (13) in (11) ˙ αb− = A0bc A†tc αt + Lu Lt J θ+ Abs0 αs + Lu Ls V s (4)

    where Ab = Abs0 A0bc , At = 0 Atc . Kinematic constraints are fully obtained from (1) and (4).

3.3 Velocity and Acceleration Analysis for Bases and Tips Spatial velocities of the joints with respect to an inertial frame can be obtained by using velocity propagation. V = φH θ˙ + φσb Vb−

αc = A†tc αt + A˙ †tc (Lt J θ˙ + Ls V s )

(5)

Furthermore, tip velocities of arms can be calculated.

(14)

where Lu = A0bc A˙ †tc + A˙ 0bc A†tc . Plugging (10) in (14), now we can obtain the equation for base ¨ θ, ˙ αs and V s . accelerations in terms of θ, αb− = La J θ¨ + Lb J θ˙ + Lc αs + Ld V s + Le (15) La Lb Lc Ld Le

= = = = =

A0bc A†tc Lt Lu Lt + A0bc A†tc L˙ t Abs0 + A0bc A†tc Ls Lu Ls + A0bc A†tc L˙ s A0bc A†tc Lt (σt φa + at ) = La (σt φa + at )

3.4 Pseudo Joints The Jacobian is a linear operator transforming joint space to a 6-dimensional task space. Under the assumption of no singularity, the Jacobian can be full rank in R6 only if the number of joints is

at least 6. Rank deficiency in the Jacobian means that there are certain directions which the tip of the arm cannot go. Finding these directions, one can reduce the size of the Jacobian to ensure that it will be full rank. A common approach is to use SVD but it may introduce discontinuity. i

h

−1 θ˙dep = E dep (−E ind θ˙ind + V r )

(17)

If E dep is not full rank, it means that the choice of S a is wrong. Similar equation to (17) can be written for accelerations too. In order to do that, let us first take the time derivative of (3). αt = Atc αc + A˙ tc V c

k+1

(18)

i

h

iτ = k+1

k

τpseudo

From(10) and (18)

i

θk

Atc αc + A˙ tc V c = Lt J θ¨ + L˙ t J θ˙ + Ls αs + L˙ s Vs + Lt σt φa + Lt at (19)

Pseudo Jo

int i O

iτ

k+1

k

τpseudo

g Using (12) in (19), and pre-multiplying it with A tc

i

Ok

key

Link k

† r ¨ g ˙ ˙ g ˙ A tc Lt J θ = Atc (Atc Atc Lt − Lt )J θ + α

(20)

where † g ˙ ˙ αr = −A tc (Ls αs + (Ls − Atc Atc Ls )Vs + Lt σt φa + Lt at )

Let us apply the operator S a to separate the dependent and independent variables.

keyway

† a−1 a ¨ g g ˙ A S θ = A tc Lt J S tc (Atc Atc Lt − −1 L˙ t )J S a S a θ˙ + αr

Fig. 6. Pseudo joint As an easy to implement and computationally efficient alternative approach, we can assume as if there were extra joints. Later we have to calculate the torques to keep those joints at zero angle at all times. These torques in reality correspond to internal torques where the joints were assumed. Fig.6 displays the idea.

Assigning   † a−1 g ˙ ˙ A = C dep C ind tc (Atc Atc Lt − Lt )J S Finally we get −1 θ¨dep = E dep (−E ind θ¨ind + C r θ˙ind + αr +

C dep E dep V r ) −1

3.5 Dependent Joint Variables

s g V r = −A tc Ls V

where C r = C ind − C dep E dep E ind . −1

˙ Now we will investigate the feasible sets of θ. Combining (3) and (8), and pre-multiplying both g sides of the result by A tc , we get: r ˙ g A tc Lt J θ = V

(21)

(16)

To uniquely determine the dependent joint velocities in terms of the independent ones, the elements ˙ are reordered in a way that of the joint space, θ, these elements are grouped into independent and dependent subspaces using the matrix S a formed by rearranging the rows of an identity matrix. As the choice of these subspaces is not unique, one can determine his/her choice based on preference and the structure. a−1 a ˙ g A S θ =Vr tc Lt J S

With the following definitions, (17) is obtained.  dep   dep ind  θ˙ a−1 a˙ g A L J S = S θ = E E tc t θ˙ind

3.6 Dynamic Constraints Dynamic constraints have a dual relationship with the kinematic constraints. Dual of (1) is     F b− I I = Fb (22) F b+ where Fb is a spatial force whose torque component represents the frictional term that slows down rolling wheel. Dual of (4) is    s  T T  F b− F = (23) Ab At Ft Fc where 1

   fs 0  2 fs  0 s c f f    F =  3 fs  F = K (Mf αf + bf ) K =  0  4 I fs

3.7 Equation of Motion From the acceleration propagation we can obtain all the accelerations as: α = φ(H θ¨ + a + σb αb− )

(24)

Now we can substitute (15) in (24) α = φ((H + σb La J )θ¨ + σb Lb J θ˙ + σb Lc αs + σb Ld V s + σb Le + a)

(25)

General force equation is F = φT (M α + b) + φT σtT Ft

(26)

Base forces can be written as: F b+ =

σbT F

(27)

Now the equation of motion can be obtained as θ¨ = M−1 (T − C θ˙ − La αs − Lb V s − D+ BbT Fb − BtT Fta )

(35)

3.8 Tip Forces Taking the time derivative of (16) and substituting (35) in that, we obtain Fta = Ga T + Gb Fb + Gc θ˙ + Gd αs + Ge V s + Gf (36)

Using (27) and (26) Fb+ = σbT φT (M α + b) + φTt,b Ft

(28)

We can write two dynamic constraints using (22) ATbs0 Fb− = F s

M=(H + σb La J )T φT M φ(H + σb La J ) C=(H + σb La J )T φT M φσb Lb J La=(H + σb La J )T φT M φσb Lc Lb=(H + σb La J )T φT M φσb Ld D=(H + σb La J )T φT (M φ(σb Le + a) + b)+ (A†tc Lt J )T F c g Bb=−La J Bt = A tc Lt J

AT0bc Fb− + ATtc Ft = F c (29)

g using (29) and annihilator of Atc , A tc T

g Ft = A†tc (F c − AT0bc Fb− ) + A tc Fta T

(30)

Using (22), let us substitute (28) in (30). T Ft = A†tc F c

+ (φσb A0bc A†tc )T (M α + b) + T † T g (φt,b A0bc A†tc )T Ft + A tc Fta − (A0bc Atc ) Fb

Finally we get Ft = (φσb La )T (M α + b) + (A†tc Lt )T F c + T T g (A (31) tc Lt ) Fta − L Fb a

Let us substitute (31) in (26)

Ga=(Bt M−1 BtT )−1 Bt M−1 Gb=(Bt M−1 BtT )−1 Bt M−1 BbT ˙ −1 g ˙ g Gc=(Bt M−1 BtT )−1 (A C) tc Lt J + Atc Lt J − Bt M −1 T −1 g −1 Gd=(Bt M Bt ) (Atc Ls − Bt M La ) ˙ −1 g Ge=(Bt M−1 BtT )−1 (A Lb ) tc Ls − Bt M −1 T −1 g Gf =(Bt M Bt ) (Atc Lt (σt φa + at ) − Bt M−1 D) 3.9 Suspension Dynamics Substituting (36) in (35) θ¨ = Na T + Nb Fb + Nc θ˙ + Nd αs + Ne V s + Nf (37) Na=M−1 (I − BtT Ga ) Nb = M−1 (BbT − BtT Gb ) Nc=−M−1 (C + BtT Gc ) Nd = −M−1 (La + BtT Gd ) Ne=−M−1 (Lb + BtT Ge ) Nf = −M−1 (D + BtT Gf ) Assigning B(Fb , αs , Vs ) = Nb Fb + Nd αs + Ne Vs + Nf θ¨ = Na T + Nc θ˙ + B

Now, let us separate the joint variables subject to F = [φ(I + σb La σt φ)] (M α + b) + + dynamic effect of suspension mechanism from the others. To do that, we first reorder the equations. T T g (A (32) tc Lt σt φ) Fta − (La σt φ) Fb −1 −1 S b θ¨ = S b Na S b S b T + S b Nc S b S b θ˙ + S b B Next, substituting (25) in (32) we get the force equation in terms of the joint accelerations Here, S b is obtained by rearranging the rows of h  T T ¨ a 54 × 54 identity matrix. Separation is done as F=(I + σb La σt φ) φ M φ (H + σb La J )θ+    i follows: s s σb Lb J θ˙ + Lc α + Ld V + Le + a + b + (33)  susp   susp  T g (A†tc Lt σt φ)T F c − (La σt φ)T Fb + (A θ¨ θ˙ tc Lt σt φ) Fta S b θ¨ = ¨nosusp S b θ˙ = ˙nosusp θ θ   The torques on the links can be extracted as Na1 Na2 b b−1 T S Na S = T = H F . Hence, Na3 Na4   s s ¨ ˙ T = Mθ + C θ + La α + Lb V + D+ N c1 N c2 b b−1 (34) S Nc S = BbT Fb + BtT Fta N c3 N c4 T

(A†tc Lt σt φ)T F c

 S T = b

T susp T nosusp



 S B= b

B susp B nosusp



Ya=Pa4 − Pa3 Pa−1 Pa2 1 −1 Yc=Pc4 − Pa3 Pa1 Pc2

Suspension dynamics equation is

Let us write (42) in the form of first order ODE:

M susp θ¨susp = −d θ˙susp − k θsusp

(38)

where d and k are damper and spring constants respectively. 

Yb = Pb4 − Pa3 Pa−1 P b2 1 true −1 pseudo Yd = B − Pa3 Pa1 B

    susp  θ¨susp Na2 B nosusp = T + + Na4 B nosusp θ¨nosusp    susp  −Na1 d + Nc1 Nc2 θ˙ ˙θnosusp + −Na3 d + Nc3 Nc4    susp  −Na1 k 0 θ (39) −Na3 k 0 θnosusp

Let the following be defined as:     b−1 Na2 b−1 Na1 k 0 Pc = S Pa = S Sb Na4 Na3 k 0   −1 Na1 d + Nc1 Nc2 Pb = S b Sb Na3 d + Nc3 Nc4



    true  θ¨true Yb Yc θ˙ = + true I 0 θtrue θ˙   Ya T nosusp true + Yd 0

(43)

4. CONCLUSION This paper has introduced a 24 independent dof model for the simulation of a four-wheel fullsuspension vehicle and developed an algorithm to solve the forward dynamics problem providing the information of velocities, accelerations, forces and torques of all joints. Among these, contact forces and torques which give ability to implement control with a desired objective such as preventing roll-over can be obtained. Smooth irregular terrains have been represented by four independently mobile platforms.

Now, the equation of motion can be modified to include suspension dynamics. REFERENCES

θ¨ = Pa T nosusp + Pb θ˙ + Pc θ + B

(40)

3.10 Torques Acting on Pseudo Joints Let S c , and S d be sorting matrices. S c θ¨ = S c Pa S d S d T nosusp + S c Pb S c S c θ˙ + −1

−1

−1

S c Pc S c S c θ + S c B 

θ¨pseudo θ¨true  P b1 P b3  P c1 P c3

  Pa1 Pa2 T pseudo + Pa3 Pa4 T nosusp true   pseudo  P b2 θ˙ + P b4 θ˙true   pseudo   pseudo  P c2 θ B + (41) P c4 θtrue B true





=

By definition θ¨pseudo = 0, θ˙pseudo = 0, θpseudo = 0 T pseudo = −Pa−1 (Pa2 T 1

nosusp true

+ Pb2 θ˙true +

Pc2 θtrue + B pseudo )

3.11 Forward Dynamics Including pseudo joint torques in the equation of motion, we get: θ¨true = Ya T nosusp true + Yb θ˙true + Yc θ

true

+ Yd (42)

A. Albagul and Wahyudi. Dynamic modelling and adaptive traction control for mobile robots. International Journal of Advanced Robotic Systems, 1(3):149–154, 2004. Rajankumar Bhatt, Chin Pei Tang, Michel AbouSamah, and Krovi Venkat. A screw-theoretic analysis framework for payload manipulation by mobile manipulator collectives. In Proceedings of ASME International Mechanical Engineering Congress & Exposition, Orlanda, FL, 2005. K. M. Captain, A. B. Boghani, and D.N. Wormley. Analytical tire models for dynamic vehicle simulation. Vehicle System Dynamics., 8:1–32, 1979. Robert Holmberg and Khatib Oussama. Development and control of a holonomic mobile robot for mobile manipulation tasks. International Journal of Robotics Research, 19(11): 1066–1074, 2000. G. Rodriguez, A. Jain, and K. Kreutz. Spatial operator algebra for multibody systems dynamics. The Journal of the Astronautical Sciences, 40 (1):27–50, January-March 1992. M. Sayers and D. Han. A generic multibody vehicle for simulating handling and braking. In 1995 IAVSD Symposium, Ann Arbor, MI., 1995.