We present canonical linearized equations of motion for the Whipple bicycle model consisting of four rigid laterally symmetric ideally hinged parts: two wheels, a frame and a front assembly. The wheels are also axisymmetric and make ideal knife-edge rolling point contact with the ground level. The mass distribution and geometry are otherwise arbitrary. This conservative non-holonomic system has a seven-dimensional accessible configuration space and three velocity degrees of freedom parametrized by rates of frame lean, steer angle and rear wheel rotation. We construct the terms in the governing equations methodically for easy implementation. The equations are suitable for e.g. the study of bicycle self-stability. We derived these equations by hand in two ways and also checked them against two nonlinear dynamics simulations. In the century-old literature, several sets of equations fully agree with those here and several do not. Two benchmarks provide test cases for checking alternative formulations of the equations of motion or alternative numerical solutions. Further, the results here can also serve as a check for general purpose dynamic programs. For the benchmark bicycles, we accurately calculate the eigenvalues (the roots of the characteristic equation) and the speeds at which bicycle lean and steer are self-stable, confirming the century-old result that this conservative system can have asymptotic stability.
In 1818, Karl von Drais showed that a person riding forward on a contraption with two in-line wheels, a sitting scooter of sorts, could balance by steering the front wheel (Herlihy 2004). Later, the velocipede of the 1860s which had pedals directly driving the front wheel as on a child's tricycle, could also be balanced by rider-applied steering control. This ‘boneshaker’ had equal-size wooden wheels and a vertical (untilted) steering axis passing through the front wheel axle. By the 1890s, it was well known that essentially anyone could learn to balance a ‘safety bicycle’, which had pneumatic tyres and a chain drive. More subtly, but more importantly for balance and control, the safety bicycle also had a tilted steer axis and fork offset (bent front fork) like a modern bicycle. In 1897, French mathematician Carvallo (1899) and then, more generally, Cambridge undergraduate Whipple (1899) used rigid body dynamics equations to show in theory what was surely known in practice, that some safety bicycles could, if moving in the right speed range, balance themselves. Presently, these same two basic features of bicycle balance are clear.
A controlling rider can balance a forward-moving bicycle by turning the front wheel in the direction of an undesired lean. This moves the ground contact points under the rider, just like an inverted broom or stick can be balanced on an open hand by accelerating the support point in the direction of lean.
Some uncontrolled bicycles can balance themselves. If an appropriate typical bicycle is given a push to approximately 6 m s−1, it steadies itself and then progresses stably until its speed gets too low. The torques for the self-correcting steer motions can come from various geometric, inertial and gyroscopic features of the bicycle.
Beyond these two generalities, there is little that has been solidly accepted in the literature, perhaps owing to the lack of need. Through trial and error bicycles had evolved by 1890 to be stable enough to survive to the present day with essentially no modification. Since bicycle design has been based on tinkering rather than equations, there has been little scrutiny of bicycle analyses.
To better satisfy general curiosity about bicycle balance and perhaps contribute to the further evolution of bicycle design, we aim here to firmly settle some basic, and largely previously presented, bicycle stability science. The core of the paper is a set of easy-to-use and thoroughly checked linearized dynamics equations for the motion of a somewhat elaborate, yet well-defined bicycle model. These are given in equation (5.3) and appendix A. Future studies of bicycle stability aimed, for example, at clarifying especially point (ii) above, can be based on these equations.
Many methods can be used to derive the equations using various choices of coordinates, each leading to vastly different-looking governing equations. Even matching initial conditions between solution methods can be a challenge. However, the roots of the characteristic equation (eigenvalues) and the speed range of stability are independent of all of these differences. Thus, for example, a computer-based study of a bicycle based on any formulation can be checked for correctness and accuracy by comparing with the benchmark eigenvalues here.
The work here may also have more general use in that the bicycle balance problem is close to that for skating, walking and running in their coupling of lean, steer and balance. Furthermore, there is a dearth of non-trivial examples with precisely known solutions that can be used to check general purpose multibody dynamics simulators (such as used for machine, vehicle and robot design). This paper provides such a non-trivial benchmark system.
2. Brief literature review
Since their inception bicycles have attracted attention from the more or less well-known scientists of the day, including thermodynamicist William Rankine, the mathematicians Carlo Bourlet, Paul Appell and Emmanuel Carvallo, the meteorologist Francis Whipple, the mathematical physicist Joseph Boussinesq, and the physicist Arnold Sommerfeld working with mathematician Felix Klein and engineer Fritz Noether (brother of Emmy). A later peak in the ‘single-track vehicle’ dynamics literature began in about 1970, perhaps because digital computers eased integration of the governing equations, because of the increased popularity of large motorcycles (and attendant accidents), and because of an ecology-related bicycle boom. This latter literature includes work by dynamicists such as Neı̆mark, Fufaev, Breakwell and Kane. Starting in the mid-1970s, the literature increasingly deviates from the rigid body treatment that is our present focus.
Over the past 140 years, scores of other people have studied bicycle dynamics, either for a dissertation, a hobby or sometimes as part of a life's work on vehicles. This sparse and varied research on the dynamics of bicycles modelled as linked rigid bodies was reviewed in Hand (1988). A more general but less critical historical review, which also includes models with structural compliance and tyre models, is given in Sharp (1985).
Many bicycle analyses aimed at understanding rider control are based on qualitative dynamics discussions that are too reduced to capture the ability of an uncontrolled moving bicycle to balance itself. The paper by Jones (1970) is the best-known of these. The paper by Maunsell (1946) carefully considers several effects. Qualitative dynamics discussions can also be found in Lallement (1866), Rankine (1869), Appell (1896), Sharp (1896), Wallace (1929), Jones (1942), Den Hartog (1948), Higbie (1974), Kirshner (1980), Le Hénaff (1987), Olsen & Papadopoulos (1988), Patterson (1993), Cox (1998) and Wilson (2004).
A second class of papers does use analysis to study the dynamics. Some, appropriately for basic studies of rider control, use models with geometry and/or mass distribution that are too reduced to allow self-stability. Others, even if using a bicycle model that is sufficiently general, use rules for the control of the steer and thus skip the additional equation for uncontrolled steer dynamics. Such simple and/or steer-controlled approaches are found in Bourlet (1899), Boussinesq (1899a,b), Routh (1899), Bouasse (1910), Bower (1915), Pearsall (1922), Loı̆cjanskiı̆ & Lur'e (1934), Timoshenko & Young (1948), Haag (1955), Neı̆mark & Fufaev (1972), Lowell & McKell (1982), Getz & Marsden (1995), Fajans (2000), Åström et al. (2005) and Limebeer & Sharp (2006).
Finally, we have found approximately 30 rigid body dynamics models that have general enough geometry and mass distribution for self-stability to be possible, and which also allow uncontrolled steer dynamics. These governing equations are complex and different authors use slightly different modelling assumptions, different parametrizations and different choices of dynamic variables, and most authors did not know of most of their predecessors. Thus, only a small fraction of the 200 or more chronologically possible cross-checks have been performed in detail. Of these, a large fraction are by Hand and ourselves. The evaluations below are based on comparison with our own derivations (Papadopoulos 1987; Meijaard 2004; Schwab et al. 2005), and on comparisons made by the first six authors below, especially Hand.
Correct equations for the Whipple model are in Döhring (1955), who built on the less-general Carvallo model presented in Klein & Sommerfeld (1910), in Weir (1972), who checked Sharp (1971), in Eaton (1973), who checked Weir (1972) and Sharp (1971), in Hand (1988), who checked these papers and others, and in Mears (1988) who checked Weir and Hand. Singh & Goel (1971) use Döhring's correct equations, but we were unable to reproduce their eigenvalues. The paper by Dikarev et al. (1981) corrects an error in Neı̆mark & Fufaev (1972), independently also corrected later by Hand, and we have found no fault with it, but we have not confirmed the final equations. Based on graphical agreement of Psiaki's (1979) plots against solutions of the equations here, we expect that Psiaki's complex equations are correct, but we have not confirmed their formal equivalence to ours.
Equations of similar models are in Carvallo (1899), which is slightly simplified, Whipple (1899), which has some typographical errors, Klein & Sommerfeld (1910), which follows Carvallo and is slightly simplified, and Herfkens (1949), which has some typographical errors. We recently discovered a report by Manning (1951) that has no evident flaws but we have not checked it in detail. Sharp (1971) is correct before he eliminates tyre compliance and is the foundation for much subsequent tyre-based vehicle modelling. Van Zytveld (1975) is correct when his slightly incorrect and more general model is simplified to the Whipple model, and Weir & Zellner (1978) has a sign error. Neı̆mark & Fufaev (1972) has more substantial but still correctable errors (Dikarev et al. 1981; Hand 1988).
Others who studied complex rigid body bicycle models include Collins (1963), Singh (1964), Rice & Roland (1970), Roland & Massing (1971), Roland & Lynch (1972), Roland (1973), Rice (1974, 1976), Singh & Goel (1976), Lobas (1978), Koenen (1983), Franke et al. (1990), Lennartsson (1999), Åström et al. (2005) and Limebeer & Sharp (2006). We continue to discover more promising papers (Kondo et al. 1963; Ge 1966). Despite all these decades of careful good work, as of this writing there is no standard journal publication in English that we are confident presents complete and correct equations for the canonical Whipple model. Electronic supplementary material, appendix 1 expands on the historical review above.
3. The bicycle model
We use the Whipple bicycle consisting of four rigid bodies: a rear wheel R, a rear frame B with the rider body rigidly attached to it, a front frame H consisting of the handlebar and fork assembly and a front wheel F (figure 1). Within the constraint of overall lateral (left–right) and circular symmetry of the wheels, the shape and mass distributions are general with one caveat. A model that respects these symmetries allows non-planar (thick) wheels. We allow for such thickness in our inertial properties but, like Whipple, restrict attention to knife-edge rolling point contact, thus excluding, for example, contact with toroidal wheels. We neglect the motion of the rider relative to the frame, structural compliances and dampers, joint friction and tyre models with compliance and slip.
The model delineation is not by selecting the most important aspects for describing real bicycle behaviour. For understanding basic features of active rider control, the model here is undoubtedly unnecessarily and inappropriately complex. For example, some aspects included here have very small effects, like the non-planarity of the inertia of the real wheel, and other neglected aspects may be paramount, e.g. the rider's flexibility and control reflexes. Even for the study of uncontrolled stability, tyre deformation and frame compliance seem necessary for understanding wobble (a rapid steering oscillation). In summary, the model here includes all the sharply defined rigid body effects, while leaving out a plethora of terms that would require more subtle and less well-defined modelling.
Our bicycle design is fully characterized by 25 parameters described below. Table 1 lists the numerical values used for the numerical benchmark. Most numerical values are representative of real bicycles, but some values (e.g. wheel inertial thickness as represented by IRxx>IRyy/2) are exaggerated to guarantee a detectable role in the benchmark numerical studies. The bicycle design parameters are defined in an upright reference configuration with both wheels on the level flat ground and with zero steer angle. The reference coordinate origin is at the rear wheel contact point P. We use the conventions of vehicle dynamics (J670e, SAE 2001) with positive x pointing generally towards the front contact point, positive z pointing down and the y-axis pointing to the rider's right.
The radii of the circular wheels are rR and rF. The wheel masses are mR and mF with their centres of mass at the wheel centres. The moments of inertia of the rear and front wheels about their axles are IRyy and IFyy, respectively. The moments of inertia of the wheels about any diameter in the xz-plane are IRxx and IFxx. The wheel mass distribution need not be planar, so any positive inertia is allowed with IRyy≤2IRxx and IFyy≤2IFxx. All front wheel parameters can be different from those of the rear so, for example, it is possible to investigate separately the importance of angular momentum of the front and rear wheels.
Narrow high-pressure, high-friction tyre contact is modelled as non-slipping rolling point-contact between the ground and the knife-edge wheel perimeters. The frictionless wheel axles are orthogonal to the wheel symmetry planes and are located at the wheel centres. In the reference configuration, the front wheel ground contact Q is located at a distance w (the wheel base) in front of the rear wheel contact P. The front wheel ground contact point trails a distance c behind the point where the steer axis intersects the ground. Although c>0 for most bicycles, the equations allow a negative trail (c<0) with the wheel contact point in front of the steer axis.
The rear wheel R is connected to the rear frame assembly B (which includes the rider body) at the rear axle. The centre of mass of B, with mass mB, is located at (xB, yB=0, zB<0). The moment of inertia of the rear frame about its centre of mass is represented by a 3×3 moment of inertia matrix where all mass is symmetrically distributed relative to the xz-plane, but not necessarily on the plane. The centre of mass of the front frame assembly (fork and handlebar) H is at (xH, yH=0, zH<0) relative to the rear contact P. H has mass mH. As for the B frame, IHyy can be less than IHxx+IHzz. The moment of inertia matrices of the rear and front assemblies are(3.1)To be physical (i.e. no negative mass), the moment-of-inertia matrix must have all principal values positive and also satisfy the triangle inequalities that no one principal value is bigger than the sum of the other two. The steer axis tilt angle λ is measured from the upwards vertical, positive when tipped back as on a conventional bicycle with −π/2<λ<π/2 (all angles in radians). The steer tilt is π/2 minus the conventional ‘head angle’; a bicycle with head angle of 72° has λ=18°=π/10. The steer axis location is implicitly defined by the wheel base w, trail c and steer axis tilt angle λ.
Two non-design parameters are the downward gravitational acceleration g and the nominal forward speed v. This model, or slight simplifications of it, is a common idealization of a bicycle (electronic supplementary material, appendix 1). Motorcycle modelling is often based on an extension of this model using toroidal wheels, tyre compliance, tyre slip and frame compliance. Theories of bicycle and motorcycle control are often based on simplifications of this model or, alternatively, on simple analogous systems that do not come from reductions of this model.
(a) How many parameters describe the bicycle model?
The bicycle model here is defined completely by the 25 design parameters above (table 1). This is not a minimal description for dynamic analysis, however. For example, the inertial properties of the rear wheel R, except for the polar moment of inertia (i.e. mR and IRxx but not IRyy), can be combined with the inertial properties of the rear frame B, reducing the number of parameters by 2. The same combination can be used for the front frame, reducing the number of parameters to 25−2−2=21. The polar inertia of each wheel can be replaced with a gyrostat constant that gives its spin angular momentum in terms of forward velocity. This does not reduce the number of parameters in nonlinear modelling. But in linear modelling, the radius of the wheels is irrelevant for lean and steer geometry and their effect on angular momentum is embodied in the gyrostat constants. Eliminating wheel radii reduces the number of parameters by 2 to 21−2=19. Finally, in the linearized equations of motion, the polar (yy components) of the moments of inertia of the two frames are irrelevant, reducing the necessary number of design parameters to 19−2=17.
In their most reduced form, the linearized equations of motion (5.3) have 11 arbitrary independent matrix entries. Each entry is a complex combination of the 17 parameters just described. Still, further reduction can be obtained by inspection of the fourth-order characteristic equation (6.5). After scaling by the leading coefficient det(M), there remain four coefficients, each a polynomial in the forward speed. There are seven independent coefficients of these velocity polynomials. By reduction using suitable length and time-scales, two of these coefficients can be eliminated. Thus, the space of scaled root loci plots is only five-dimensional. For simpler comparisons, we use all 25 design parameters.
(b) How many degrees of freedom does the bicycle model have?
Since this system has non-holonomic kinematic constraints, the concept of ‘degree of freedom’ needs clarification. The holonomic (hinges and ground contacts) and non-holonomic (non-slip rolling) constraints restrict these four linked three-dimensional objects in space as follows. Start with the 24 degrees of freedom of the four rigid bodies, each with three translational and three rotational degrees of freedom in physical space (4×(3+3)=24). Then, subtract five degrees of freedom for each of the three hinges and one more for each wheel touching the ground plane: 24−3×5−2=7. Thus, before we consider the non-slipping wheel contact constraints, the accessible configuration space is seven-dimensional. The four non-holonomic rolling constraints (two for each wheel-to-ground contact) do not further restrict this accessible configuration space: kinematically allowable rolling motions can translate and steer the bicycle on the plane in arbitrary ways and also can rotate the wheels relative to the frame with no net change of overall bicycle position or orientation. For example, even though side-slip is not allowed, a bicycle can move sideways by the same motions used to parallel-park a car. Thus, the accessible configuration space for this model is seven-dimensional.
(i) Description of the seven-dimensional configuration space
This seven-dimensional configuration space can be parametrized as follows (figure 2). The location of the rear wheel contact with the ground is (xP, yP) relative to a global fixed coordinate system with origin O. The orientation of the rear frame with respect to the global reference frame O-xyz is given by a sequence of angular rotations (312 Euler angles). These rotations are shown in figure 2 with fictitious hinges, each represented as a can in the drawing, in series, mounted at the rear hub: a yaw rotation, ψ, about the z-axis, a lean rotation, ϕ, about the rotated x-axis and a pitch rotation, θB, about the rotated y-axis. Note that the pitch θB is not one of the seven configuration variables because it is determined by a three-dimensional trigonometric relation that keeps the steered front wheel on the ground. The steering angle δ is the rotation of the front handlebar frame with respect to the rear frame about the steering axis. A right turn of a forward-moving bicycle has δ>0. Finally, the rotation of the rear R and front F wheels with respect to their respective frames B and H are θR and θF. In summary, the configuration space is parametrized here with (xP, yP, ψ, ϕ, δ, θR, θF). Quantities such as wheel centre coordinates and rear frame pitch are all determined by these.
(ii) Velocity degrees of freedom
As explained previously, the accessible configuration space is seven-dimensional. However, the four non-holonomic rolling constraints (one longitudinal and one lateral for each wheel-to-ground contact) reduce the seven-dimensional accessible configuration space to 7−4=3 velocity degrees of freedom. This three-dimensional, kinematically accessible, velocity space can conveniently be parametrized by the lean rate of the rear frame, the steering rate and the rotation rate of the rear wheel R relative to the rear frame B.
4. Basic features of the model, equations and solutions
(a) The system behaviour is unambiguous
The dynamics equations for this model follow from linear and angular momentum balance applied to each part, along with the assumption that the kinematic constraint forces follow the rules of action and reaction and do no net work. These equations may be assembled into a set of ordinary differential equations, or differential-algebraic equations by various methods. One can assemble governing differential equations using the Newton–Euler rigid body equations, Lagrange equations with Lagrange multipliers for the in-ground-plane rolling-contact forces, or methods based on the principle of virtual velocities (e.g. Kane's method), etc. But the subject of mechanics is sufficiently well defined that we know that all standard methods will yield equivalent sets of governing differential equations. Therefore, a given consistent-with-the-constraints initial state (positions and velocities of all points on the frames and wheels) will always yield the same subsequent motions of the bicycle parts. Thus, while the choice of variables and the recombination of governing equations may lead to quite different-looking governing equations, any difference between dynamics predictions can only be due to errors.
(b) The system is conservative but not Hamiltonian
The only friction forces in this model are the lateral and longitudinal forces at the ground contact points. But, owing to the no-slip conditions (constraints) these friction forces are modelled by non-dissipative constraint forces. The hinges and ground contacts are all workless kinematic constraints. In uncontrolled bicycle motion, the only external applied forces are the conservative gravity forces on each part. That is, there are no dissipative forces and the system is energetically conservative; the sum of the gravitational and kinetic energies is a constant for any free motion. But the non-holonomic kinematic constraints preclude writing the governing equations in standard Hamiltonian form, so theorems of Hamiltonian mechanics do not apply. One result, surprising to some cultured in Hamiltonian systems, is that the non-dissipative bicycle equations can have asymptotic (exponential) stability (figure 4). This apparent contradiction of the stability theorems for Hamiltonian systems is because the bicycle, while conservative, is, by virtue of the non-holonomic wheel contacts, not Hamiltonian. A similar system that is conservative but has asymptotic stability is the uncontrolled skateboard (Hubbard 1979) and simpler still is the classical Chaplygin sleigh described in e.g. Ruina (1998).
(c) Symmetries in the solutions
Without explicit use of the governing equations, some features of their solutions may be inferred by symmetry.
(i) Ignorable coordinates
Some of the configuration variables do not appear in any expression for the forces, moments, potential energies or kinetic energies of any of the parts. These so-called cyclic or ignorable coordinates are: the location of the bicycle on the plane (xP, yP), the heading of the bicycle ψ and the rotations (θR, θF) of the two wheels relative to their respective frames. Thus, one can write a reduced set of dynamics equations that do not include these ignorable coordinates. The full configuration as a function of time can be found afterwards by integration of the kinematic constraint equations, as discussed at the end of appendix B. These ignorable coordinates cannot have asymptotic stability; a small perturbation of, say, the heading ψ will lead to a different ultimate heading.
(ii) Decoupling of lateral dynamics from speed dynamics
The lateral (left–right) symmetry of the bicycle design along with the lateral symmetry of the equations implies that the straight-ahead unsteered and untipped (δ=0, ϕ=0) rolling motions are necessarily solutions for any forward or backward speed v. Moreover, relative to these symmetric solutions, the longitudinal and the lateral motion must be decoupled from each other to first order (linearly decoupled) by the following argument. Owing to lateral symmetry, a perturbation to the right must cause the same change in speed as a perturbation to the left. But by linearity, the effects must be the negative of each other. Therefore, there can be no first-order change in speed due to lean. Similarly, speed change alone cannot cause lean. Thus, the linearized fore–aft equations of motion are entirely decoupled from the lateral equations of motion and a constant speed bicycle has the same linearized equations of lateral motion as a constant energy bicycle. This argument is given in more detail in the electronic supplementary material, appendix 4.
(iii) A fore–aft symmetric bicycle cannot be asymptotically self-stable
For any rigid body system with workless kinematic constraints and state-dependent forces, any solution q(t) implies a solution q(−t). Thus, any bicycle motion is also a solution of the equations when moving backwards, with all particle trajectories being traced at identical speeds in the reverse direction. Thus, a bicycle that is exponentially stable in balance when moving forwards at speed v>0 must be exponentially unstable when moving at −v (backwards at the same speed). Consider a fore–aft symmetric bicycle. Such a bicycle has a vertical central steering axis (or a horizontal steering axis) and has a handlebar assembly, front mass distribution and front wheel that mirrors that of the rear assembly. If such a bicycle has exponentially decaying solutions in one direction, it must have exponentially growing solutions in the opposite direction owing to time reversal. By symmetry, it must therefore also have exponentially growing solutions in the (supposedly stable) original direction. Thus, such a bicycle cannot have exponentially decaying solutions in one direction without also having exponentially growing solutions in the same direction, and thus cannot be asymptotically self-stable. Such a symmetric bicycle might, however, have the oscillatory (neutral) stability of the type Hamiltonian systems can have.
(d) The nonlinear equations have no simple expression
In contrast with the linear equations we present below, there seems to be no reasonably compact expression of the full nonlinear equations of motion for this model. The kinematic loop, from rear wheel contact to front wheel contact, determines the rear frame pitch through a complicated equation (Psiaki 1979) which can be expressed as a quartic polynomial equation. Therefore, there is no simple expression for rear frame pitch for large lean and steer angles. Thus, the writing of nonlinear governing differential equations in a standard form that various researchers can check against alternative derivations is a challenge that is not addressed here, and might never be addressed. However, when viewed as a collection of equations, one for each part, and a collection of constraint equations, a large set of separately comprehensible equations may be assembled. An algorithmic derivation of nonlinear equations using such an assembly, suitable for numerical calculation and benchmark comparison, is presented in Basu-Mandal et al. (2007) where a complete set of no-hands circular motions are also presented.
5. Linearized equations of motion
Here, we present a set of linearized differential equations for the bicycle model, slightly perturbed from upright straight-ahead motion, in a canonical form. To aid in organizing the equations, we include applied lean and steer torques which are later set to zero for study of uncontrolled motion.
(a) Derivation of governing equations
Mostly correct derivations and presentations of the equations of motion for a relatively general bicycle model, although not necessarily expressed in the canonical form of equation (5.3), are found in Carvallo (1899), Whipple (1899), Klein & Sommerfeld (1910), Herfkens (1949), Döhring (1953, 1955), Sharp (1971), Weir (1972), Eaton (1973) and Van Zytveld (1975). Dikarev et al. (1981) have a derivation of equations equivalent to equation (5.3), based on correcting the errors in Neı̆mark & Fufaev (1972) as does Hand (1988), which just predates Mears (1988). Papadopoulos (1987) and Meijaard (2004) also have derivations that were generated in preparation for this paper.
The derivations above are generally long, leading to equations with layers of nested definitions. This is at least part of the reason for the lack of cross-checking in the literature. A minimal derivation of the equations using angular momentum balance about various axes, based on Papadopoulos (1987), is given in appendix B. Note that this derivation, as well as all of the linearized equations from the literature, is not based on a systematic linearization of full nonlinear differential equations. Thus far, systematic linearizations have not achieved analytical expressions for the linearized equation coefficients in terms of the 25 bicycle parameters. However, part of the validation process described later includes numerical comparison with full nonlinear simulations, and also comparison with numerical values of the linearized equation coefficients as determined by these same nonlinear programs.
(b) Forcing terms
For numerical benchmark purposes, where eigenvalues are paramount, we neglect control forces or other forcing (except gravity, which is always included). However, the forcing terms help to organize the equations. Moreover, forcing terms are needed for study of disturbances and rider control.
In addition to the gravity forces, consider an arbitrary distribution of forces Fi acting at various points on the bicycle. Their net effect is to contribute to the forces of constraint (the ground reaction forces, and the action–reaction pairs between the parts at the hinges) and to contribute to the accelerations . Three generalized forces can be defined by writing the power of the applied forces, kept at their current values, associated with arbitrary perturbations of the velocities that are consistent with the hinge assembly and ground wheel contact constraints. This ‘virtual’ power necessarily factors into a sum of three terms(5.1)because the perturbations of the velocities Δvi of all material points are necessarily linear combinations of the perturbations of the generalized velocities . The generalized forces are thus each linear combinations of the components of the various applied forces Fi.
The generalized forces are energetically conjugate to the generalized velocities. The generalized forces can be visualized by considering special loadings each of which contributes to only one generalized force when the bicycle is in the reference configuration. In this way of thinking
is the propulsive ‘force’, expressed as an equivalent moment on the rear wheel. In practice, pedal torques or a forward push on the bicycle contribute to and not to Tϕ and Tδ.
Tϕ is the right lean torque, summed over all the forces on the bicycle, about the line between the wheel ground contacts. A force, perpendicular on the rear frame located directly above the rear contact point contributes only to Tϕ. A sideways wind gust, or a parent holding a beginning rider upright, contributes mainly to Tϕ.
Tδ is an action–reaction steering torque. A torque causing a clockwise (looking down) action to the handlebar assembly H along the steer axis and an equal and opposite reaction torque on the rear frame contributes only to Tδ. In simple modelling, Tδ would be the torque that a rider applies to the handlebars. Precise description of how general lateral forces contribute to Tδ depends on the projection implicit in equation (5.1). Some lateral forces make no contribution to Tδ, namely those acting at points on either frame which do not move when an at-rest bicycle is steered but not leaned. Lateral forces applied to the rear frame directly above the rear contact point make no contribution to Tδ. Nor do forces applied to the front frame if applied on the line connecting the front contact point with the point where the steer axis intersects the vertical line through the rear contact point. Lateral forces at ground level, but off the two lines just described, contribute only to Tδ. Lateral forces acting at the wheel contact points make no contribution to any of the generalized forces.
Just as for a pendulum, finite vertical forces (additional to gravity) change the coefficients in the linearized equations of motion but do not contribute to the forcing terms (e.g. a magnet under a pendulum changes the effective g in but not the 0). Similarly, propulsive forces also change the coefficients but have no first-order effect on the lateral forcing. Thus, the equations presented here only apply for small propulsive and small additional vertical forces.
(c) The first linear equation: with no forcing, forward speed is constant
A solution of both the linearized and the full nonlinear equations is straight-ahead δ=0 upright ϕ=0 motion at any constant forward speed . The governing equations here describe the evolution of small perturbations from this reference solution. As explained above and in more detail in the electronic supplementary material, appendix 4, lateral symmetry of the system, combined with the linearity in the equations, precludes any first-order coupling between the forward motion and the lean and steer. Therefore, the first linearized equation of motion is simply obtained from two-dimensional (xz-plane) mechanics as(5.2)where mT is total bicycle mass (appendix A). In cases with no propulsive force, the nominal forward speed is therefore constant (to first order).
(d) Lean and steer equations
The linearized equations of motion for the two remaining degrees of freedom, the lean angle ϕ and the steer angle δ, are two coupled second-order, constant coefficient, ordinary differential equations. Any such set of equations can be linearly combined to get an equivalent set. We define the canonical form below by insisting that the right-hand sides of the two equations consist only of Tϕ and Tδ, respectively. The first of these two equations is the lean equation and the second is the steer equation. Mechanical systems have linearized equations of the form . For the bicycle model, these equations can be written, with velocity and gravity explicit, as (Papadopoulos 1987)(5.3)where the time-varying quantities are q=[ϕ, δ]T and f=[Tϕ, Tδ]T. The constant entries in matrices M, C1, K0 and K2 are defined in terms of the 25 design parameters in appendix A. Briefly, M is a symmetric mass matrix which gives the kinetic energy of the bicycle system at zero forward speed as . The damping-like (there is no real damping) matrix C=vC1 is linear in the forward speed v and captures skew-symmetric gyroscopic torques due to steer and lean rates. C1 also captures inertial reactions due to steer rate: part from the lateral acceleration due to rear wheel yaw rate (path curvature) and part from system yaw acceleration. The stiffness matrix K is the sum of two parts: a velocity-independent symmetric part gK0 proportional to the gravitational acceleration, which can be used to calculate changes in potential energy with qT[gK0]q/2, and a part v2K2, which is quadratic in the forward speed and is due to gyroscopic and centrifugal effects. The matrix subscripts match the exponents of the v multipliers.
6. Benchmark model and solutions
To facilitate comparisons with other formulations we have defined a benchmark bicycle with all parameter values given in table 1. The parameter values were chosen to minimize the possibility of fortuitous cancellation that could occur if used in an incorrect model. We also wanted numbers that could easily be described precisely. In the benchmark bicycle, the two wheels are different in all properties and no two angles, masses or distances match. A second simpler benchmark is in the electronic supplementary material, appendix 5.
(a) Coefficients of the linearized equations of motion
Substitution of the values of the design parameters for the benchmark bicycle from table 1 in the expressions from appendix A results in the following values for the entries in the matrices in the equations of motion (5.3):(6.1)(6.2)(6.3)(6.4)To serve as a precise benchmark, the coefficients are given with 14 decimal places (trailing zeros suppressed) above and elsewhere. Many-digit agreement with results obtained by other means provides near certainty that there is also an underlying mathematical agreement, even if that agreement is not apparent analytically.
(b) Linearized stability, eigenvalues for comparison
The characteristic equation eigenvalues are independent of coordinate choice and equation form; any non-singular change of variables yields equations with the same eigenvalues. Thus, eigenvalues serve as convenient benchmark results, permitting comparison between different analyses. The eigenvalues are calculated by assuming an exponential solution of the form q=q0 exp(λt) for the homogeneous equations (f=0 in equation (5.3)). This leads to the characteristic polynomial equation,(6.5)which is quartic in λ. After substitution of the expressions from appendix A, the coefficients in this quartic polynomial become complicated expressions of the 25 design parameters, gravity and speed v. The zeros λ of the characteristic polynomial for a range of forward speeds are shown in figure 3. Eigenvalues with a positive real part correspond to unstable motions, whereas eigenvalues with a negative real part correspond to asymptotically stable motions for the corresponding mode. Imaginary eigenvalues correspond to oscillatory motions. As mentioned in §4c, the time reversability of this system leads to a symmetry evident in the characteristic equation (6.5): if (v, λ) is a solution then (−v, −λ) is also a solution. This means that figure 3 is point symmetric about the origin as revealed in fig. 9 of Åström et al. (2005).
This fourth-order system has four distinct eigenmodes except at special parameter values associated with multiple roots. A complex (oscillatory) eigenvalue pair is associated with a pair of complex eigenmodes. At high enough speeds, the two modes most significant for stability are traditionally called the capsize mode and weave mode. The capsize mode corresponds to a real eigenvalue with eigenvector dominated by lean: when unstable, a capsizing bicycle leans progressively into a tightening spiral with steer and lean both increasing proportionally as it falls over. The weave mode is an oscillatory motion in which the bicycle steers sinuously about the headed direction with a slight phase lag relative to leaning. The third eigenvalue is large, real and negative. It corresponds to the castering mode, which is dominated by steer in which the front ground contact follows a tractrix-like pursuit trajectory, like the straightening of a swivel wheel under the front of a grocery cart.
At near-zero speeds, typically 0<v<0.5 m s−1, there are two pairs of real eigenvalues. Each pair consists of a positive and a negative eigenvalue and corresponds to an inverted pendulum-like falling of the bicycle. The positive root in each pair corresponds to falling, whereas the negative root corresponds to the time reversal of this falling (i.e. rising). For the benchmark bicycle here, when speed is increased to vd≈0.7 m s−1, two real eigenvalues coalesce and then split to form a complex conjugate pair; this is where the oscillatory weave motion emerges. At first, this motion is unstable but at vw≈4.3 m s−1, the weave speed, these eigenvalues cross the imaginary axis in a Hopf bifurcation (Strogatz 1994) and this mode becomes stable. At a higher speed, the capsize eigenvalue crosses the origin in a pitchfork bifurcation at vc≈6.0 m s−1, the capsize speed, and the bicycle becomes mildly unstable (Basu-Mandal et al. 2007). The speed range for which the uncontrolled bicycle shows asymptotically stable behaviour, with all eigenvalues having negative real parts, is vw<v<vc. For comparison by future researchers, benchmark eigenvalues are presented at various forward speeds in table 2.
7. Validation of the linearized equations of motion
The linearized equations of motion here, equation (5.3) with the coefficients as presented in appendix A, have been derived by pencil and paper in two ways (Papadopoulos 1987; Meijaard 2004), and agree exactly with some of the past literature (§2). We have also checked equation coefficients via the linearization capability of two general nonlinear dynamics simulation programs described below. Comparisons with the work here using nonlinear simulations have also been performed by A. Lennartsson (2006, personal communication, based on Lennartsson 1999) and Basu-Mandal et al. (2007). Finally, in the self-stable speed range, steering and lean transients can be measured on a physical bicycle with narrow high-pressure tyres. Kooijman et al. (2007) measured the mass and geometry parameters on a riderless bicycle and found good comparison between the experimentally measured eigenvalues and the eigenvalues predicted by the formulae here.
(a) Equations of motion derived with the numeric program Spacar
Spacar is a program system for dynamic simulation of multibody systems developed by Van der Werff (1977), Jonker (1988), Jonker & Meijaard (1990), Meijaard (1991), Schwab (2002) and Schwab & Meijaard (2003). Spacar is based on finite element principles laid out by Besseling (1964). Spacar handles systems of rigid and flexible bodies connected by various joints in both open and closed kinematic loops, and where parts may have rolling contact. Spacar generates numerically, and solves, full nonlinear dynamics equations using minimal coordinates (constraints are eliminated). The Spacar model used in this paper uses the rigid body, point mass, hinge and rolling wheel contact features of the program (Schwab & Meijaard 1999, 2003). Spacar can also find the numeric coefficients for the linearized equations of motion based on a semi-analytic linearization of the nonlinear equations. As determined by Spacar, the entries in the matrices of the linearized equations of motion (5.3) agree to 14 digits with the values presented in §6a (see electronic supplementary material, appendix 2 for more about the Spacar model).
(b) Equations of motion derived with the symbolic program AutoSim
We also derived the nonlinear governing equations using the multibody dynamics program AutoSim, v. 2.80 (Sayers 1991a,b). AutoSim is a Lisp (Steele 1990) program mostly based on Kane's (1968) approach. It consists of function definitions and data structures allowing the generation of symbolic equations of motion of rigid body systems. AutoSim works best for systems of objects connected with prismatic and revolute joints arranged with the topology of a tree (no loops). AutoSim generates equations in the form(7.1)Here, q are the generalized coordinates; u are the generalized velocities; S is the kinematic matrix that relates the rates of the generalized coordinates to the generalized speeds; M is the system mass matrix; and Q contains all force terms and velocity-dependent inertia terms. Additional constraints are added for closed kinematic loops, special joints and non-holonomic constraints. For example, the closed loop holonomic constraint for both bicycle wheels touching the ground cannot be solved simply in symbolic form for the dependent coordinates, requiring the solution of a quartic polynomial equation. An iterative numerical solution for this constraint was used, destroying the purely symbolic nature of the equations.
In general, the standard linearization procedure of AutoSim is not applicable for systems with closed loops. However, in the present case, the dependent pitch angle is zero to first order and no difficulties arise. The final AutoSim-based linearization output consists of a MATLAB script file that numerically calculates the matrices of the linearized equations. The entries in the matrices of the linearized equations of motion (5.3) as determined by the program AutoSim agree to 14 digits with the values presented in §6a. More details about the AutoSim verification are given in the electronic supplementary material, appendix 3.
8. Energy conservation and asymptotic stability
When an uncontrolled bicycle is within its stable speed range, lean and steer perturbations die away in a seemingly damped fashion. However, the system has no true damping and conserves energy. The energy in the lean and steer oscillations is transferred to the forward speed rather than being dissipated. As the forward speed is affected only to second order, linearized equations do not capture this energy conservation. Therefore, a nonlinear dynamic analysis with Spacar was performed on the benchmark bicycle model to demonstrate the flow of energy from lateral perturbations into forward speed. The initial conditions at t=0 are the upright reference position (ϕ, δ, θR)=(0, 0, 0) at a forward speed of v=4.6 m s−1, which is within the stable speed range of the linearized analysis, and a lean rate of . In the full nonlinear equations, the final upright forward speed is augmented from the initial speed by an amount determined by the energy in the lateral perturbation. In this case, the speed-up was approximately 0.022 m s−1. Figure 4 shows this small increase in the forward speed v, while the lateral motions die out, as expected. Figure 4 also shows that the period for the lean and steer oscillations is approximately T0=1.60 s, which compares well with the 1.622 s from the linearized stability analysis. The lack of agreement in the second decimal place is from finite amplitude effects, not numerical accuracy issues. When the initial lateral velocity is decreased by a factor of 10, the period of motion matches the linear prediction to four digits. The steering motion has a small phase lag relative to the lean motion , visible in the solution in figure 4.
9. Conclusions, discussion and future work
We have presented reliable equations for a well-delineated model for more deeply studying controlled and uncontrolled stability of a bicycle. The equations of motion, equation (5.3) with appendix A, are buttressed by a variety of historical and modern simulation comparisons and, we feel, can be used with confidence. They can also be used as a check for those who derive their own equations by comparison to
the tabulated eigenvalues, or
the speed range of self-stability for the benchmark parameters.
This paper firms up Carvallo's discovery in 1897 that asymptotic self-stability of an uncontrolled bicycle is explicable with a sufficiently complex conservative rigid body dynamics model. It only narrowly answers the question ‘how does an uncontrolled bicycle stay up?’ by showing that it follows from the equations. A simple explanation does not seem possible because the lean and steer are coupled by a combination of several effects including gyroscopic precession, lateral ground reaction forces at the front wheel ground contact point trailing behind the steering axis, gravity and inertial reactions from the front assembly having centre-of-mass offset from the steer axis, and from effects associated with the moment of inertia matrix of the front assembly.
The equations here can be the basis for future work addressing how bicycle self-stability does and does not depend on the bicycle design parameters. For example, we hope to dispel some bicycle mythology about the need for mechanical trail or gyroscopic wheels for bicycle self-stability.
Technical and editorial comments from Karl Åström, Anindya Chatterjee, Andrew Dressel, Neil Getz, Richard Klein, Anders Lennartsson, David Limebeer, Mark Psiaki, Keith Seffen, Alessandro Saccon, Manoj Srinivasan and five anonymous reviewers have improved the paper. Claudine Pouret at the French Academy of Sciences provided information about the Prix Fourneyron. J.P.M. was supported by the Engineering and Physical Sciences research Council (EPSRC) of UK A.R. and J.M.P. were supported by an NSF presidential Young Investigator Award and A.R. further partially supported by NSF biomechanics and robotics grants.