Multibody Syst Dyn (2006) 15:325–350 DOI 10.1007/s11044-006-9013-7

Fuzzy control for equilibrium and roll-angle tracking of an unmanned bicycle Chih-Keng Chen · Thanh-Son Dao

Received: 28 December 2004 / Accepted: 13 March 2006  C Springer Science + Business Media B.V. 2006

Abstract This study presents steady turning motion and roll-angle tracking controls for an unmanned bicycle. The equations of motion describing the dynamics of a bicycle are developed using Lagrange’s equations for quasi-coordinates. Pure rolling without slipping constraints between the ground and two wheels are also considered in this model. These constraints introduce four holonomic and four non-holonomic constraint equations to the model. For the developed bicycle dynamics, one PID and one fuzzy controller that create steering torque are derived to recover the balance of the bicycle from a near-fall state. Furthermore, another fuzzy controller is added for controlling the bicycle to a desired roll angle which leads to its steady circular motion. The bicycle can track a given roll angle while maintaining its balance. The effectiveness of the control schemes is proved by simulation results. Keywords Bicycle dynamics · Bicycle balancing · Non-holonomic constraint · Multibody · Fuzzy control

1. Introduction Balancing a bicycle by human control is possible for almost everyone, yet stabilizing an unmanned bicycle proves to be rather demanding. Due to the challenges in fully understanding their dynamics and stabilization, bicycles have been attracting considerable attention from a number of researchers in the fields of physics, automation and control. One of the first publications of the full non-linear and its linearized equations of motion for an uncontrolled bicycle was proposed by Whipple [18]. His work was thoroughly verified by Hand [11] and Weir [19]. Another single-track vehicle model for motorcycles was studied comprehensively by Sharp et al. [15, 16] on numerous problems regarding dynamics, stability and control.

C.-K. Chen () · T.-S. Dao Department of Mechanical and Automation Engineering, Da-Yeh University, 112 Shan-Jiau Rd., Changhua, Taiwan 515 ROC e-mail: [email protected] Springer

326

Multibody Syst Dyn (2006) 15:325–350

To understand the nature of the dynamics and stability of a bicycle, Jones [12] did a variety of investigations on his bicycle. He pointed out that, in order to balance a ridden bicycle, an enough centrifugal force could be generated to correct its fall by steering the fork into the direction of the fall. This theory was well formalized mathematically by Timoshenko et al. [17] and is certainly confirmed by our bicycle riding experience. Schwab et al. [14] developed linearized equations of motion for a bicycle as a benchmark. In their study, the results obtained by pencil-and-paper, the numerical multibody dynamics program SPACAR and the symbolic software AutoSim were compared. Control efforts for stabilizing unmanned bicycles have also been addressed in previous studies. Beznos et al. [5] modeled a bicycle with gyroscopes that enabled the vehicle to stabilize itself in an autonomous motion along a straight line as well as along a curve. In their study, the stabilization unit consisted of two coupled gyroscopes spinning in opposite directions. Han et al. [10] derived a simple kinematic and dynamic formulation of an unmanned electric bicycle. The controllability of the stabilization problem was also investigated and a control algorithm for self-stabilization of the vehicle with bounded wheel speed and steering angle using nonlinear control based on sliding patch and stuck phenomena was proposed. Yavin [20] dealt with the stabilization and control of a riderless bicycle by a pedaling torque, a directional torque and a rotor mounted on the crossbar that generated a tilting torque. Getz et al. [9] derived a controller using steering and rear-wheel torque to recover the balance of their bicycle from a near fall as well as to converge to a time-parameterized path in the ground plane. In another study [8], Getz applied internal equilibrium control to the problem of path-tracking with balance for the bicycle. From the internal dynamics of the bicycle, an internal equilibrium manifold, a submanifold of the state-space, was constructed. Recently, some researchers have even added a balance mass to their electric bicycle system [13] or utilized fuzzy and intelligent control [7] for the stabilization problem. Most previous studies have dealt with simplified mathematical bicycle models, which were subsequently used to implement simulations, analysis and experiments. However, due to their simplicity, the mathematical models were incapable of presenting all of the dynamic motions of the system in certain situations. The nonlinear model itself should be used for the precise controller design and simulation. On the other hand, in our opinion, what makes bicycles different from other two-wheeled vehicle systems (like motorcycles) is that as investigating a bicycle system, we can eliminate the weave and wobble modes [15, 16] from tyre-road interaction. The relatively light weight and low speed motion of bicycles validate our assumption of pure-rolling wheels. By such observations, the complicated tyre model can be replaced by non-slipping kinematic constraints within a reasonable operational range. With that in mind, we approach the problem by modeling the bicycle as system with eleven generalized coordinates in three-dimensional space using Lagrange’s equations for quasi-coordinates. The constraint conditions of the two wheels are also considered. The derived equations of motion are then used to develop controllers, including PID and fuzzy controllers. With appropriate control laws, our concern is not only to stabilize the vehicle but also to control it to track a desired roll angle to make a steady turning motion, which was also suggested by Getz [8] and Berriah et al. [4]. They used much simpler models and different approaches however. This paper is organized as follows: In Section 2, the equations of motion of the bicycle model are briefly described. Four holonomic and four non-holonomic constraints are also added to describe the non-slipping motions of the two wheels. Section 3 deals with the search for the circular steady motions of the bicycle. Control schemes combining PID and fuzzy controllers to regulate the unmanned bicycle to a desired steady motion or a desired rollangle are discussed in detail in Section 4. Simulation results are also shown and analyzed. In Section 5 are some concluding remarks. Springer

Multibody Syst Dyn (2006) 15:325–350

327

Fig. 1 Schematic of a bicycle

2. Dynamics of bicycle model In this section, the compact equations of motion with nine coordinate variables are developed to describe the dynamics of the bicycle. The equations discussed here are developed by Lagrange’s equations for quasi-coordinates. The pure rolling constraints between the wheels and the ground are also considered as parts of the model of the bicycle system. 2.1. Coordinate systems The schematic bicycle model is shown in Figure 1. Let the uppercase letters A, B, D and F represent the vehicle body, the rear wheel, the front wheel and the fork, respectively, while the lowercase ones a, b, d, and f are used to designate the center of mass of each part. Reference point c is between the saddle and the vehicle body; e is a point between the vehicle body Springer

