It has been shown recently that the shimmy motion of towed wheels can be predicted in a wide range of parameters by means of the so-called memory effect of tyres. This delay effect is related to the existence of a travelling-wave-like motion of the tyre points in contact with the ground relative to the wheel. This study shows that the dynamics within the small-scale contact patch can have an essential effect on the global dynamics of a four-wheeled automobile on a large scale. The stability charts identify narrow parameter regions of increased fuel consumption and tyre noise with the help of the delay models that are effective tools in dynamical problems through multiple scales.
In engineering, the prediction of the dynamical behaviour of designed systems often requires the use of many degrees of freedom (d.f.) models leading to large systems of ordinary differential equations (ODEs) or partial differential equations (PDEs). The design of these systems needs approximate analytical models and solutions, too, in order to study the effects of certain parameters, to explore new physical phenomena of the systems, and also to check the correctness and precision of the corresponding numerical solutions. When analytical solutions to multiple scale problems are to be constructed, delay differential equations (DDEs) can be used efficiently. These DDEs can still capture the infinite dimensional structure of the corresponding phase spaces while using fewer parameters and state variables than do large ODE or PDE systems. Another advantage of DDEs is that they use time only as an independent variable, so their analytical study can use tools and methods originated directly in the theory of ODEs.
There are several examples representing the efficiency of delay models. As one of the most characteristic ones, consider the recently introduced delayed model of ghost (or phantom) traffic jams [1,2]. It has been shown by Orosz & Stépán  that the large-scale dynamics of cars following each other on a road is essentially affected by the delay in the drivers' reflexes; actually, this is the reason why the possible Hopf bifurcations are all subcritical in these cases. If we wanted to describe this system with the correct modelling of the small-scale nerve systems of the drivers, then the resulting system of coupled PDEs and ODEs has a ‘large’ state space with many (infinite dimensional) state variables and a large number of system parameters. If, however, the whole small-scale dynamic behaviour is substituted by a single delay parameter of a driver's reflex, and a certain gain with which the travelling-wave-like signal arrives because of the lag between the driver seeing something and acting accordingly, then the large-scale dynamics of the traffic can be captured, and its intricate dynamics can still be explored with a few relevant parameters in an infinite dimensional phase space. Another similar recent example is the gene regulatory networks  where, again, the time delay is the most important parameter that characterizes small-scale dynamics and also retains the infinite dimensional character of the phase space for the large-scale dynamics.
The same philosophy is in use in the case of the so-called delayed tyre model when the intricate dynamics of the rolling tyre is to be explored, and new generation tyres are to be designed. Clearly, mechanical modelling of tyres requires different levels of simplifications and approximations of the real complex structure in order to have tyre models  that can be handled by the available mathematical and/or computational tools. The dynamics of a vehicle strongly depends on the tyre behaviour, which involves many key nonlinearities from essential geometrical ones to physical ones such as the non-smooth friction force characteristics. All these make tyre dynamics one of the most challenging research tasks in vehicle dynamics.
Accurate and detailed mechanical modelling of the tyres are required for numerical experiments to study specific motions and vibrations of vehicles. Such tyre models are presented elsewhere [6,7,8,9]. One type of ‘introductory problem’ to the field is the fascinating phenomenon of the shimmy motion of a single towed wheel. The results gained first for simple wheel shimmy showed that many vibration issues of vehicles originate in the dynamics of the tyres. These vibration issues include the reduction in rolling resistance and the noise reduction in the rolling tyre , which are currently important aspects of research and development that require sophisticated tyre modelling .
Although some studies on shimmy were carried out during the first half of the twentieth century , the elimination of shimmy is an ongoing problem that requires the special attention of engineers even today. One of the first detailed scientific reports on this phenomenon was written by von Schlippe & Dietrich , who introduced the basic ideas of the so-called stretched string tyre model and the delay effect. Nevertheless, the mathematical model and its analysis suffered from unjustified approximations because of the lack of the rigorous mathematical theory of DDEs that began to be developed later.
Pacejka  developed a quasi-static approximation for lateral tyre forces in order to construct a manageable and still realistic model for shimmying wheels. Moreover, Pacejka made his creep-force idea ready for use in other industrial problems by means of his ‘magic formula’, which was also implemented successfully in anti-lock brake systems of vehicles. Using Pacejka's tyre model, the linear stability analyses of motorcycles and aircraft main landing gear were also studied accurately by Sharp et al.  and Besselink , respectively.
In spite of these early successes, some unexpected shimmy phenomena have been experienced with low probability but with large damage costs in vehicles, such as self-excited, quasi-periodic and chaotic oscillations induced by large perturbations. The study of these ‘mystical’ vibrations required the nonlinear analysis of shimmy. It was shown, for example, by Takács et al.  that geometrical nonlinearities even with rigid wheels can lead to such separated bistable parameter domains, in which the stable stationary rolling of a towed wheel coexists with large-amplitude shimmy motion, a phenomenon that originates in a subcritical Hopf bifurcation. A practically more relevant model of an aircraft nose landing gear with reduced geometrical nonlinearities but with an elastic tyre model was investigated by Thota et al. , who explained that quasi-periodic vibration is the interaction of the torsional and bending vibration modes of the suspension system.
Stépán  demonstrated quasi-periodic vibration in a single rigid-body d.f. mechanical model, in which the lateral tyre deformation in the contact patch was described without any stationary approximation for the deformation shape. Namely, the memory effect of the tyre–ground contact was taken into account, which was neglected in the simplified tyre models. The existence of these quasi-periodic vibrations was validated later in laboratory experiments by Takács . The exact solution of the tyre–ground contact problem also predicted some new linearly unstable parameter domains for wheels towed with long casters. Takács & Stépán  showed that the non-smooth friction force characteristics in the contact patch limit the amplitude of the emerging self-excited vibrations in the corresponding parameter domains. The resulting so-called ‘micro-shimmy’ of the towed wheel was detected by means of measuring the increased tyre temperature and the increased rolling resistance.
All these preliminary results for towed wheels motivated us to implement a delayed tyre description for vehicle models. In this paper, the instantaneous lateral tyre deformations in thin contact patches are calculated with mathematical rigour, whereas a simple so-called bicycle model of a passenger car is used with the simplest so-called brush-type tyre. The linear stability of the straight stationary motion of a vehicle is analysed and new unstable parameter domains are identified, which are candidates for micro-shimmy leading to increased vehicle noise and fuel consumption. The results demonstrate the use of delay models in multi-scale problems, because the small-scale dynamics in the contact patch has an essential effect on the intricate large-scale dynamics of the vehicle.
2. Mechanical model
The mechanical model in question is shown in figure 1. This classical model of a four-wheeled vehicle is constructed by means of the following assumptions. The vertical dynamics of the vehicle is neglected, and the vehicle moves in the ground-fixed (X,Y) coordinate system. Because the linear stability of the straight stationary rolling of the vehicle will be analysed, the so-called bicycle model can be used. This simplified mechanical model eliminates the lateral dimension of the vehicle that has no effect on the lateral stability in the case of a symmetric car body and wheel set-up. Thus, the car body is substituted by a perfectly rigid beam between the front and the rear wheels. The distance between the wheels, namely the axle base, is 2l, and the mass of the vehicle and mass moment of inertia of the vehicle with respect to the z-axis at the centre of gravity C are denoted by m and JC, respectively. The location of the centre of gravity C along the beam is given by the eccentricity e. In order to study the straight stationary motion of the vehicle, the steering angle is set to be zero (δ=0). The (x,y) coordinate system is fixed to the vehicle, and the coordinates x and y will be called the longitudinal and lateral directions, respectively.
In order to describe the planar motion of the vehicle, one can choose three general coordinates: the position of the centre of gravity of the vehicle XC, Y C and the rotation angle ψ of the vehicle around the z-axis (figure 1b). In this study, we consider that the longitudinal velocity v of the vehicle is prescribed and constant; in other words, the vehicle is driven with a constant speed in the vehicle-fixed x-direction. This means that all the points of the beam have the same longitudinal velocity, including the centre of gravity, i.e. vCx=v. The corresponding kinematic constraint reads 2.1where dot denotes differentiation with respect to time.
The equation of motion of this non-holonomic system is determined with the help of the Appell–Gibbs equations. In order to do this (see ), the general velocities are eliminated by intuitively chosen pseudo-velocities, which guarantee the satisfaction of the kinematic constraint (2.1). A possible choice for the pseudo-velocities βk, k=1,2 is 2.2where β1=vCy is the lateral velocity of the centre of gravity in the vehicle-fixed y-direction, and β2 is the angular velocity of the vehicle. Using these definitions, the three general velocities can be expressed uniquely by the two pseudo-velocities in the form 2.3
(a) Tyre model
The front and rear wheels of the bicycle model are described by means of the simplest so-called brush model (figure 2), which assumes that the (linearly) elastic lateral deformation of the tyre contact patch is formed by the lateral displacements of distinct elastic tyre tread elements attached to the rigid central disc of the wheel. The lateral deformation of the front and rear tyres is given by qF(xF,t) and qR(xR,t), respectively; the longitudinal deformations are negligible in this model. Note that the new general coordinates qF(xF,.) and qR(xR,.) are distributed along the vehicle-fixed longitudinal directions, which makes the state space of the system infinite dimensional.
The lengths of the contact patches are 2a for both the front and the rear tyres. The distributed lateral forces acting on the tyre can be calculated by means of the specific lateral stiffness k (N m−2); the damping effect of the tyre is assumed to be negligible here.
In this study, only the rolling of the tyre is investigated, which means that the tyre particles that are in contact with the ground have zero velocities vP=0. These additional kinematic constraints are given by the PDEs 2.4for xF∈[−a,a], xR∈[−a,a] and . Prime denotes differentiation with respect to the space coordinates. The corresponding boundary conditions of the PDEs are 2.5which come from the zero deformation of the leading edges in the two contact patches in accordance with the brush model (figure 2). The transition of the tyre particles in the contact patch is given by their longitudinal velocities relative to the vehicle, 2.6
(b) Governing equations
Without the detailed derivation of the Appell–Gibbs equations (which is given in the electronic supplementary material), we present only the equations of motion for the bicycle model, 2.7which are coupled to the equations of the general velocities (2.3), to the PDEs (2.4) of the rolling condition, including their Dirichlet boundary conditions (2.5), and the longitudinal transition velocities (2.6). On the right-hand sides of (2.7), clearly, the resultant force and moment of the lateral forces acting on the wheel appear with respect to the centre of gravity of the vehicle.
3. Memory effect of contact patches
In the case of rolling, the tyre tread particles in contact with the ground have zero velocity, which means that their positions are constant in the ground-fixed (X,Y) coordinate system during the time of contact. This implies that the tyre lateral deformations travel backwards relative to the vehicle, which can be expressed in the form 3.1for xF∈[−a,a], xR∈[−a,a]. The time delay τF(xF) is the duration needed by a tyre particle of the front tyre to travel backwards relative to the car body from the leading edge at xF=a of the contact patch to its actual position xF at the time instant t; the meaning of the delay τR(xR) is the same with respect to the rear wheel.
The position of a particle in the contact patch is given by 3.2with xF∈[−a,a] for the front tyre, and by 3.3with xR∈[−a,a] for the rear tyre. Using these expressions and the boundary conditions (2.5) in (3.1), the travelling-wave-like solutions of (2.4) assume the forms 3.4and 3.5
While the time delays τF and τR cannot be expressed from these equations in closed forms, the derivatives of the local space coordinates with respect to the time delays can be calculated as 3.6and 3.7By means of these expressions, the integrals in (2.7) with respect to the variables xF and xR can be carried out by the substitution of τF and τR, respectively, and, after the substitution of (3.4) and (3.5), the governing equations are transformed to retarded functional differential equations (RFDEs) and the PDEs (2.4) become decoupled from these governing equations, whereas the RFDEs remain coupled to the ODEs (2.3). All this means that the system is now described by the actual values of the pseudo-velocities β1, β2 and the actual and delayed values of the general coordinates XC, Y C, ψ, whereas the distributed general coordinates qF(x,.) and qR(x,.) are eliminated. Nevertheless, the system will retain its infinite dimensional character, because the RFDE has infinite dimensional state space. Owing to the complexity of the equations, the detailed derivations of the linearized RFDEs are presented only in §4.
4. Linearization at straight stationary rolling
The straight stationary rolling of the vehicle corresponds to the solution 4.1for xF∈[−a,a] and for xR∈[−a,a]. The small-perturbed motion of the vehicle can be described by the linearized form of the governing equations (2.7), (2.3) and (2.4) with (2.6) at the trivial solution.
In the case of small-amplitude oscillations, the first pseudo-velocity simplifies to . This suggests that it is worth dropping the pseudo-velocity notation and using only the general coordinates as state variables, which means that we return to the second-order representation of the equations of motion in order to have them in the most dense mathematical form. Thus, the linearized type of the governing equations forms a system of coupled integro-/partial/ordinary differential equation, 4.2with the boundary conditions (2.5).
Because the derivatives of the space coordinates with respect to the delays are the same for the front and rear tyres in the case of small oscillations, that is, (3.6) and (3.7) result simply in 4.5one can use (4.3) in the first two equations of the linearized governing equations (4.2) when doing the integration by substitution with respect to the new variable τ. The result is a two-dimensional, second-order linear system of RFDEs, 4.6where the general forces are 4.7These two RFDEs uniquely describe the small perturbed motion of the vehicle at its straight stationary running.
5. Dimensionless parameters
By means of a new dimensionless integration variable and dimensionless time 5.1respectively, (4.6) can be transformed to dimensionless form. In order to do this, introduce also the dimensionless constrained longitudinal running speed V , the dimensionless half axle base L, the dimensionless eccentricity E and the dimensionless frequency ratio F as 5.2where 5.3are theoretical natural angular frequencies that would be valid if each d.f. of the non-cruising linearized car model were blocked separately. Figure 3 illustrates the sub-models, which would vibrate by ωI and ωII, respectively.
Finally, using the dimensionless (small) displacement of the centre of gravity 5.4the governing equation for small perturbations is given by two scalar second-order RFDEs, which can be written in matrix form as 5.5
6. Self-excited vibrations
The substitution of the trial solution 6.1leads to the characteristic function 6.2where μ represents the parameter vector μ=[V L E F]T.
The characteristic function at λ=0 can be calculated by its continuous extension as , 6.3which means that, independently of the system parameters, the characteristic root λ=0 always satisfies the characteristic equation. This corresponds to the fact that lateral perturbations of the vehicle cause the deflection of the vehicle from the original direction of the motion.
If there are additional characteristic exponents located in the imaginary axis with non-zero imaginary parts, Hopf bifurcation may occur. In this case, a pair of pure imaginary complex conjugate characteristic exponents λ1,2=±iω satisfy the characteristic equation D(λ1,2,μcr)=0. The corresponding stability boundaries can be determined in the parameter space by means of the D-subdivision method, namely the characteristic exponent λ=iω is substituted into (6.2), and the characteristic function is separated to real and imaginary parts. The solutions of ReD(iω,μcr)=0 and ImD(iω,μcr)=0 characterize the stability boundary curves (or surfaces, or hyper-surfaces) of the straight stationary motion parametrized by the angular frequency ω of the emerging self-excited vibrations. The stable parameter regions among the intricate structure of the boundary curves can be selected, for example, by the methods presented previously [23,24]. In the case of E=0 and for ω=2jπ, the stability boundaries are characterized by 6.4while for and , the asymptotes of the stability boundaries are characterized by 6.5
The linear stability chart of the straight stationary motion of the vehicle is shown in figure 4a, where the stability boundaries are plotted in the plane of the dimensionless speed V and the dimensionless eccentricity E for ω∈[0,8π]. The other dimensionless parameters are realistically chosen by means of the parameters of a medium-sized car (2l=2.5 m, 2a=0.1 m, JC≈ml2; see  for details), i.e. L=25 and F=1. The stable domains are shaded; the encircled numbers show the number of the unstable characteristic roots. The enlarged part of the stability chart is shown in figure 4b, where the intersection of the stability boundaries can be observed.
In order to show the effect of a possible positive/negative eccentricity, the rightmost characteristic exponents are shown in figure 4c for the parameter points E+ and E− of figure 4a. The corresponding lateral vibration modes of the emerging self-excited oscillations together with the characteristic lateral displacements of the front and rear tyre tread elements are also illustrated in figure 4d. As can be observed, if the centre of gravity of the car is closer to the front wheel, the tyre deformations are in phase and there is no nodal point between the wheels in the lateral mode shape. On the contrary, if the centre of gravity is closer to the rear wheel, then the tyre deformations are in the opposite phase, and there is a nodal point between the two wheels, which has zero lateral speed during the corresponding self-excited vibration. Both cases represent, however, the types of snaking motion of the vehicle body in the case of unstable stationary running.
The effect of the dimensionless frequency ratio F is represented in figure 5, where the stability boundaries are plotted for F=0.8 and for F=1.2. The location of the unstable domains can be characterized by formula (6.4) and by F/(2jπ), j=1,2,3….
In this research, the stability of the straight stationary running of a four-wheeled vehicle is investigated when the elasticity of the tyres is taken into account. While the vehicle is described by the simplest bicycle model, and a simple brush model is used for the tyre, the mechanical model introduced here takes into account the instantaneous deformed shape of the tyre tread elements during the motion. This means that a type of memory effect of the tyre–ground contact patch is taken into account, or, in other words, a delayed tyre model is attached to the classical bicycle model of the vehicle.
Even for this simplified dimensional non-holonomic model, the governing equations can be constructed as a set of integro-/partial/ordinary differential equations in an infinite dimensional phase space. This system can be radically simplified when a travelling-wave-like solution of the PDE part of the system is found, which is directly related to the memory effect in the contact patches. The simplest form of the governing equations corresponds to that of a 2 d.f. delayed oscillator.
As is expected in the case of delayed oscillators, many unexplored unstable regions can be identified in the parameter space of the prescribed running speed and the eccentricity of the vehicle's centre of gravity. The essential effects of further mechanical and geometrical parameters of the vehicle are also studied and explained.
These new unstable parameter regions were identified for undamped cases only. It was shown in our recent study  on wheel shimmy that the internal damping of the tyre does reduce and can even eliminate some of the unstable parameter domains. The bicycle model of this study has not been investigated from this view point, but similar behaviour is expected: the unstable domains shrink by increasing damping. However, it must be noted that analogous unstable parameter regions were found in the case of single towed wheels with elastic tyres where the experiments proved the existence of so-called ‘micro-shimmy’  that leads to increased fuel consumption and increased noise level. Similar effects are expected for the case of the four-wheeled vehicles studied here.
In a more general framework, the applied theoretical approach and the corresponding mathematical analysis are examples of the study of multi-scale dynamical problems in which the small-scale dynamics is that in the tiny contact patch and the large-scale dynamics is that of the vehicle itself. Clearly, the existence of the identified unstable parameter regions shows that the small-scale dynamics is relevant for the whole system and it cannot be substituted simply by its approximate stationary solutions. The approach that models the coupling of the two dynamics via DDEs is a powerful tool to study these systems—and this is represented by the applied delayed tyre model in vehicle dynamics.
This research work was supported by the Hungarian National Science Foundation under grant nos. OTKA K101714 and OTKA PD105442. This research is also connected to the scientific programme of the ‘Development of quality-oriented and harmonized R+D+I strategy and functional model at BME’ project. This project is supported by the New Széchenyi Plan (project no. TÁMOP-4.2.1/B-09/1/KMR-2010-0002).
One contribution of 17 to a Theme Issue ‘A celebration of mechanics: from nano to macro’.
- © 2013 The Author(s) Published by the Royal Society. All rights reserved.