Abstract
A classical problem of nonholonomic system dynamics—the motion of a sphere on the inside of a rough vertical cylinder—is extended to rolling friction. The case study is modelled in independent coordinates. Due to the nonholonomic constraints imposed on the sphere, the governing equations arise as a set of differential-algebraic equations. The results of numerical simulations show the transition of the sphere from a sinusoid path on the vertical cylinder surface to a fall with slip. The physics of the ‘paradoxical’ motion is explained in detail.
1. Introduction
The motion of a sphere on the inside of a rough vertical cylinder is one of the classical problems of nonholonomic system dynamics Citation1. The mathematical description of this seemingly trivial case is not simple, however, mainly because of the three-dimensional motion description and the way the nonholonomic constraints due to the slip-less rolling conditions are involved in the dynamic model. Incidentally, modelling of nonholonomic systems has always been recognized as one of the more complex tasks of analytical dynamics, and the situation is made even more difficult by a large variety of proposed formulations for nonholonomic systems Citation1-4. Many of the methods are of historical value only, either old-fashioned or cumbersome for computer applications. Therefore, the projective approach to the modelling and simulation of nonholonomic systems is worth mentioning Citation5,Citation6. Due to the geometrical interpretations, the projection method is conceptually simple and of intuitive nature, and the proposed compact matrix formulation is especially well suited to computer implementation.
The curiosity referred to the considered motion, predicted analytically in Citation1, is the ‘paradox’ that the sphere, when rolled horizontally on the inside of a rough vertical cylinder, does not descend constantly but recurrently goes down and up, moving against gravity when going up. More strictly, when looking on the motion from above (along the cylinder axis), the sphere mass centre moves on a circle with a constant angular speed, and the side-view is so that the mass centre describes a sinusoid path on the cylinder side surface. A similar phenomenon refers to a ‘paradoxical’ motion of a basketball which, when rolling inside the basket rim, often goes up and falls out of the basket. The motion features are confirmed in the sequel through numerical simulations, and an attempt to explain the physics of the ‘paradoxical’ motion is undertaken. To make the simulation more realistic, rolling resistance (often referred to as rolling friction) is considered.
The problem at hand can be modelled using different sets of state variables Citation1. In this paper, we present a formulation in independent coordinates, leading to governing equations in the form of differential-algebraic equations (DAEs) supplemented by nonholonomic constraints. The framework for the mathematical modelling is a projection method developed in Citation6.
2 Mathematical modelling
2.1 Sphere dynamic equations
Consider a homogeneous sphere of mass m and radius a which rolls without slip on the inside of a rough vertical cylinder of radius b ( ). Let us treat in this section the sphere as unconstrained, whose state is described by n = 6 coordinates and n velocities. A possible choice of the state variables is p = [xc yc zc φ θ ψ] T and v = [[xdot] C [ydot] C ż C ωξ ωη ωζ] T , where xc, yc and zc are the coordinates of mass centre C in the inertial reference frame Oxyz, φ, θ and ψ are the classical Euler angles, that is the rotation, nutation and precession angles that describe the angular orientation of a body-fixed reference frame Cζηζ (not shown in ) with respect to Oxyz, and ωζ, ωη and ωη are of the sphere angular velocity components in Cζηζ. The respective dynamic equations of the unconstrained sphere are then Newton – Euler equations which, for a homogeneous sphere featured by J Cζ = J Cη = J Cζ = J = 2ma 2/5 and J Cζη = J Cζζ = J Cηζ = 0, simplify to
2.2 Holonomic constraint and the dynamic equations in independent coordinates
The sphere is subject to m h = 1 holonomic and m nh = 2 nonholonomic constraints. The holonomic constraint follows from the condition that the distance of mass centre C from the axis Oz remains c = b – a, and its equation expressed in the dependent coordinates p (given in implicit form Φ(p) = 0 Citation6,Citation7) is as follows:
For the sphere subject solely to holonomic constraint Equation(3), which can be interpreted as the sphere moving on the ideally smooth inside surface of the cylinder (no friction between the sphere and the cylinder), the dynamic
Equationequation (1)
becomes
The position of the constrained sphere can be described in k = n – m
h = 5 independent coordinates q, and a convenient set is q = [γ zc
φ θ ψ]
T
, where γ is defined in
. The holonomic constraint, given implicitly in Equation(3), can then be expressed explicitly as p = g(q) Citation6,Citation7, and one has Φh(g(q)) ≡ 0. By differentiating p = g(q) with respect to time, the velocity form of the explicit constraint equation follows as v = D(q) [qdot], where D = B
−1(∂g/∂q) and, according to Equation(2)
, B
−1 = diag(I,R). The implicit form C
h
v = 0 of the holonomic constraint at the velocity level is satisfied identically when expressed in terms of [qdot], i.e. C
h
D[qdot] ≡ 0. Thus, it follows that
The holonomic constraint equations, given explicitly at the position, velocity and acceleration levels, respectively, p = g(q), v = D [qdot] and [vdot] = D[qddot] + [Ddot][qdot], are defined by
The orthogonal projection, i.e. the remaining
Equationequations of (9), then yields a relation for determination of λ
h, i.e.
2.3 Nonholonomic constraints and the governing equations
The m nh = 2 nonholonomic constraints follow from the non-slip condition—the condition of vanishing velocity of that sphere point which is in contact with the cylinder. Expressing the condition in the vertical and horizontal directions, the nonholonomic constraint equations can conveniently be introduced as (see also Citation1):
Imposing the nonholonomic constraints Equation(17) on the sphere moving on the inside of the vertical cylinder, the dynamic
Equationequations (10)
become
due to the generalized reaction force
related to q. The Lagrange multipliers λ
nh = [λnh1 λnh2]
T
are the physical constraint reaction forces in the horizontal and vertical directions (
), which are the forces generated by static friction between the sphere and the cylinder.
The governing equations of the sphere rolling without slip on the inside of the vertical cylinder are then 2k + m nh = 12 differential-algebraic equations (DAEs) in k = 5 independent coordinates q, k = 5 velocities ν = [pdot], and m nh = 2 Lagrange multipliers λ nh:
An inconvenience associated with handling DAEs Equation(20) is that they are inherently unstable Citation9, which is usually referred to as a constraint violation problem Citation10. More strictly, due to numerical truncation errors, the original constraint conditions Equation(17)
may be violated by the solution q(t) and [qdot](t) to DAEs Equation(20)
. The problem, especially significant for long simulation times, can be resolved by applying the constraint violation elimination method developed in Citation11. The method (not reported in this contribution) was applied in the numerical simulations, which are reported in the sequel.
2.4 Involvement of rolling resistance
To make the analysis more realistic, rolling resistance can be included as well. The holonomic constraint Equation(3) becomes nonideal in this case, and the relevant constraint reaction λh, when still normal to the side cylinder surface, is moved a distance f (coefficient of rolling resistance) forward in the direction of the velocity of the sphere mass centre C. The consequence is that λh yields a moment about C on the sphere, whose components along the O1, O2 and O3 axes are M
C1 = 0, M
C2 = −λh1
fż/V and
, where
. The effect of the moments is that the generalized applied force vector [fbar]∗ in the governing
Equationequations (20)
must be replaced by [fbar]∗ as follows
3 Numerical simulations
The numerical simulations were carried out for a steel ball of radius α = 2.5 cm, and thus m = 0.517 kg and J = 1.293 × 10−4 kg m2. The cylinder radius was taken as b = 22.5 cm (c = 20 cm). The initial state of motion was rolling of the ball without slip in the horizontal direction (direction of axis 2) with speed V 0 = 4 m s−1, i.e.:
Some features of the obtained results are worthy of mention. First, the high-frequency variations of both φ and θ angles, which are related to the high values of angular velocity of the ball, are intriguing. The physics of the ‘paradoxical’ (sinusoidal) motion of the sphere can then best be explained by analyzing the variations of the angular velocity components and the gyroscopic effects (which cause the unusual motion of the sphere) in the O123 reference frame. By denoting ω = [ωx ωy ωz]
T
and ω′ = [ω1 ω2 ω3]
T
, the second kinematic relation Equation(15) can be written as ω = R′ω′, where R′(γ) is the matrix of rotational transformation between the O123 and Oxyz reference frames. The dynamic equations of the sphere rotational motion about C, expressed in the Oxyz reference frame axes in terms of ω components, are
, where G = diag(J,J,J) and n
C
= [−λnh2
a sin γ λnh2
a cos γ λnh1
a]
T
. The dynamic equations of the rotational motion expressed in the O123 frame axes in terms of ω' components can be obtained as Citation6
, i.e.:
To make the simulation more realistic, rolling resistance was introduced in the way described in Section 2.4, and the coefficient of rolling resistance used was f = 0.02a. The results of numerical simulations are shown in
. The mass centre path on the cylinder surface is still sinusoidal, but, due to the rolling resistance, its mean level descends, and the amplitude grows up. The resistance also causes the ball speed V and its angular velocity ω to decrease in time, and the amplitudes of variations of φ and θ angles tend to grow up. The decrease in speed , which is mainly due to the decrease in
(ż is relatively small and varies according to the sinusoidal motion), results then in the decrease in λh. The diminishing value of the normal force will always cause a situation when |λnh2|, the friction force between the sphere and cylinder, exceeds the value μλh, where μ is the static friction coefficient. Such a situation, for μ = 0.35, is illustrated in . In reality, this will lead to a transition from the slip-less rolling to a downward skidding of the sphere on the cylinder surface. Evidently, the tendency to the transition to fall is dependent on the value of μ and the initial value of V; see .
4 Conclusions
The paradoxical motion of a sphere on the inside of a rough vertical cylinder, predicted analytically in Citation1, has been verified by numerical simulations. The causes for the sphere being pushed up and down during its sinusoidal motion on the cylinder inner side are the gyroscopic torques and the resultant vertical friction force between the sphere and cylinder (a reaction of nonholonomic constraint). The involvement of rolling resistance does not markedly change the character of the motion. Of greater importance for the amplitude of vertical position variations and for the physical realizability of the motion (for keeping the no-slip rolling) is the initial horizontal velocity, . Increased values of V
0 yield diminished vertical position variations, while smaller values of V
0 require larger values of friction coefficient to keep the no-slip conditions and delay the transition to downward fall. The motion features can be validated by a simple experiment, with a mouse ball rolling inside a glass cylindrical vase, for example. The motion has some features in common with the motion of a basketball on the basket rim, though the latter problem should rather be modelled as the motion of a ball on a torus, and the stick-slip phenomena in the contact point seem to play an important role for the features of this motion.
References
- Neimark , J. I. and Fufaev , N. A. 1972 . Dynamics of Nonholonomic Systems , Vol. 33 , Providence, RI : American Mathematical Society . Translations of Mathematical Monographs
- Hamel , G. 1949 . Theoretische Mechanik, Eine einheitliche Einführung in die gesamte Mechanik , Berlin : Springer .
- Papastavridis , J. G. 1998 . A panoramic overview of the principles and equations of motion of advanced engineering dynamics . Applied Mechanics Review , 51 : 239 – 265 .
- Ostovskaya , S. and Angeles , J. 1998 . Nonholonomic systems revisited within the framework of analytical mechanics . Applied Mechanics Review , 51 : 415 – 433 .
- Arczewski , K. and Blajer , W. 1996 . A unified approach to the modelling of holonomic and nonholonomic mechanical systems . Mathematical Modelling of Systems , 2 : 157 – 174 .
- Blajer , W. 2001 . A geometrical interpretation and uniform matrix formulation of multibody system dynamics . Zeitschrift für Angewandte Mathematik und Mechanik , 81 : 247 – 259 .
- Schiehlen , W. 1997 . Multibody system dynamics, roots and perspectives . Multibody System Dynamics , 1 : 149 – 188 .
- Campbell , S. L. and Meyer , C. D. Jr . 1979 . Generalized Inverses of Linear Transformations , London : Pitman .
- Ascher , U. M. and Petzold , L. R. 1998 . Computer Methods for Ordinary Differential Equations and Differential-Algebraic Equations , Philadelphia, PA : SIAM .
- Eich-Soellner , E. and Führer , C. 1998 . Numerical Methods in Multibody Dynamics , Stuttgart : B.G. Teubner .
- Blajer , W. 2002 . Elimination of constraint violation and accuracy aspects in numerical simulation of multibody systems . Multibody System Dynamics , 7 : 265 – 284 .