328

Multibody Syst Dyn (2006) 15:325–350

and the front fork; o and s are the contact points between the ground and the rear and front wheels, respectively. There are three SAE-standard coordinate systems used in the bicycle model: (1) an inertial frame o (I, J, K) fixed on the ground, (2) a reference frame c (ic , jc , kc ) mounted on the model at point c, and (3) a frame e (ie , je , ke ) placed on the front fork at point e. The coordinate e is obtained by rotating about c a rake and a steering angle as shown in Figure 1. In this paper the dynamics of the bicycle is described by the motion of the reference point c. Six coordinates are used to signify the positions and orientations of point c. The other three coordinate variables are the rotating angles of the front and rear wheels and the steering angle of the front fork. According to the foregoing definitions, the generalized coordinates can be written as q=[X

Y

ψ

Z

φ

θ

δ

φr

φ f ]T ,

(1)

where (X, Y, Z) are the position parameters and the relative rotation from o to c is determined  T by the 3-2-1 Euler angles Φ = ψ φ θ which are the rotations in the order about K, J and L axes successively; δ is the steering angle; φ f and φr are the rotating angles of the front and rear wheels, respectively. The velocity vector is  u = vx

vy

vz

ωx

ωy

ωz

δ˙

φ˙r

φ˙ f

T

,

(2)

the components of which are quasi-velocities (generalized speeds). Let R and Rce be the rotation matrices from o to c and from c to e , respectively. These two matrices are as follows ⎡ ⎤ cos ψ cos φ sin ψ cos φ − sin φ ⎢ ⎥ cos ψ cos θ + sin ψ sin θ sin φ cos φ sin θ ⎦, R = ⎣ − sin ψ cos θ + cos ψ sin θ sin φ sin ψ sin θ + cos ψ sin φ cos θ

− cos ψ sin θ + sin ψ sin φ cos θ

cos θ cos φ (3)

and ⎡

cos δ

sin δ

0

⎤⎡

cos ε

0

− sin ε

⎤

⎢ ⎥⎢ ⎥ 1 0 ⎦ Rce = ⎣ − sin δ cos δ 0 ⎦ ⎣ 0 0 0 1 sin ε 0 cos ε ⎡ ⎤ cos δ cos ε sin δ − cos δ sin ε ⎢ ⎥ sin δ sin ε ⎦ , = ⎣ − sin δ cos ε cos δ sin ε

0

cos ε

where ε is the inclined angle of the front fork as shown in Figure 1. 2.2. Dynamics of bicycle model For simplicity, the position vector, the velocity and angular velocity of body M in frame n will be denoted by rnM , vnM and ω nM , respectively. One can write the position vector of point c in o as roc = [ X Y Z ]T . Let the position vectors of the centers of mass of the Springer

Multibody Syst Dyn (2006) 15:325–350

329

vehicle body, rear wheel and point e relative to point c as ρa , ρb , and ρe , respectively. That is, ρac = [xa 0 z a ]T , ρcb = [xb 0 z b ]T and ρce = [xe 0 z e ]T in c . In a similar fashion, the position vectors 0ρ f and ρd of the centers of mass of the fork and the front wheel relative to point e are expressed in e as ρef = [x f 0 z f ]T and ρed = [xd 0 z d ]T . The angular velocities of the vehicle body and the rear wheel are written in c as ω cA = [ωx ω y ωz ]T ω cB = ω cA + ω cB/A = [ωx ω y − φ˙r ωz ]T,

(4)

where ω cB/A = [ 0 −φ˙r 0 ]T is the angular velocity of the rear wheel relative to the vehicle body. The vehicle body’s angular velocities ωz , ω y and ωx are referred to as yaw, pitch and roll rates, respectively. The angular velocities are related to the time rates of the three Euler angles by the formula ˙ ω A = SΦ,

(5)

where ⎡

− sin φ

⎢ S = ⎣ cos φ sin θ cos θ cos φ

0 cos θ − sin θ

1

⎤

⎥ 0⎦. 0

The angular velocities of the fork and the front wheel can be expressed in e as ω eF = Rce ω cA + ω eF/A , ω eD = ω eF + ω eD/F = Rce ω cA + ω eF/A + ω eD/F ,

(6)

˙ T is the angular velocity of the fork relative to the vehicle body and where ω eF/A = [0 0 δ] e T ω D/F = [0 −φ˙ f 0] is the angular velocity of the front wheel relative to the fork. Denote the velocity vector of point c by vc , one has vcc = [vx v y vz ]T . The velocities of all parts of the vehicle are expressed as vac vcb vcd vcf

= vcc + ω cA × ρac , = vcc + ω cA × ρcb ,

e  T ω F × ρed , = vcc + ω cA × ρce + Rce

T ω eF × ρef . = vcc + ω cA × ρce + Rce

(7)

The kinetic energy of all parts can be written in terms of their centers of mass as 1 c T c 1 c T m a va va + ω A I A ω cA , 2 2 1 c T c 1 c T TB = m b vb vb + ω B I B ω cB , 2 2 1 T 1 T TF = m f vcf vcf + ω eF I F ω eF , 2 2 1 c T c 1 e T TD = m d vd vd + ω D I D ω eD , 2 2 TA =

(8)

Springer

330

Multibody Syst Dyn (2006) 15:325–350

where T A , TB , TF and TD are the kinetic energy of the bicycle body, the rear wheel, the fork and the front wheel, respectively; I A , I B , I F and I D are the inertia matrices of each corresponding body. The total kinetic energy is obtained by summing all the kinetic energy of all the parts, thereby resulting in T = T A + T B + TF + T D =

1 T u Ju, 2

(9)

where J is the inertia matrix of the system. The potential energy has the form V = mgh,

(10)

 V = −g3T m a rao + m b rob + m f rof + m d rod ,

(11)

or,

where g3 = [ 0 0 g ]T . The generalized velocities q˙ are related to the quasi-velocities u by u = Yq, ˙

or q˙ = Wu,

(12)

where Y and W are the 9 × 9 transform matrices defined by ⎡

and q˙ = [ X˙

R

0

0

⎢ Y=⎣0

S

⎥ 0⎦,

0

0

I

Y˙

ψ˙

Z˙

⎡

⎤

RT

⎢ W = Y−1 = ⎣ 0 0 φ˙

θ˙

δ˙

0

0

S−1

⎥ 0⎦,

0

φ˙r

⎤ (13)

I φ˙ f ]T .

Lagrange’s equations for quasi-coordinates [3] can be formulated as d dt



∂T ∂u

+

∂T ∂T ∂V T Δ− W+ W = Unc , ∂u ∂q ∂q

(14)

˙ and Unc = WT Qnc are the nonconservative forces. where ∂∂uT = uT J, dtd ( ∂∂uT ) = u˙ T J + uT J, ∂u ∂u ˙ )W = ( dtd ∂u − ∂u )W, where ∂u = [ ∂q , . . . , ∂q ]. The coefficient matrix Δ is Δ = (Y − ∂u ∂q ∂ q˙ ∂q ∂q 1 9 One can rewrite Equations (14) in the standard form of differential equations as ˙ − ΔT Ju + WT Ju˙ = −Ju

∂T ∂q



T − WT

∂V ∂q

T + Unc ,

(15)

or simply, Ju˙ = Q. Springer

(16)

Multibody Syst Dyn (2006) 15:325–350

331

Fig. 2 Schematic of rear wheel

2.3. Constraint conditions for wheels By considering that the front and rear wheels of the bicycle roll without slipping, eight constraint equations, including four holonomic and four non-holonomic ones, can be found. 2.3.1. Rear wheel Figure 2 shows the schematic of the rear wheel. Let Rr = [−rr sin(αr ) 0 rr cos(αr )]T be the position vector of contact point o’ relative to the center of mass b of the rear wheel. The position of the contact point o is ro = rb + Rr .

(17)

Express Equations (17) in o and note that the contact point o is on the ground. The K component of roo is zero. This leads to a holonomic constraint Z − sin φ(xb − rr sin αr ) + cos φ cos θ(z b + rr cos αr ) = 0,

(18)

where αr , which is used to define the contact point due to the orientation change, is the included angle between Rr and kc . The symbol rr denotes the radius of the rear wheel. Furthermore, denote the intersection vector between the ground and the rear wheel planes by pro = [x y 0]T , since it is on the I–J plane. By the observation that pro is perpendicular to jc , their dot product is zero. That is

 R pro · jc = x (cos ψ sin φ sin θ − sin ψ cos θ) + y (sin ψ sin φ sin θ + cos ψ cos θ) (19) = 0,

where jc = [ 0 1 0 ]T . Equation (19) gives ⎧ ⎫ ⎧ ⎫ ⎪ ⎨x ⎪ ⎬ ⎪ ⎨ − (sin ψ sin φ sin θ + cos ψ cos θ) ⎪ ⎬ cos ψ sin φ sin θ − sin ψ cos θ pro = y = . ⎪ ⎪ ⎩ ⎪ ⎭ ⎪ ⎩ ⎭ 0 0

(20) Springer

332

Multibody Syst Dyn (2006) 15:325–350

Similarly, pr and Rr are perpendicular. By equating pr · Rr (with pr given in Equation (20)) to zero yields another holonomic constraint including the extra parameter αr cos φ cos θ sin αr − sin φ cos αr = 0.

(21)

To find the non-holonomic constraint equations, we write the velocity of contact point o’ as vo = vb + vo /b , ˙ r . Vector vo can be expressed in c as where vb = vc + ω A × ρb and vo /b = R

c ˙r + ω cB × Rrc = [vo x vco = vcb + R rel

vo y

vo  z ] T .

(22)

Assume that the wheel rolls without slipping, that is vco = 0. In o , vo can be written as voo = RT vco = [vo x 

vo y 

vo z  ]T = 0.

(23)

where vo x  = vo x (cos ψ cos φ) + vo y (cos ψ sin φ sin θ − sin ψ cos θ ) +vo z (cos ψ sin φ cos θ + sin ψ sin θ) = 0,

(24)

vo y  = vo x (sin ψ cos φ) + vo y (sin ψ sin φ sin θ + cos ψ cos θ ) +vo z (sin ψ sin φ cos θ − cos ψ sin θ) = 0, vo z  = −vo x sin φ + vo y cos φ sin θ + vo z cos φ cos θ = 0.

(25) (26)

The constraint vo z  = 0 can be proved to be an integrable equation, which can be obtained by differentiating Equations (18) and (21) (see [6]). Therefore, Equation (26) is trivial since it is only a velocity form of the holonomic constraint. Thus, Equations (24) and (25) represent two non-holonomic constraints for the system. 2.3.2. Front wheel The constraint conditions of the front wheel (see Figure 3) can be obtained by following the same procedures as those used with the rear wheel. Vector R f = [ −r f sin(α f ) 0 r f cos(α f ) ]T is used to designate the position of the center of mass d of the wheel relative to the contact point s. The position of the contact point s is rs = rd + R f .

(27)

Writing Equation (27) in o and equating the K component of ros to zero gives a holonomic constraint Z − xe sin φ + z e cos φ cos θ + z d (− sin φ sin ε + cos φ cos θ cos ε) +xd (− sin φ cos δ cos ε + cos φ sin θ sin δ − cos φ cos θ cos δ sin ε) Springer

Multibody Syst Dyn (2006) 15:325–350

333

Fig. 3 Schematic of front wheel

+r f sin α f (− cos φ sin θ sin δ + cos φ cos θ cos δ sin ε + sin φ cos δ cos ε) +r f cos α f (− sin φ sin ε + cos φ cos θ cos ε) = 0,

(28)

where α f is the included angle between R f and ke , r f is the radius of the front wheel. Angle α f is used to define the front contact point due to the orientation change of the front wheel. The intersection vector between the ground and the front wheel planes, p f , can be obtained by equating p f je to zero, that gives ⎧  ⎫ cos δ (sin ψ sin φ sin θ + cos ψ cos θ) − sin ψ cos φ sin δ cos ε ⎪ ⎧ ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪− ⎪ ⎪ ⎪ ⎨x ⎪ ⎬ ⎪ ⎨ ⎬ + sin δ sin ε (sin ψ sin φ cos θ − cos ψ sin θ ) o p f = y = cos δ (cos ψ sin φ sin θ − sin ψ cos θ) − cos ψ cos φ sin δ cos ε . (29) ⎪ ⎪ ⎩ ⎪ ⎭ ⎪ ⎪ ⎪ ⎪ ⎪ + sin δ sin ε (cos ψ sin φ cos θ + sin ψ sin θ ) 0 ⎪ ⎪ ⎪ ⎪ ⎩ ⎭ 0 Another holonomic constraint can be found by equating p f R f to zero. That gives (cos φ sin θ sin δ − cos φ cos θ cos δ sin ε − sin φ cos δ cos ε) cos α f + (cos φ cos θ cos ε − sin φ sin ε) sin α f = 0.

(30)

To find the non-holonomic constraints of the front wheel, we express the velocity of the contact point s as vs = ve + vd/e + vs/d , ˙ f . Vector vs can be written in c where ve = vc + ω A × ρe , vd/e = ω F × ρd , and vs/d = R as

e 



e  T T ˙e T vcs = vcc + ω cA ρce + Rce ω F × ρed + Rce R f rel + Rce ω D × Ref = [vsx vsy vsz ]T. (31) Springer

334

Multibody Syst Dyn (2006) 15:325–350

With the assumption of non-slip rolling of the wheel and zero-velocity motion of the contact point s, one has vcs = 0. Furthermore, in o , Equation (31) becomes vos = RT vcs = [vsx  vsy  vsz  ]T , or vsx  = vsx cos ψ cos φ + vsy (cos ψ sin φ sin θ − sin ψ cos θ ) +vsz (cos ψ sin φ cos θ + sin ψ sin θ) = 0,

(32)

vsy  = vsx sin ψ cos φ + vsy (sin ψ sin φ sin θ + cos ψ cos θ) +vsz (sin ψ sin φ cos θ − cos ψ sin θ ) = 0, vsz  = −vsx sin φ + vsy cos φ sin θ + vsz cos φ cos θ = 0.

(33) (34)

It is also proved that Equation (34) can be obtained by differentiating Equations (28) and (30) (see [6]). Hence, Equations (32) and (33) give two non-holonomic constraints. 2.4. Equations of motion with constraints As a result, eight constraint equations are obtained from the no slipping conditions. Among these, there are four holonomic constraints (Equations (18), (21), (28) and (30)) and four non-holonomic constraints (Equations (24), (25), (32) and (33)). To derive the constraint equations, two algebraic variables, αr and α f , are introduced. Hence, the generalized coordinates and velocity vectors are expanded to qe = [X U = [vx

Y vy

Z vz

ψ φ θ δ ω x ω y ωz

φr φ f αr α f ]T , ˙δ φ˙r φ˙ f α˙ r α˙ f ]T .

(35)

Differentiating all the holonomic constraints to yield their velocity forms and then combining with the non-holonomic constraints leads to constraint equations as follows BU = 0,

(36)

where B is an 8 × 11 matrix, referred to as the constraint Jacobian matrix. The equations of motion with constraint conditions appear to be one set of DAE (differential-algebraic equation) as ⎧ ⎪ ⎨ q˙ e = We U ˙ = Qe + τ − BT λ , (37) Je U ⎪ ⎩ BU = 0 where ⎡

W

⎢ We = ⎣ 0 0 Springer

0

0

⎤

1

⎥ 0⎦

0

1

Multibody Syst Dyn (2006) 15:325–350

335

is an 11 × 11 matrix, ⎡

J

⎢ Je = ⎣ 0 0

0

0

⎤

0

⎥ 0⎦

0

0

is an 11 × 11 generalized mass matrix, ⎧ ⎫ ⎪ ⎨Q⎪ ⎬ Qe = 0 ⎪ ⎩ ⎪ ⎭ 0 is vector of applied forces, λ represents the eight Lagrange multipliers or constraint forces coupled to the system by the 8 × 11 constraint Jacobian matrix B, and τ is the generalized non-conservative force vector. Remark 2.1. In the above DAE (37), eleven generalized coordinates and eleven velocity variables are used for the system. The system is constrained by the eight constraints as shown in Equation (36) (four holonomic and four non-holonomic constraints) and by such there are only three degrees of freedom in the velocity space. Due to the four non-holonomic constraints, there are four kinematic coordinates in the system. In other words, the independent states of the system can be described by three generalized speeds together with seven generalized coordinates. Remark 2.2. By the differentiation of the velocity constraints in Equation (36) and the rearrangement of the terms in Equation (37), one can have the governing equations of motion in the form      ˙ Qe + τ Je B T U = . (38) ˙ B 0 λ −BU Numerical algorithms such as coordinate reduction and embedding methods [1] can be used to solve the DAE (38). Some numerical stabilization schemes such as Baumgarte method and post-stabilization technique [2] can help to ensure the convergence of the solution and increase the integration accuracy.

3. Steady circular motion of bicycle Steady motion or dynamic equilibrium is an essential and useful concept in dynamics. For a bicycle system, steady turning motion is necessary to understand its interesting dynamic behaviors and to control the vehicle to perform a given circular motion. To calculate the dynamic equilibriums of the bicycle system, one can utilize the constraint equations to compute the state variables when the system is at equilibrium. The dynamic equation in (37) indicates that ˙ = Qe + τ − BT λ. Je U

(39) Springer

336

Multibody Syst Dyn (2006) 15:325–350 Table 1 Simulation parameters Name

Value

Name

ma md ρa ρd ρf g

11.05(kg) 3.92(kg) (0.1296, 0,0.285) (0, 0,0.601) (0.017, 0,0.1083) 9.80665(m/s2 )

(a) mb 2.09 (kg) mf 4.04 (kg) ρb (−0.365, 0,0.503) ρe (−0.789, 0,0.078) r f = rr 0.325(m) ε 15◦ Moment of inertia (kg-m2 ) (b) ⎡ ⎤ 0.407 0 −0.068 ⎣ IA = 1.934 0 ⎦ 1.558 ⎡ ⎤ 0.421 0 −0.025 ⎢ ⎥ IF = ⎣ 0.384 0 ⎦

Vehicle body

Front fork ⎡ Wheels

Value

⎢ IB = ⎣

0.109

0

0

0.218

0

0.041 ⎡ 0.204 ⎥ ⎢ ⎦; I D = ⎣ ⎤

0.109

0

0

0.408

0

⎤ ⎥ ⎦

0.204

Practice shows that only the steering torque is needed to control the bicycle to an equilibrium point. Therefore, in the force vector τ , only the seventh element, the steering torque τ s , being nonzero is necessary; that is τ = [0, . . . , 0, τs , 0, . . . , 0]T . In order to eliminate the Lagrange multipliers, Equation (39) is multiplied by the matrix T which is a (n − m) × n orthogonal complement of matrix B, so that TBT = 0 (see [1]). Equation (39) thus becomes ˙ = TQe + Tτ . TJe U

(40)

Matrix T can locally reduce the constrained system (Equation (39)) into an unconstrained one. Therefore the number of equations is reduced to n − m = 11 − 8 = 3 which is the same as its number of degrees of freedom. At the steady motion situation, the right-hand side of the locally unconstrained system (Equation (40)) should become vanished, i.e. TQe + Tτ = 0, while the constraint equations are still satisfied. Note that the constraint equations are independent of the variables X, Y , ψ, φ r and φ f and at a steady turning state, δ˙ = α˙ r = α˙ f = θ˙ = 0. There are three equations in TQe + Tτ = 0 from which one can search for the speed vx and the steering torque τ s if the roll angle δ is given. This gives the equilibrium to maintain the bicycle at the steady turning motion. The search results are selectively shown in Table 2. Figure 4 shows the relationship between vx and the steering angle δ when the bicycle is at steady motion with a constant roll angle θ. This figure indicates that to maintain a constant θ angle, angle δ must increase when vx decreases. One can see that, at high speeds, the change in the slope of the curves is much milder than it is at low speeds. However, one should note that for larger roll angle, say θ ≥ 30◦ , the relation may not be practical since the no slipping assumption at the contact points could be violated and the contact point between the tyre and ground could be shifted to the sides of the tyres. Springer

Multibody Syst Dyn (2006) 15:325–350

337

Table 2 Selected equilibrium points θ (deg) States

10

15

20

25

X˙ (km/h) Y˙ (km/h) Z˙ (km/h)

NAC NAC 0 10.283 0.750 −0.167 0.581 0 0 0 −9.051 −9.390 0 0 NAC NAC −0.813 NAC −0.175 10 14.739 NAC NAC −0.178 11.787 −0.0931 4.9418

NAC NAC 0 9.552 1.168 −0.373 0.906 0 0 0 −8.772 −9.574 0 0 NAC NAC −0.796 NAC −0.329 15 22.542 NAC NAC −0.341 7.615 −0.0745 2.9514

NAC NAC 0 8.371 1.568 −0.657 1.291 0 0 0 −8.289 −10.148 0 0 NAC NAC −0.774 NAC −0.424 20 31.157 NAC NAC −0.452 1.450 −0.0201 1.8459

NAC NAC 0 6.611 2.065 −1.065 1.798 0 0 0 −7.608 −11.085 0 0 NAC NAC −0.747 NAC −0.562 25 41.731 NAC NAC −0.366 −7.355 −0.0010 1.0745

vx (km/h) v y (km/h) vz (km/h) ψ˙ (rad/s) φ˙ (rad/s) θ˙ (rad/s) δ˙ (rad/s) φ˙r (rad/s) φ˙ f (rad/s) α˙ r (rad/s) α f (rad/s) X (m) Y (m) Z (m) ψ (deg) φ (deg) θ (deg) δ (deg) φr (deg) φ f (deg) αr (deg) α f (deg) τs (N-m) Rc (m)

NAC stands for ‘Not a constant’ Rc is the radius of the circular trajectory of reference point C on the ground plane

The equilibrium for the unmanned bicycle system shows the region of the steady-state motion, which is useful for controller design. The following section shows how an unmanned bicycle can be balanced and controlled to a desired equilibrium state.

4. Controller design With the developed equations of motion, one can design a control strategy to balance the unmanned bicycle system for stable straight-running motion, or furthermore, to control the vehicle to a desired steady turning motion which is one of the dynamic equilibrium states. In this section, these issues are discussed and the corresponding control schemes are given. The fuzzy controllers with the properties that formulate human knowledge by rule bases from human expert demonstration and experiences are used. This study follows the human experiences in riding a bicycle and proposes a cascaded fuzzy structure as shown later to control the bicycle. Springer

338

Multibody Syst Dyn (2006) 15:325–350

Fig. 4 Relationship between vx and δ with constant θ angles

Fig. 5 Balancing control scheme

4.1. Balancing control For balancing a bicycle in a straight-running motion, a control scheme including fuzzy and PID controllers is developed. The block diagram of the closed-loop system is shown in Figure 5. The PID controller for the inner loop is τs (k) = k p eδ (k) + kd eδ (k) + ki

k 

eδ (i)

i=1

where k p , ki and kd are the three corresponding parameters for the proportional, integral and derivative parts of the PID controller, eδ = δref − δ denotes the error, eδ (k) = eδ (k) − eδ (k − 1) is the difference between the successive error signals. In this study, (k p , ki , kd ) are chosen as (1000, 1, 4000). The PID controller controls the bicycle by generating a torque τ s to steer the fork into the direction of fall. Springer

Multibody Syst Dyn (2006) 15:325–350

339

In the outer loop is a fuzzy controller whose inputs are the roll angle θ and its change θ . The value of the controller output is the reference steering angle δ ref . The simulation parameters are listed in Table 1. These numerical values are measured from an actual bicycle system. To balance an unmanned bicycle, one needs to adjust the reference angle δ ref according to θ and θ until both θ and θ become zero. With this in mind, one can design a fuzzy controller in the outer loop to feed appropriate δ ref into the inner loop to stabilize the system. Let θ and θ be the inputs of the fuzzy controller, as shown in Figure 5. The output variable of the fuzzy controller is δ ref . Next, linguistic quantification is used to specify a set of rules to describe the strategy for stabilizing the system. Assume that one is facing the heading direction of the bicycle; thus, θ is negative when the bicycle leans to the left and positive when it leans to the right. In considering the bicycle in the following three typical situations, the corresponding rules are discussed:

1. If θ is negative large (NL) and θ is NL, then δ ref is NL. This rule quantifies the situation wherein the bicycle is at a significant roll-angle to the left and it is moving away from the upright position. Hence, it is clear that one should steer the fork to the left at a large angle to counteract the movement so that the bicycle moves towards the upright position. 2. If θ is zero (Z) and θ is Z, then δ ref is Z. This rule quantifies the situation wherein the bicycle is very near the upright position. In this case, one need not steer the fork at all, because the bicycle is already in its proper position. 3. If θ is positive large (PL) and θ is PL, then δ ref is PL. In this case, the bicycle is at a significant roll-angle to the right and it is falling to the right (the bicycle is moving in the direction of positive θ). Therefore, it is clear that a large positive angle (steering to the right) should be applied to the fork to correct the fall.

Using the above approach, one could continue to formulate rules for all possible cases, as shown in Table 3. Notice that the body of the table lists the linguistic-numeric consequents of the rules, and the left column and top row of the table contain the linguistic-numeric premise terms. For this controller, with two inputs and seven linguistic values for each of these, there are at most 72 = 49 possible rules.

Table 3 Rule-base for balancing fuzzy controller θ

NL

NM

NS

Z

PS

PM

PL

θ NL NM NS Z PS PM PL

NL NL NL NM NM NS Z

NL NL NM NM NS Z PS

NL NM NM NS Z PS PM

NM NM NS Z PS PM PM

NM NS Z PS PM PM PL

NS Z PS PM PM PL PL

Z PS PM PM PL PL PL Springer

340

Multibody Syst Dyn (2006) 15:325–350

For all fuzzy controllers in this study, triangular membership functions are used. A triangular membership can be presented by the formula ⎧ 0 x c

(41)

By using Equation (41), the membership functions for input θ are defined as follows: θ θ θ θ

is NL: triangle (x : −2, −π/2, −0.7) and θ is PL: triangle (x : 0.7, π/2, 2); is NM: triangle (x : −0.9, −0.6, −0.3) and θ is PM: triangle (x : 0.3, 0.6, 0.9); is NS: triangle (x : −0.4, −0.2, 0) and θ is PS: triangle (x : 0, 0.2, 0.4); is Z: triangle (x : −0.15, 0, 0.15) .

Similarly, the membership functions for θ are defined as follows: θ θ θ θ

is NL: triangle (x : −∞, −10, −6.7) and θ is PL: triangle (x : 6.7, 10, +∞); is NM: triangle (x : −10, −6.7, −3.3) and θ is PM: triangle (x : 3.3, 6.7, 10); is NS: triangle (x : −5, −2.5, 0) and θ is PS: triangle (x : 0, 2.5, 5); is Z: triangle (x : −1, 0, 1) .

The membership functions for δ ref are as follows: δ ref is NL: triangle (x : −1.3, −1, −0.7) and δ ref is PL: triangle (x : 0.7, 1, 1.3); δ ref is NM: triangle (x : −0.9, −0.6, −0.3) and δ ref is PM: triangle (x : 0.3, 0.6, 0.9); δ ref is NS: triangle (x : −0.4, −0.2, 0) and δ ref is PS: triangle (x : 0, 0.2, 0.4); δref is Z: triangle (x : −0.15, 0, 0.15) . These membership functions are plotted in Figure 6. Various numerical tests are performed to evaluate the effectiveness of the designed fuzzy controller. Figure 7 shows the time histories of θ and δ with the same initial roll and steering angles, θ = 30 ◦ and δ = −15 ◦ , but different initial speeds, vx = 20 (km/h), vx = 30 (km/h) and vx = 50 (km/h). One can see that the controlled roll and steering angles asymptotically approach zero. Figure 7 indicates that, with a higher initial forward speed, the bicycle can move toward a straight-running status faster since a larger raising torque is produced to lift the bicycle to the upright position faster. Figure 8 shows the trajectories of the reference points of the bicycle on the ground plane in three simulation cases, wherein the paths asymptotically approach three different straight lines. Remark 4.1. In simulations, a constant initial forward speed is assumed. However, this forward speed is crucial in the bicycle balance control since it can affect the magnitudes of the raising torques generated. The raising torques, which are produced by the steering action, can be used to control the lean (or roll) angle. The energy for the rolling motion is converted into the system by draining the energy from the forward motion due to the change of steering angle. That is one of the interesting properties of non-holonomic systems. At a very low forward speed, this raising torque becomes inadequate to overcome the gravity and to lift the bicycle upward fast enough. Therefore, controlling the forward speed is critical for the bicycle balancing control. The steering torque alone can not do. Springer

Multibody Syst Dyn (2006) 15:325–350

341

Fig. 6 Membership functions for fuzzy system FISδ

Fig. 7 Time histories for roll and steering angles with different initial velocities Springer

342

Multibody Syst Dyn (2006) 15:325–350

Fig. 8 Bicycle trajectories with different initial velocities

Fig. 9 Control scheme with additional angle δ a

4.2. Equilibrium and roll-angle tracking controls An interesting result occurs when an appropriate additional angle δ a is added to the output of the previous fuzzy controller, as shown in Figure 9. The bicycle immediately falls to its equilibrium state for its corresponding speed and thus leads to a circular motion. Figure 10 shows a simulation case in which a constant angle δ a = 10◦ is applied. One can see that the roll-angle approaches a significantly large constant value other than zero. The value of the resulting θ angle is dependent on the value of the additional angle δ a . The circular trajectory of the bicycle in this case is shown in Figure 10b. With this result, it appears that the bicycle can be controlled to an equilibrium point by applying an appropriate additional angle δ a to the output of the fuzzy controller FISδ . This result also suggests that one can design another controller to track a desired roll-angle or to control the bicycle to a desired steady motion. Springer

Multibody Syst Dyn (2006) 15:325–350

343

Fig. 10 Resulting roll-angle (a) and circular trajectory of reference point c on ground plane (b) with δ a = 10◦

Fig. 11 Equilibrium control scheme

Springer

Multibody Syst Dyn (2006) 15:325–350

Fig. 12 Membership functions for FISδ a

344

Springer

Multibody Syst Dyn (2006) 15:325–350

345

Fig. 13 An equilibrium point with θ = 30◦ (a) and circular trajectory of reference point c on ground plane (b)

The control scheme following this idea is shown in Figure 11. Another control loop is added to the existing one. This control loop uses a fuzzy controller (FISδa ), whose inputs are eθ (k) = θref (k) − θ (k), the difference between the desired and the actual roll-angles, and its change, eθ = eθ (k) − eθ (k − 1), to generate appropriate δ a values. Linguistic quantification used to specify a set of rules for this controller is characterized by the following three typical situations: Springer

346

Multibody Syst Dyn (2006) 15:325–350

Fig. 14 Roll-angle time history (a) and resulting trajectory of reference point c (b) on ground plane

1. If eθ is NL and eθ is NL, then δ a is PL. This rule quantifies the situation wherein the desired roll-angle is much smaller than the actual one and the bicycle is falling to the right at a significant rate. Hence, one should steer the fork to the right at a large positive angle to make the bicycle lean to the left. 2. If eθ is Z and eθ is Z, then δ a is Z. This rule quantifies the situation wherein the bicycle is already in its proper position. No control effort is needed. 3. If eθ is PL and eθ is PL, then δ a is NL. This rule quantifies the situation wherein the desired roll-angle is much larger than the actual one and the bicycle is falling to the left at a significant rate. Therefore, one should steer the fork to the left at a large angle to make the bicycle lean to the right. In a similar fashion, the complete rule-base is constructed as listed in Table 4. Notice that the body of the table lists the linguistic-numeric consequents of the rules, and the left column Springer

Multibody Syst Dyn (2006) 15:325–350

347

Fig. 15 Ground contact paths of two wheels for tracking a varying θ angle

and top row of the table contain the linguistic-numeric premise terms. For this controller, with two inputs and seven linguistic values for each of these, there are at most 72 = 49 possible rules. The membership functions for input eθ are defined in the form of Equation (41) as follows: eθ eθ eθ eθ

is NL: triangle (x : −2, −1.1, −0.35) and eθ is PL: triangle (x : 0.35, 1.1, 2); is NM: triangle (x : −0.9, −0.5, −0.1) and eθ is PM: triangle (x : 0.1, 0.5, 0.9); is NS: triangle (x : −0.2, −0.07, 0) and eθ is PS: triangle (x : 0, 0.07, 0.2); is Z: triangle (x : −0.01, 0, 0.01).

Similarly, the membership functions for input eθ are defined as follows: eθ is NL: triangle (x : −∞, −10, −6.7) and eθ is PL: triangle (x : 6.7, 10, +∞) ; eθ is NM: triangle (x : −9, −6, −2) and eθ is PM: triangle (x : 2, 6, 9) ; Table 4 Rule-base for equilibrium fuzzy controller eθ

NL

NM

NS

Z

PS

PM

PL

eθ NL NM NS Z PS PM PL

PL PL PL PM PM PS Z

PL PL PM PM PS Z NS

PL PM PM PS Z NS NM

PM PM PS Z NS NM NM

PM PS Z NS NM NM NL

PS Z NS NM NM NL NL

Z NS NM NM NL NL NL Springer

348

Multibody Syst Dyn (2006) 15:325–350

Fig. 16 Roll-angle time response (a), tracking error (b) and resulting trajectory of reference point c on ground plane (c) for tracking θ (t) = 20◦ sin (0.5π t) Springer

Multibody Syst Dyn (2006) 15:325–350

349

eθ is NS: triangle (x : −4, −2, 0) and eθ is PS: triangle (x : 0, 2, 4) ; eθ is Z: triangle (x : −0.5, 0, 0.5) . The membership functions for output δ a are as follows: δa δa δa δa

is NL: triangle (x : −5.7, −3, −1.5) and δ a is PL: triangle (x : 1.5, 3, 5.7); is NM: triangle (x : −3, −1.5, −0.6) and δ a is PM: triangle (x : 0.6, 1.5, 3); is NS: triangle (x : −1, −0.4, 0) and δ a is PS: triangle (x : 0, 0.4, 1); is Z: triangle (x : −0.03, 0, 0.03).

These membership functions are plotted in Figure 12. To reduce the steering vibration when the bicycle is already close to its proper position, the intervals of the membership functions of the control inputs and output around zero are chosen to be rather small. With this control strategy, one can control the unmanned bicycle to a desired equilibrium point. The effectiveness of this control scheme is demonstrated by a variety of simulation cases, one of which is shown in Figure 13. In this simulation, the bicycle is controlled to an equilibrium point with the roll-angle θ = 30◦ , with which the corresponding motion of the reference point c results in a circular path on the ground plane, as shown in Figure 13b. Apparently, with this control scheme, the bicycle can be controlled to track a time-varying roll-angle other than a constant one. Figure 14 shows a simulation in which the bicycle is controlled to follow a varying roll-angle. As shown in Figure 14a, the desired roll-angle is initially zero, after 1 second the bicycle is controlled to a roll-angle θ = 20 ◦ and is subsequently switched to θ = −20◦ at 5 sec. The bicycle turns to the upright position after 9 sec. The resulting trajectory of the reference point c on the ground plane is shown in Figure 14b. The paths of the contact points by the two wheels are shown in Figure 15. Another simulation, in which the bicycle is controlled to follow a desired sinusoidal rollangle θ (t) = 20o sin (0.5π t) , is shown in Figure 16. Figure 16b shows the tracking error. The resulting trajectory of the reference point c on the ground plane is shown in Figure 16c. The result indicates the effectiveness of the proposed control scheme for controlling the unmanned bicycle to follow a time-varying roll-angle.

5. Conclusion In this study, a dynamic model of an unmanned bicycle has been developed by using Lagrange’s equations for quasi-coordinates. By considering the contact relationship between the two wheels and the ground plane, the constraint conditions, including four holonomic and four non-holonomic ones, have been derived. According to this mathematical model, the fuzzy and PID controllers have been designed to stabilize the bicycle in its straight-running motion. The simulation results indicate that the bicycle can be controlled to the upright position faster for higher initial speeds vx . This result is compatible with our bicycle-riding experience. By adding an additional angle δ a , the bicycle turns out to fall to a circular equilibrium point, which is shown in Figure 10. The results obtained by adding δ a inspired us to design another fuzzy controller to control the bicycle to a desired equilibrium point or to track a desired roll-angle. These have also been confirmed by numerical tests. Future work will continue with ground-path following problem. By considering the distance from the reference point c to the given path, the vehicle is controlled to follow the heading direction the path. Once this is done, obstacle avoidance, path-planning and associated problems will certainly be interesting problems. Springer

350

Multibody Syst Dyn (2006) 15:325–350

Acknowledgements The authors would like to thank the National Science Council of the Republic of China, Taiwan, for financially supporting this research under project number NSC 93-2213-E-212-037. The authors also wish to express appreciations to the valuable comments and suggestions from the reviewers.

References 1. Amirouche, F.M.L.: Computational Methods in Multibody Dynamics. Prentice-Hall, New Jersey (1992) 2. Ascher, U.M., Chin, H., Reich, S.: Stabilization of DAEs and invariant manifolds. Numerische Mathematik 67(2), 131–149 (1994) 3. Baruh, H.: Analytical Dynamics. McGraw-Hill, Singapore (1999) 4. Berriah, S., Lachiver, G.: Control of equilibrium and trajectory of a remotely controlled bicycle. In Proc. of the IEEE Canadian Conference on Electrical and Computer Engineering (in French) 2, 1014–1019 (1999) 5. Beznos, A.V., Formal’ sky, A.M., Gurfinkel, E.V., Jicharev, D.N., Lensky, A.V., Savitsky, K.V., Tchesalin, L.S.: Control of autonomous motion of two-wheeled bicycle with gyroscopic stabilization. In Proc. of the 1988 IEEE International Conference on Robotics & Automation, Leuven, Belgium 3, 2670–2675 (1988) 6. Chen, C.K., Dao, T.S., Yang, C.K.: Turning dynamics and equilibrium of two-wheeled vehicles. J of Mechanical Science and Technology (Special Edition) 19(1), 377–387 (2005) 7. Chen, P.H.: A scheme of fuzzy training and learning applied to elebike control system. In Proc. of the Ninth IEEE International Conference on Fuzzy Systems 2, 810–816 (2000) 8. Getz, N.H.: Internal equilibrium control of a bicycle. In Proc. of the 34th IEEE Conference on Decision and Control. New Orleans, 4285–4287 (1995) 9. Getz, N.H., Marsden, J.E.: Control for an autonomous bicycle. In Proc. of the IEEE International Conference on Robotics and Automation 2, 1397–1402 (1995) 10. Han, S., Han, J., Ham, W.: Control algorithm for stabilization of tilt angle of unmanned electric bicycle. Transaction on Control, Automation and Systems Engineering 3, 176–180 (2001) 11. Hand, R.S.: Comparisons and Stability Analysis of Linearized Equations of Motion for a Basic Bicycle Model. MSc Thesis, Cornell University (1988) 12. Jones, D.E.H.: The stability of the bicycle. Physics Today American Institute of Physics 23(4), 34–40 (1970) 13. Lee, S., Ham, W.: Self stabilizing strategy in tracking control of unmanned electric bicycle with mass balance. in Proc. of the 2002 IEEE/RSJ International Conference on Intelligent Robots and Systems 3, 2200–2205 (2002) 14. Schwab, A.L., Meijaard, J.P., Papadopoulos, J.M.: Benchmark results on the linearized equations of motion of an uncontrolled bicycle, J of Mechanical Science and Technology (Special Edition) 19(1), 292–304 (2005) 15. Sharp, R.S., Limebeer, D.J.N.: A motorcycle model for stability and control analysis. Multibody System Dynamics 6(2), 123–142 (2001) 16. Sharp, R.S., Evangelou, S., Limebeer, D.J.N.: Advances in the modelling of motorcycle dynamics. Multibody System Dynamics 12(3), 251–283 (2004) 17. Timoshenko, S., Young, D.H.: Advanced Dynamics, McGraw-Hill, New York pp. 239–240 (1948) 18. Whipple, F.J.W.: The stability of the motion of a bicycle. The Quarterly Journal of Pure and Applied Mathematics 30(120), 312–348 (1899) 19. Weir, D.H., Zellner, J.W.: Lateral-directional motorcycle dynamics and rider control. Trans. SAE 1364– 1388 (1979) 20. Yavin, Y.: Stabilization and control of the motion of an autonomous bicycle by using a rotor for tilt moment. Compute. Methods Appl. Mech. Engrg. 178, 233–243 (1999)

Springer