Nonlinear dynamics of a basketball rolling around the rim

We analyse the dynamics of a basketball which rolls around the rim of a basketball hoop. The rolling steady motions are determined, and we investigate falling, slipping, and instability. The qualitative behaviour of the global dynamics is analysed and the possible trajectories are categorised. We investigate the effect of initial conditions which cause the basketball to fall inside or outside the basket or to remain on the rim.

times, and the observer cannot predict whether it falls inside or outside the basket. Analysis of this type of motion is the motivation of the research work in this paper.
The basketball throw was analysed in [8,19] focusing on the appropriate shooting angle and velocity. Later, the outcome of the different shooting initial conditions was investigated numerically [11,27] including the effects of bounces of the ball on the rim and of air resistance. The different shooting strategies were investigated analytically by Huston and Grau [16]. Several numerical simulations were published by Okubo and Hubbard [20,21,23], focusing on the ball-rim interactions and rebounds of the ball. Other works deal with the identification of the parameters of the basketball itself from measurements [2,22]. The dynamics of the basketball on the rim is similar to that of the golf ball rolling on the edge of the hole [14,15]. The problem of a sphere rolling on different types of fixed [6,7,13,28] or rotating [5,17] surfaces contains many interesting dynamical details related nonholonomic constraints conserved quantities. The dynamics of rolling disk on a flat surface [24] has some structural similarities to that of the basketball, as well.
In the present paper, the main goal is to analyse and explain the long-term rolling of the basketball around the rim by considering the nonlinear dynamics of a rigid body model of the ball-rim system. Previously, Liu et al. [18] discovered the existence of a twoparametric family of steady motions in a similar model. The stability of these steady motions was investigated by the authors in [12].
In the present paper, we provide a thorough analysis of the ball rolling on the rim. The effects of slipping and falling from the rim are included. The steady motions, the global dynamics of the phase space of the system and the physical limitations of the long-term motion are analysed in detail.
The paper is organised as follows. In Sect. 2, the mechanical model is presented and the equations of motion are derived for the basketball rolling and slipping on the rim. In Sect. 3, the steady motions of the rolling ball are found and the physical restrictions of these solutions caused by instability, slipping, falling and limited kinetic energy are determined. In Sect. 4, the global dynamics of the rolling ball is analysed based on the symmetries of the system, and the possible trajectories are categorised. Section 5 focuses on the longterm realisable rolling motions.

Configuration
The basketball is modelled as a rigid sphere of radius r , and the rim is modelled as a rigid torus with major radius R and minor radius a. A sketch of the system can be seen in Fig. 1, and the parameters of the model can be found in Table 1.
During the analysis, we focus on describing the dynamics of the ball which is already in permanent contact with the rim-including both rolling or slipping. The possible collisions preceding this state are not covered by the calculations, but the conditions of separation from the rim are checked from the normal contact force.
Let us consider the orthonormal basis (e 1 , e 2 , e 3 ) fixed to the rim, where the basis vector e 3 is perpendicular to the middle plane of the rim. ) From the geometric centre O of the rim, the location of the contact point C is given by where α denotes the angle around the circumference of the rim and β denotes the angle around the crosssection of the rim (see Fig. 1). Note that at β = 0, the contact point is located in the middle plane of the rim in the inner side. At the contact point C, we define another orthonormal basis (n 1 , n 2 , n 3 ) in the following way: Let n 1 and n 3 be the unit vectors tangent to the parameter lines of α and β of the torus, ∂r OC /∂α ∂r OC /∂α = e 1 cos α + e 2 sin α, (2) n 3 = ∂r OC /∂β ∂r OC /∂β = sin β(sin α e 1 − cos α e 2 ) and let n 2 = n 3 × n 1 = cos β(cos α e 2 − sin α e 1 ) + sin β e 3 (4) be the unit vector normal to the rim surface at the contact point C. The connection between the bases (e 1 , e 2 , e 3 ) and (n 1 , n 2 , n 3 ) can be given in the form n 1 = T e 1 , n 2 = T e 2 , n 3 = T e 3 , where T(α, β) denotes the proper orthogonal tensor defined by (2)- (4). Note, that this tensor can be composed from two rotations in the form T = Q(β, cos αe 1 + sin αe 2 )Q(α, e 3 where Q(θ, e) = I + (1 − cos θ) · e ⊗ e + sin θ · skwt(e) (7) is the rotation tensor about the unit vector e by the angle θ (see [25, p. 224]). In (7), I denotes the identity tensor, ⊗ is the dyadic product of vectors, and skwt(.) gives the skew-symmetric tensor associated to an (axial) vector. The angles α and β can be considered as the first two Euler angles parametrizing the rotation between the bases (e 1 , e 2 , e 3 ) and (n 1 , n 2 , n 3 ). A third rotation is not needed because n 1 must remain in the plane spanned by e 1 and e 2 .
The configuration of the ball is prescribed by the location of contact point C, which is given by α and β, and the orientation of the ball. The orientation can be described by  The values of R, a, r, ρ and m are based on the official basketball rules [1]. The value of the friction coefficient μ is a typical value based on measurements in [2]. These default values are used throughout this paper where P and Q are arbitrary points of the ball in the current configuration, P and Q are the same material points in a given initial configuration, and R is an arbitrary proper orthogonal tensor. The tensor R could be expressed by scalar variables such as Euler angles (see [10, p. 140]), but such parametrization is not necessary for the present calculation. Then, the configuration of the ball is given by the triplet where Q denotes the configuration space composed from the torus S 1 ×S 1 and the special orthogonal group SO(3) of three-dimensional rotations.

Velocity state
During the motion of the ball, the basis (e 1 , e 2 , e 3 ) is fixed, but (n 1 , n 2 , n 3 ) is rotating as the contact point C changes. By differentiating (2)-(4), the time derivatives of these basis vectors are given in the forṁ where × denotes the cross product of vectors, the dot denotes differentiation with respect to time, and = ax(ṪT −1 ) =β·n 1 +α sin β·n 2 +α cos β·n 3 . (11) The function ax(.) denotes the mapping between vectors and skew-symmetric tensors in three dimensions, which satisfies ax(X) × z = Xz for any skewsymmetric tensor X and vector z. The function ax(.) is the inverse of the mapping skwt(.) mentioned above (see [25, p. 497]). Note, that the vector is not the angular velocity of the ball, but it can be considered as the angular velocity of a virtual rigid body co-rotating with the basis (n 1 , n 2 , n 3 ).
The location of the centre G of the ball is given by where the parameter ρ := a +r is the distance between G and the middle circle of the rim (see the left panel of Fig. 1). The velocity of the ball at G can be obtained by the time derivative of (12) considering (10)- (11), The angular velocity of the ball can be defined by the rotation tensor R, The velocity state of the ball is now expressed by the derivativesq = (α,β,Ṙ). Instead, we introduce the vector of five quasi-velocities (see [10, p. 217]), which are chosen to be the components of the angular velocity of the ball, ω = ω 1 n 1 + ω 2 n 2 + ω 3 n 3 (16) and the components of the velocity of the contact point, The ball is rolling on the rim when v C = 0. The derivativesq = (α,β,Ṙ) can be expressed from quasivelocities; from (14), we geṫ and from the rigid body reduction formula v G = v C + ω × r CG with r CG = r n 2 , we can expresṡ The physical meaning of the quasi-velocities (can be seen in Fig. 2) is important because these variables will be used trough the analysis of the paper: -The component ω 3 is called circular angular velocity. In the case of pure circular rolling (all components of s is zero except for ω 3 ), (19)- (20) show that the ball is rolling around the rim with β = const, which is also called the toroidal direction of the torus. -The component ω 1 is called transversal angular velocity. In the case of pure transversal rolling, - (20) show that the ball is rolling around a transversal cross-section of rim with α = const, which is also called the poloidal direction of the torus. -The component ω 2 is called orthogonal angular velocity. In the case of pure orthogonal rotation (called spinning in basketball), (19)- (20) show that the points C and G are fixed and the ball is selfrotating around the axis of the normal contact direction. -The components u 1 and u 3 are the slipping velocities between the ball and the torus, which are the components of the velocity of the contact point C if the ball is slipping on the rim. These components are zero in the case of rolling.
In summary, the configuration of the ball is described by (9) and the velocity state of the ball is described by (15). Thus, the state of the ball is determined by z = (q, s) in the state space Z z. By the usage of quasi-velocities, it is not needed to express the equation of motion of the system as a set of second-order differential equations for q, but instead, it is expressed as a set of first-order differential equations for z. The expression forq is already given in (18)- (20), and the expression forṡ is obtained from the dynamical equations presented below.

Newton-Euler equations
For more complicated systems expressed with the quasi-velocities, the governing equations can be derived algorithmically by using, e.g. the Appell-Gibbs equations (see [10, p. 254]) or the Boltzmann-Hamel equations (see [10, p. 226]). However, for our system, the simple Newton-Euler equations of rigid bodies give results effectively.
Let us consider a rigid body with a mass m and a mass moment of inertia tensor J computed at the centre of gravity G of the body. When the external force system consists of only concentrated forces F 1 . . . F n act- Fig. 2 Description of the velocity state of the ball by quasivelocities. The angular velocity can be separated into transversal (ω 1 ), orthogonal (ω 2 ) and circular (ω 3 ) components by using the coordinate system (n 1 , n 2 , n 3 ). The velocity of the contact point is given by the slipping velocities u 1 and u 3 . In the figure, the components of the angular velocity are denoted by double arrows, and the components of the slipping velocity are denoted by single solid arrows. In the cases when all quasi-velocities are zero expect one of the angular velocity components, the motion of the ball is symbolised by thin dashed arrows and grey silhouettes ing at P 1 . . . P n , respectively, then the Newton-Euler equations in an inertial reference frame can be written in the form where a G is the acceleration of the centre of gravity G and ε is the angular acceleration of the body.
In the case of a basketball, the mass distribution is spherically symmetric. Thus, the moment of inertia has the form J = jmr 2 I independently from the coordinate system, where I is the identity matrix, and j is the dimensionless mass moment of inertia. We assume j = 2/3 because the mass of the basketball is distributed approximately on a thin spherical shell. The spherical symmetry makes gyroscopic term ω × (J · ω) vanish from (21).
The two external forces are gravity force F G acting at G, and the contact force F C acting at C, which are given by F G = −mg e 3 = −mg sin β n 2 − mg cos β n 3 , (22) respectively, where g is the gravitational acceleration, F 2 is the normal contact force and F 1 , F 3 are the components of the friction force. Then, Newton-Euler Eq. (21) of the ball become Note, that Eqs. (24) and (25) are expressed in an inertial reference frame fixed to the rim, while the rotating coordinate system (n 1 , n 2 , n 3 ) is used to express the coordinates of the quantities. Therefore, the derivatives of the basis vectors (10) have to be considered at the calculations.
The acceleration and angular acceleration can be computed from (13) and (14), whereα andβ are given by (19)- (20). In (24)-(25), we neglect dissipation effects such as the torques acting at C (drilling friction, rolling resistance) and the air resistance. Our purpose is to explore and understand the fundamental structure of the dynamics. However, the consequence of dissipation effects are considered in Sect. 5.3.
By taking into account (22) Our eight unknown variables are the derivatives of the five quasi-velocities inṡ and the three components of the contact force F C in (23). The missing two scalar equations are provided by the friction model at the contact point.
Let us apply the Coulomb friction model, where we assume that the coefficients of static and dynamic friction are identical and equal to μ. When the normal contact is ensured, the tangential contact state of the ball can be rolling or slipping. In the rolling case, the velocity v C of the contact point is zero, which is equivalent to In the slipping case, the friction forces are given by Now, Eqs. (24)- (25) and (28) determine all components ofṡ and F C in the rolling case; while Eqs. (24)- (25) and (29) determine these variables in the slipping case.

Equation of motion for rolling
In the rolling case, the rolling constraint (28) directly givesu 1 = 0 andu 3 = 0. Then, the expressions of the contact forces are given by and the derivatives of the angular velocities becomė Now, the dynamics of the variable setŷ = (q, s) is fully determined by (18)- (20), (33)-(35) andu 1 = u 3 = 0. The dynamics of u 1 and u 3 are trivial, and thus, it can be separated from the system. Moreover, the cyclic variables R and α do not appear on the right-hand side of the evolution equations of the other variables.
(The system is symmetric with respect to the change of orientation R of the ball and with respect to the angle α around the rim.) That is, by considerinġ from (20) in the rolling case u 3 = 0, the dynamics of the variables β, ω 1 , ω 2 , and ω 3 is fully determined by the set (33)-(36) of four first-order ordinary differential equations (ODEs). Formally, the rolling dynamics can be written aṡ where is the variable set for rolling in the four-dimensional reduced rolling state space x ∈ X ∼ = S 1 ×R 3 . Note that once the solution for x is known from (37), the solution for the cyclic variables R and α can be obtained from (18)- (19). A rolling state is physically valid if is satisfied. Moreover, the contact condition can be checked to avoid the separation of the ball from the rim, as well. Condition (39) provides a stricter condition than (40). The properties of rolling trajectories satisfying these conditions are analysed later in Sect. 5.

Equations of motion for slipping
In the slipping case, direct calculation from Eqs. (24)- (25) shows that the expression of the normal contact force is given by and the derivatives of the quasi-velocities are ⎡ whereα andβ are given by (19)- (20). The variables R and α do not appear on the right-hand side of the equations, again. Thus, (20) and (42) can be formally written aṡ where is the variable set for slipping in the six-dimensional reduced slipping state space y ∈ Y ∼ = S 1 × R 5 . The condition u 1 = u 3 = 0 determines the fourdimensional subspace X ⊂ Y of the rolling dynamics (37). On this set, the differential Eq. (43) is not defined due to the discontinuity of F 1 and F 3 (see (29)). During the slipping motion, (40) has to be satisfied again. In the case F 2 < 0, the slipping dynamics (43) is not valid, but the ball separates from the rim. After the separation, the motion of the ball is free-fall with possible impacts with the rim.
It can be checked by direct calculation from (43) that in the slipping case, no equilibria y = y 0 exists which satisfies g(y 0 ) = 0.

Steady motions of the rolling ball
In this section, we find and analyse the equilibria of the rolling Eq. (37) in the form x = x 0 satisfying f(x 0 ) = 0. In such states, the state variables x = (β, ω 1 , ω 2 , ω 3 ) are constant in time, but the cyclic variables R and α are changing according to (18) and (19). Thus, the equilibrium point x 0 of the reduced rolling state space X x correspond to special but non-stationary trajectories in the full state space Z z = (q, s). Physically, these solutions correspond to steady motions when the ball is rolling around the rim with a uniform speed at a constant height β (see Fig. 4) while the components of the angular velocity are constant measured in the rotating coordinate system (see Fig. 2). Throughout this paper, the phrase steady motion refers simultaneously to this mechanical meaning and the corresponding equilibrium point of the of differential equationẋ = f(x) in the reduced phase space.
These steady motions of the system were found first in [18]. In this section, we present a throughout analysis of the steady motions, including stability analysis and restrictions by other physical effects. In the next two sections, we show that these solutions play a fundamental role in the global dynamics of the ball and the realisation of the long-time rolling motion on the rim.

Determining of the steady motions
We look for the steady motions of the rolling ball from (37) in the form where x 0 ∈ X and f(x 0 ) = 0. We can find two typical families of equilibria which we call trivial and nontrivial steady motions.

Proposition 1 (Trivial steady motions) A state of the ball given by
is an equilibrium point of (33)-(36) for any ω 20 ∈ R. These solutions are called trivial steady motions.
These two formulae do not correspond to discrete equilibrium points but to two one-parametric families of equilibria in the form (46). The states x + 0 correspond to the case when the ball is placed on the top of the rim with an arbitrary orthogonal angular velocity ω 20 .
The states x − 0 correspond to the case when the ball is 'hanging' on the bottom of the rim, which is physically impossible due to the requirement of the non-negative normal force F 2 . More interesting, non-trivial steady motions exist for the rolling ball: Proposition 2 (Non-trivial steady motions) The nontrivial equilibria of the system (33)-(36) are given bȳ where Proof In (36), the conditionβ = 0 leads to that is, the transversal angular velocity is zero (see Fig. 2). Moreover, ω 10 = 0 satisfies the conditionṡ ω 2 = 0 andω 3 = 0 automatically in (34)-(35). Finally, Formulae (47)-(48) lead to a two-parametric family of steady motions, which form a two-dimensional surface in the four-dimensional state space X . The physical meaning of this family of these steady motions is the following: If we choose an angle β 0 around the minor circle of the rim and choose a rolling angular velocity ω 30 around the rim, as well, then the formula (48) gives the appropriate orthogonal angular velocity ω 20 (β 0 , ω 30 ) which is necessary to maintain this state The expression in (48) diverges if β 0 = ±π/2 or ω 30 = 0 (see the thick solid lines in Fig. 3). The states β 0 = ±π/2 correspond to the top and the bottom position on the ball around the minor radius of the rim (see Fig. 4). The condition ω 30 = 0 corresponds to the state when the ball is not moving around the circumference of the rim. Although ω 20 in (48) diverges along any of these singularity lines, the singularities vanish at the intersection of these lines (denoted by circles in Fig. 3), where we get ω 20 = 0. The point x 0 (π/2, 0) = [π/2, 0, 0, 0] T corresponds to the case when the ball is just placed on the top of the rim without rotation. Atx 0 (−π/2, 0) = [−π/2, 0, 0, 0] T , the ball would 'hang' at the bottom of the rim without rotation. These special solutions coincide with the trivial steady motions x + 0 (0) and x − 0 (0). Considering the periodicity of the angle β, β = −π and β = π correspond to the same physical state. Then, we can see according to Fig. 3, that the singularity lines cos β 0 = 0 and ω 30 = 0 divide the family of non-trivial steady motions into four branches.
There are several physical limitations on the practical occurrence of the steady motions. We analyse three main limitations: instability of the steady motion, slipping of the ball (including falling from the rim) and the limited kinetic energy of the ball.
In the case of the trivial steady motions, there exist simple, general limitations. The state x + 0 is always unstable, which can be proved in a similar way as will be done in Sect. 3.2. (It is interesting, that independently from ω 20 , one of the corresponding eigenvalues has always the positive value √ (g/ρ)/(1 + j).) In the case of x − 0 , the ball is always falling down from the rim, because these motions would require negative normal contact force (F 2 < 0). In the case of the non-trivial steady motions, the situation is more complicated. In the subsequent subsections, we analyse the three limitations (instability, slipping-falling, limited kinetic energy) and determine the subset of the points (β 0 , ω 30 ) for that the non-trivial steady motions are physically realisable.

Limitation from instability
For the linear stability analysis, let us expand (33)-(36) around a chosen (non-trivial) steady motion x 0 (β 0 , ω 30 ), where A(β 0 , ω 30 ) is a 4 × 4 matrix and O 2 denotes the higher order terms. The matrix A can be obtained by direct calculation, and formally, we get where the elements A i j depend on the values β 0 , ω 30 of the given equilibrium point. The four eigenvalues of the system are where Hence, each steady motion has a double zero eigenvalue. It can be checked that the eigenvectors corresponding to this zero eigenvalues span exactly the tangent space ofx 0 (β 0 , ω 30 ) at the given equilibrium. This means that small perturbations in this tangent space just push the system into a different steady motion of the family of the steady motions. Hence, the equilibria are neutrally stable in these directions.
If h(β 0 , ω 30 ) > 0 then the eigenvalues λ 3,4 are real, given by λ 3,4 = ± √ h(β 0 , ω 30 ). Then the equilibrium is a degenerate saddle, and it is unstable. This kind of motion cannot be realized physically, because small perturbations make the system leave the vicinity of the equilibrium.
If h(β 0 , ω 30 ) < 0 then the eigenvalues λ 3,4 are pure imaginary, given by These pure imaginary eigenvalues imply periodic oscillations in the linearized system. As the total mechanical energy of the ball is conserved, the equilibrium is expected to be neutrally stable (stable but not asymptotically stable), and small perturbations create periodic oscillations of the state of the ball. This conjecture will be proven in Proposition 13 from the properties of the global dynamics, and these equilibria prove to be degenerate centre.
The stable and unstable regions in the (β 0 , ω 30 ) plane can be seen in Fig. 5. The typical values of the diagram can be computed analytically using (53). For small values of the circular angular velocity ω 30 , the steady motions are typically unstable in the 'inner' region |β 0 | < π/2 of the rim (see Fig. 4).
Assume that R > ρ (1+ j), which is satisfied for the parameter values of a standard basketball and rim (see Table 1). Then, we can find the following two typical values of the diagram in Fig. 5.
-For large circular angular velocities ω 30 , the unstable region reduces and it tends to (It can be proved that (55) has a single real root in −1 < c β < 1.) The points with β 0 =β are denoted by black dots in Fig. 5. It can be proved thatβ < β * . In the region labelled 'U', the steady motions are unstable. In this case, small perturbations pull the system away from the steady motion. With the data shown in

Limitation from slipping and falling
As it mentioned in (39), the magnitude of tangential contact force is limited by the Coulomb friction model μF 2 ≥ F 2 1 + F 2 3 . By substituting (47) into (30)-(32), the slipping condition (39) becomes If (56) is satisfied then the rolling steady motion is realizable without slipping, which is denoted by the shaded region in Fig. 6 for μ = 1. If (56) is violated then the given steady motion is not realizable and the ball starts slipping immediately. For small values of the circular angular velocity ω 30 , the realizable region tends to |β 0 − π/2| <β = arctan μ. This interval of β is on upper side of the rim (see Fig. 4). For large values of ω 30 , the region tends to |β 0 | <β, which is related to the inner side of the rim (see Fig. 4).

Fig. 6
Boundaries showing the effect of slipping and falling of the non-trivial steady motions (47). In the shaded region, the rolling steady motions occur for μ = 1. Outside this region, these steady motions are not realisable due to slipping. In the figure, the boundary curves are depicted for several values of the friction coefficient μ. In the limit case μ → 0, the region without slipping shrinks to a line. In the limit case μ → ∞, we get the condition of falling of the ball from the rim. With the data shown in Table 1 In the limit case μ → ∞, the expression (56) simplifies to Condition (57) coincides to the requirement (40) of the non-negative normal contact force, which can be calculated from (31) at the steady motion. This curve is shown in Fig. 6, as well. At a steady motion, the ball can remain in contact with the rim only if (57) is satisfied, otherwise, the ball falls from the rim. It can be seen in Fig. 6 that for a finite friction coefficient μ, the slipping condition always provides a stricter condition than the falling condition.

Limitation of maximal kinetic energy
If a steady motion is in the vicinity of the singularity lines β 0 = ±π/2 and ω 30 = 0, the absolute value of the orthogonal angular velocity ω 20 tends to infinity (see (48) and Fig. 3). To restrict ourselves to physically real-isable cases, an upper estimation of the kinetic energy T of the ball can be given by where v throw is the typical maximum speed of the basketball throw. Since the kinetic energy of the basketball is both converted to potential energy and lost due to collisions with the rim and the backboard, (58) is an upper estimation of the energy level of physically realisable steady motions. By substituting (16)- (17) and (47) The boundary curves of these conditions can be seen in Fig. 7 for different values of v throw . Inside this closed curves, the kinetic energy of the steady motion does not exceed the energy of a throw with a speed v throw . Outside the curves, the steady motion would require too much kinetic energy. By reducing the energy level through v throw , the curves in Fig. 7 shrink to the two points which were denoted by circles in Fig. 3. It was mentioned that these states correspond to placing the ball statically onto the top of the rim or 'hanging' it to the bottom of the rim.

Realizable steady motions
In the previous subsections, we analysed different factors which restrict the occurrence of the steady motions of the ball. The trivial steady motions (46) are practically non-realisable because x + 0 is always unstable and x − 0 is unphysical. As for the non-trivial steady motions, let us consider the intersection of the conditions of stability (the shaded region in Fig. 5), the condition of rolling with μ = 1 (Fig. 6) and the maximum kinetic energy corresponding to the shot speed v throw = 10 m/s (Fig. 7). In Fig. 8, we can see the physically realisable steady motions that satisfy all these conditions. For such large throw velocity, we get three typical regions, which are denoted by A, B and C. It can be checked from Fig. 4 that region A corresponds to the inside-above part of the rim, region B corresponds to the inside-below part of the rim and region C is related to the outside-above part of the rim. In real basketball games, we see the long-time rolling of the ball in region A only, which can be explained from the properties of these regions. On one hand, region A has the largest area, that is, it can be reached from a wider range of initial conditions. On the other hand, a lower level of kinetic energy is required for the steady motions in region A. It can be computed numerically that for these parameters, region C exists for v shot 7m/s, region B exists for v shot 5m/s while the existence of region A requires only v shot 1m/s. In these regions in Fig. 8, we denoted by "+" and "−" the sign of ω 20 of the corresponding steady motion. It can be seen that in region A, both signs occur, and on the dotted line, ω 20 = 0. That is, there are steady motions in this region with very low values of the orthogonal angular velocity, which is a further reason why the solutions of region A occur in practice.
In this section, we analysed the steady motions of the ball, where the angular velocity components of the ball remains constant. A long-time rolling solution is not necessarily an exact steady motion, but it can probably be a solution in the vicinity of a steady motion. This is the topic of the following two sections. In Sect. 4, we analyse the global dynamics of the rolling ball with the

Symmetries and global dynamics of the rolling ball
During this section, we analyze the global dynamics of the rolling system (33)-(36). One of the main goals is to determine the possible types of solutions and their relations in the phase space. In this section, the slipping condition (39) and the falling condition (40) are not considered. These effects are considered in the next section.

Symmetries
We can speak about symmetries of the system from different points of view. In this section, we present the important symmetric properties of the rolling equations and use them to obtain information of the global dynamics of the system.
The dynamics of the four phase variables ω 1 , ω 2 , ω 3 and β are invariant under the change of the other variables of the mechanical system. This symmetry has been taken into consideration at the derivation of the differential Eqs. (33)-(36), and the further variables-the angle α and the orientation of the ballare excluded from the phase space.
As the mechanical model does not contain dissipation, the total energy of the system is conserved, that isĖ = 0. This conservative property can be directly checked from the Eqs. (33)-(36). For a given energy level E 0 , the equation E = E 0 determines a smooth 3 dimensional hypersurface in the 4 dimensional phase space. The form of Eqs. (33)-(36) shows that the phase space is invariant under the transformation ω 2 → −ω 2 , ω 3 → −ω 3 , and the direction of the trajectories is unchanged. That is, the phase space has a symmetry under the reflection to the set ω 2 = ω 3 = 0.
It can be shown from Eqs. (33)-(36) that the transformation ω 1 → −ω 1 results in a same phase plane but the opposite direction of the trajectories. We can say that the phase space has a time-reverse symmetry under the reflection to the plane ω 1 = 0.

General periodic solutions
The time-reverse symmetry related to the plane ω 1 = 0 is an essential property of the system. In this subsection, we use this property to prove that periodic solutions are typical in the phase space.

Proposition 3 Consider a solution x(t)
: R → R 4 of the system (37) and let us denote its components by ω 1 (t), ω 2 (t), ω 3 (t) and β(t). Assume that ω 1 (t 1 ) = 0 andω 1 (t 1 ) = 0 for some t 1 . Then, for any t ∈ R, Proof Consider the integration of the solution from the initial condition x(t 1 ) in both forward and backward time. Due to the time-reverse symmetry, the two half-trajectories can be mapped onto each other by the transformation ω 1 → −ω 1 .
That is, if a trajectory intersects the plane ω 1 = 0 then the each piece of the trajectory determines the other. Consequently, if the function ω 1 (t) twice intersects the ω 1 = 0 plane transversely then the solution is periodic: Consider a solution x(t) : R → R 4 of the system (37) and let us denote its components by Then, x(t) is a periodic function with a period 2(t 2 − t 1 ).
Proof It can be shown from (61) that ω 1 (t 1 + k(t 2 − t 1 )) = 0 for any k ∈ Z. That is, there are infinitely many zeroes of ω 1 (t), and these zeroes divide the time line t ∈ R into uniform intervals with a length t 2 − t 1 . Moreover, it follows from (61) that the solution is the same on every second of these intervals.

Non-periodic solutions
In the previous subsection, we showed that all solutions are periodic which has at least two transverse intersection with the plane ω 1 = 0. Non-periodic solutions can exist either if there is a maximum of one zero of ω 1 (t) or, the intersection is not transverse. If the intersection is not transverse at ω 1 (t 1 ) = 0 theṅ ω 1 (t 1 ) = 0. By substituting these values into (33)-(36), we getω 1 =ω 2 =ω 3 =β = 0, which corresponds to the equilibrium solution (45). That is, these trajectories are connected to the equilibrium points, which can be the homoclinic or heteroclinic trajectories of the unstable equilibrium points.
The other possibility for a non-periodic solution is if ω 1 (t) has a single zero or it does not have a zero at all. In both cases, there exists a half-line in the timeline for which ω 1 (t) does not change sign. Consider the case when ω(t) > 0 for t > t 1 , which is assumed in the following proposition. The cases with ω(t) < 0 and t < t 1 can be stated and proved analogously.

Proposition 5 Consider a solution x(t)
: R → R 4 of the system (37) and let us denote its components by ω 1 (t), ω 2 (t), ω 3 (t) and β(t). Assume that there exist t 1 ∈ R such that ω 1 (t) > 0 for any t > t 1 . Then, (33)-(36) leads to Proof Assume that the solution has a given energy level E = E 0 . Then, from (60), the value of ω 1 is limited by where the lower limit comes from the assumption of the Proposition. That is, ω 1 (t) is finite and positive all along t ∈ (t 1 , ∞), and thus, β(t) increases strictly monotonically according to (36). Consequently, the independent variable of (34)-(35) can be changed from t to β, which leads to (62).
Note that even if (62) provides a solution for ω 2 and ω 3 , we must substitute it into (33) and (36) and check if ω 1 is indeed positive for t > t 1 . Otherwise, the solution of (62) is not valid.
Proof By substituting ω 2 = ω 3 ≡ 0 into (33)-(36), we get the systeṁ which is analogous to the dynamics of a single pendulum. It can be shown from (60) that ω 1 (t) does not change sign for E > mgr, which satisfies the conditions of Proposition 5.
System (62) is a linear, homogeneous but nonautonomous system. It can be written into the form dz/dβ = A(β) · z where z = [ω 2 , ω 3 ] T and A(β) is a 2 × 2 matrix periodically depending on β. The solution of this type of system can be determined by  Table 1 is denoted by a circle the Floquet theorem (see [9, p. 55]). First, the monodromy matrix M should be determined such that z(β +2π) = M ·z(β). The eigenvalues of the matrix M (the Floquet multipliers) determine the stability of the solutions. We used numerical simulations to determine the monodromy matrix M and the Floquet multipliers for the physically relevant parameter range of (62). This leads to the following result: Proposition 7 Consider the parameter range 0 < ρ/R < 1 − a/R and 0 < j ≤ 2/3 of the parameters of the system (62). In this parameter range, there is a subset such that the corresponding Floquet multipliers μ 1 and μ 2 of the system satisfy |μ 1 | = |μ 2 | = 1. In this case, the trivial solutions of (62) is neutrally stable and the system exhibits quasi-periodic oscillations. For the complement set of the parameters in this region, the Floquet multipliers satisfy |μ 1 | < 1 < |μ 2 |, the trivial solution of (62) is unstable and the system exhibits exponentially growing oscillating solutions.
The stable and unstable regions can be seen in Fig. 9. In the stable region, the quasi-periodic solution for ω 2 and ω 3 is bounded. Then, for a sufficiently large initial value of ω 1 , the function ω 1 (t) does not change sign (see 33). Thus, the conditions of Proposition 5 are satisfied and the system (33)-(36) exhibits quasi-periodic oscillations, as well.
In the unstable region, the solution for ω 2 and ω 3 is divergent, which makesω 1 and ω 1 divergent due to (33). Then, the condition of Proposition 5 is violated and the solution is not valid for the original system (33)-(36).
It can be seen in Fig. 9 that the parameter values of a standard ball is in the unstable region. Stability would be possible either if the ball is very small (ρ/R is small), or, if the mass of the ball is concentrated close to the centre of the ball ( j is small). For a realistic basketball, we can state the following proposition: Proposition 8 Assume that the parameters ρ, R and j are located in the unstable region denoted in Fig. 9.

Then, the trivial solution of (62) is the only solution which corresponds to a valid solution of the original system (33)-(36).
That is, for these parameter values, the only nonperiodic solutions of the rolling basketball are the solutions of the system (64).

Global dynamics
To summarize the results of this section, let us state the following theorem, which is based on the results presented above According to the data in Table 2, we can try to imagine the structure of the phase space. For the realistic parameter values near to the data of a standard ball, the quasi-periodic solutions do not exist. Then, the only 4D family of the solutions is the set of general periodic solutions. That is, most of phase space is filled by periodic solutions. The only 3D family of solutions is the family of homoclinic or heteroclinic connections of unstable non-trivial equilibria. This family of connections is located on 3D surfaces, which separate the different basins of the periodic solutions. The whole structure of the 4D phase space is thus, similar to that of a planar conservative system with centers, saddles and the heteroclinic or homoclinic connections of saddles. Note that according to the simulations, it seems that all connections are homoclinic, but this has not been proved yet.
As most of the phase space is filled by the 3parametric family of periodic solutions, the levels sets of these solutions should be parameterized-at least, locally-by three independent conservative quantities. The total energy (60) can be chosen to one of these quantities. It seems that we cannot find other conservative quantities expressed in a simple algebraic form. However, we can find the missing two conserved quantities as solutions of differential equations.
Assume that a quantity related to linear or angular momentum is conserved along the motion. In our system, these momenta can be expressed as linear combination of the angular velocities ω 1 , ω 2 and ω 3 . That is, consider the quantity Let us determine the time evolution of (65) from (33)-(36) by direct calculation, which leads to the following proposition:

Proposition 9
The function (65) is conserved along the trajectories of the rolling dynamics (33)-(36) if c 1 (β) ≡ 0, and the functions c 2 (β) and c 3 (β) satisfy Note that (66) has a similar structure to that of (62) with periodic coefficients. From the analysis described above, we can check that for physically relevant parameters, the solutions of (66) are not periodic but exponentially growing oscillating solutions. Surprisingly, this kind of solutions provides the missing two conserved quantities:  Proof As (66) if a linear differential equation, its solutions can be expressed in the form where (β) is a fundamental matrix of (66) (see [9, p. 53]) and c ∈ R 2 . By substituting (67) into (65), we get By choosing two linearly independent vectors c I and c I I , (68) generates two independent conserved quantities P I (β, ω 2 , ω 3 ) and P I I (β, ω 2 , ω 3 ). All other conserved quantity in the form (65) can be expressed as the linear combination of P I and P I I .
Note, that P I and P I I do not depend on ω 1 , thus, the level sets of these quantities are parallel to the ω 1 direction in the phase space. The fundamental matrix (β) depends on j, R, ρ and β, but it cannot be expressed algebraically in terms of these quantities. Thus, it is hard to find the physical meaning of this conserved quantities. We can say based on the form (65) that P I and P I I are somehow related to conservation of linear and angular momenta. Together with the total mechanical energy (60), the global dynamics of the basketball is organized by three conserved quantities. Note that in degenerate limit cases, it is possible to obtain the analytical form of conserved quantities. The following propositions can be proved by direct calculation: Proposition 11 Assume the system (33)-(36) in the limit case R → ∞. Then, the torus shape of the rim becomes a cylinder, and the quantities are conserved along the trajectories of the system. Alternatively, we can obtain the quantities and the conservation of (69) and (70) are equivalent.
In (70), arctan 2 (y, x) denotes the two-parametric arcus tangent function defined by arctan 2 (y, x) = 2 arctan y which gives arctan 2 (y, x) = arctan(y/x) for x > 0. Physically,P I is the scalar product of the angular momentum and angular velocity in the plane of ω 2 and ω 3 , andP I I is somewhat related to the direction of these vectors measured from β. Proposition 12 Assume the system (33)-(36) in the limit case R → 0. Then, the torus shape of the rim becomes a sphere, and the quantities P I = ω 2 , are conserved along the trajectories of the system.
Physically, P I I is the angular momentum of the ball about the vertical axis. In the non-degenerate case of the rim (a torus), the conserved quantities might be somehow combined. The presence of the different frequencies of β with an irrational ratio in (69)and (72) is supposed to cause that the coefficients of conserved quantities (65) do not have a closed algebraic form. Finally, let us state an important consequence of Theorem 1, which can be used in the next section. In Subsection 3.2, we showed that the stable non-trivial equilibria are centers in the linearised system. Now, we can state that these equilibria are surrounded by periodic solutions.

Proposition 13
Consider a non-trivial equilibrium pointx 0 of the system (33)-(36), which is stable with a pair of pure imaginary eigenvalues. Then, there exist a neighbourhood ofx 0 in which all solutions are periodic. Consequently, these equilibria are nonlinear centers and these are neutrally stable.
Proof By checking the properties of the different types of solutions in Theorem 1, it can be showed that the possible solutions surroundingx 0 are either other stable equilibria or general periodic solutions. Even if quasi-periodic solutions exist, they cannot be in a small neighborhood of an equilibrium point: These trajectories require a finite, sufficiently large value of ω 1 , while the equilibria satisfies ω 1 = 0.
During this section, we assumed that all solutions of the rolling ball is realisable without slipping or falling from the rim. However, the analysis of this theoretical system lead to results which are used in the next section to determine the possible long-term rolling solutions in the presence of different mechanical effects.

Long-term realisable rolling motion and falling from the rim
Consider a realisable steady motion of the rolling ball from the regions depicted in Fig. 8. It was shown in Proposition 13 that the steady motions are surrounded by periodic solutions. During these motions, the state variables of the ball (the angle β and the components of the angular velocity) are oscillating around the steady motion. We expect that for small amplitude oscillations of this kind, the oscillating solutions are realisable in the sense of the restricting conditions of slipping and falling. By increasing the amplitude of the oscillations, the slipping and falling of the ball is anticipated, where the long-term motions are no more possible. In this section, our goal is to explore and characterise the regions of initial conditions, where long-term rolling motion of the ball is realisable. Moreover, it can be important to decide whether the motion results in a basketball score when reaching the boundary of these regions.

Effect of falling
First, we assume that the ball can fall from the rim according to (40), but the slipping condition (39) is not considered. This assumption is a reasonable approximation for large values of the friction coefficient μ. At μ → ∞, the slipping condition (39) tends to the falling condition (40).
In a given state (β, ω 1 , ω 2 , ω 3 ) of the rolling ball, the rolling motion can be preserved without falling if the normal contact force satisfies The expression does not depend on ω 2 , thus, F 2 = 0 is a cylinder-like 3D surface in the 4D state space X (β, ω 1 , ω 2 , ω 3 ). From the formula (73), we can have the following conclusions: -In the region −π < β < −π/2 (see below-outside region III in Fig. 4), the normal force F 2 is always negative, thus, the rolling ball cannot enter this region without falling from the rim. -In the below-inside region IV, falling is avoided for large values of the circular angular velocity ω 3 . -In the above-outside region II, the motion without falling is possible if both the circulation ω 3 and rotation ω 1 have small values. -The increase of the transversal angular velocity ω 1 causes falling from the rim in all cases even in the above-inside region I (see Fig. 4).
When we consider a solution beginning with a rolling state, the ball can avoid falling from the rim if throughout the motion, the normal contact force is positive. For a periodic rolling motion with a time period t, this motion is possible without falling if Let us consider the motion of the rolling ball for varying initial conditions. From this point of view,we can find two typical mechanisms for changing the sign of F min 2 . 1. When we consider a family of initial conditions corresponding to structurally similar-periodicsolutions, the value of F min 2 changes smoothly. Then, we can find some special trajectories with F min 2 = 0. The initial conditions of this property form a surface in the phase. We call this surface a deadlock-type boundary of falling. 2. When there is a surface in the phase space that separates regions of structurally different solutions, there can be a sudden change in value of F min 2 while the initial conditions change smoothly. If such boundary modifies the sign of F min 2 , we call it a separatrix-type boundary of falling.
On one side of such surfaces, we can find realisable long-term (periodic) rolling solutions with F min 2 > 0. On the other side, the ball falls down from the rim and the periodic motion is not possible.
When an initial condition is exactly at a deadlocktype boundary surface, the periodic rolling motion is still possible, but the normal contact force F 2 goes down to zero for a moment. As the normal force in (73) depends smoothly on the state variables, the time derivativeḞ 2 should be zero, as well, at this 'deadlock' point.
It was shown in Sect. 4 that the separatrices of the different regions of the phase space consists of the homoclinic or heteroclinic trajectories of non-trivial saddletype equilibria. That is, an initial condition exactly at a separatrix-type boundary surface should tend exponentially to such an equilibrium point.
These scenarios are visualized in Fig. 10 by using only two dimensions. In the four-dimensional phase space of the system, the centre and saddles are located in a two-dimensional surface, the separatrices of the saddles form three-dimensional surfaces, and the condition F 2 = 0 is also a three-dimensional surface.

The algorithm
For the systematic detection of the boundary surfaces in our four-dimensional phase space X (β, ω 1 , ω 2 , ω 3 ), no analytical tool seems to be available. Instead, we perform numerical simulations to find the regions of long-time rolling solutions. The main steps of the simulation are the following: 1. A steady motion of the rolling ball is selected. 2. The initial conditions are varied in a finite vicinity of the steady motion. 3. A numerical simulation is performed from each initial condition. 4. The long-term rolling motion or the falling are detected in each simulation. 5. The boundaries of the resulting regions are checked to be either the deadlock-type or the separatrix-type boundary surface.
The throughout exploration of the phase is beyond the scope of this paper. We selected some steady motions which are realisable from the point of view of stability, falling and slipping. From the region A of Fig. 8, we chose three typical steady motions, which can be seen in Figs. 11, 12 and 13. We take two dimensional sections from the phase space: From the variables of the initial condition, the orthogonal and circular angular velocities (ω 2 and ω 3 ) are fixed and the other two variables (ω 1 and β) are varied. In each diagram, a grid of 360 × 240 points were used as initial conditions for the simulations. The rolling simulations were implemented by Julia language by using a built-in adaptive Runge-Kutta solver [26]. During each simulation, the long-time (periodic) rolling solution can be detected by the repetition of the ω 1 = 0 value twice without falling (see Proposition 4). If falling (negative normal contact force) is detected during the simulation then from a basic calculation of free-fall, we can determine that the ball falls down inside ('in','score') or outside ('out','no score') from the rim.
When the boundaries are found between the longterm (periodic) rolling and the falling (in or out) regions then we can decide from the presented two types: If near the boundary, the trajectories satisfy F min 2 ≈ 0 then we have found a deadlock-type boundary (see above). Or, if the time period of the trajectory tends to infinity then we have found a separatrix-type boundary. The results can be seen in the top-left panels of Figs. 11, 12 and 13. (The remaining panels of these figures are produced from

The results
In Fig. 11, we can see the case which seems to be the most typical in the area A of Fig. 8. The region of periodic rolling motion is bounded by two surfaces, both are separatrix-type boundaries (denoted by thick continuous lines). In the four-dimensional space, these curves are three-dimensional surfaces which surround a tubular region of realisable periodic rolling. Note that by crossing one of these boundaries, the ball falls in (score), and by crossing the other boundary, the ball falls out (no score).
In Fig. 12, the structure of the diagram remains similar, but now, one of the boundaries becomes a deadlocktype boundary of falling (denoted by a dashed line). In Fig. 13, a new region seems to be appear surrounded entirely by a deadlock-type boundary surface. It can be checked from a different section that the two regions of periodic rolling are connected in the four-dimensional phase space.
In all cases, we proved from the simulations that around the steady motions, we can find solutions where the state variables are oscillating periodically and the ball can remain on the rim permanently. These results are modified when the effects of slipping and dissipation are considered.

Effect of slipping
Theoretically, it is straightforward to add slipping to the simulations. During the simulation, we should check the slipping condition (39) and it slipping occurs, we should switch the solution of the slipping differential Eq. (42). However, slipping and Coulomb friction in three dimensions is challenging when we want to simulate both accurately and effectively. (Note that each graph in Figs. 11, 12 and 13 were produced by running the simulation about ≈ 10 5 times).
The main issue during the simulation is the nonsmooth behaviour of the system (42) at u 1 = u 3 = 0, which occurs every time the system switches between rolling to slipping or vice versa. If we can find a realvalued event function in the form E(y) : Y → R on the phase space Y such that E(y) = 0 at the discontinuity, we can use effective event-driven simulation techniques available in many software packages. This methods can be used in many nonsmooth physical systems including the mechanical problems with planar Coulomb friction.
However, the discontinuity u 1 = u 3 = 0 of the rolling basketball cannot be written in the form of simple event functions, which is a usual issue in spatial friction problems. There is discontinuity in the system if the two conditions u 1 = 0 and u 3 = 0 are satisfied at the same time.
Recent studies show [3,4] that in the vicinity of discontinuity set u 1 = u 3 = 0, the trajectories tends to some characteristic directions called limit directions.  By utilizing this property of the vector field, we can create an effective and robust algorithm for the simulation. We present only the essential amount of mathematical background, the details can be found in [3,4].

Nonsmooth dynamics of the basketball
The phase space X (β, ω 1 , ω 2 , ω 3 ) of the rolling ball is a four-dimensional space embedded into the sixdimensional phase space Y (β, ω 1 , ω 2 , ω 3 , u 1 , u 3 ) of the slipping ball. In Y , the set X is a codimension-2 discontinuity set. Our main problem that the slipping vector field (42) is not defined in this set.
However, we can take the limit of the vector field (42) from different directions. At a chosen point x ∈ X , the set of directions in Y normal to X can be characterised by an angle ϕ = arctan(u 3 , u 1 ) ∈ [0, 2π ]. (The situation is similar to the case when an 1D line is embedded into the 3D space and we can measure the direction around the line by an angle.) Consider the last two components of (42). By taking the limit of a vector field at a point x ∈ X from the direction ϕ ∈ [0, 2π), we geṫ where and the star superscript denotes the operation of the directional limit. Let us transform the slipping velocities into polar coordinates in the form where u is the magnitude of the slipping velocity of the ball. The angle ϕ was introduced by a geometric meaning in the phase space, but it can be also imagined physically as the angle of the slipping velocity in the normal contact plane of the ball and the rim. Note that the discontinuity u 1 = u 3 = 0 corresponds to u = 0 in polar coordinates. From the time derivation of (77), we geṫ u =u 1 cos ϕ +u 3 sin ϕ, uφ = −u 1 sin ϕ +u 3 cos ϕ.
Let us take the directional limit of (78) at the discontinuity according to (75). Then, we geṫ From (80), two angles are obtained, which we call limit directions (see [4]). The limit directions give the possible direction of the slipping velocities exactly at the moment of the transition between slipping and rolling. Moreover, it can be checked by substituting back the values ϕ 1 and ϕ 2 into (80) that the limit direction ϕ 2 is always attracting (u * < 0), the trajectories follow this direction at the transition from slipping to rolling. The another limit direction can be either attracting (u * < 0) or repelling (u * > 0), depending on the values of B 1 (x), B 3 (x) and C(x). It can be shown by comparing (76) and (30)-(32) that the boundary between the attracting and repelling cases coincides with the slipping condition (39).

The algorithm
These results give important information about the structure of the trajectories in the plane u 1 − u 3 of the phase space. If both limit directions are attracting then the rolling motion is realisable and the adjacent trajectories tend to the discontinuity along the two limit directions as asymptotes. If the limit direction ϕ 1 becomes repelling then the rolling motion is not realisable and the trajectories leave the discontinuity along the limit direction ϕ 1 .
How can we use these results to create an effective simulation of the slipping-rolling transitions?
1. The simulation of the transition from rolling to slipping is easier. Beginning from a rolling state, we should check the slipping condition (39) with a usual event-driven technique. When slipping is detected at a point x, the dynamics is switched to the slipping differential Eq. (42). The repelling limit direction ϕ 1 can be calculated from (81). Then, the state u 1 = u 3 = 0 should be slightly modified to u 1 = δ cos ϕ 1 , u 3 = δ sin ϕ 1 where δ is a small positive number. Then, the slipping simulation is launched to the appropriate direction. (See the left panel of Fig. 14.) 2. Before a slipping solution turns into rolling, the trajectories approach the set u 1 = u 3 = 0 along one of the limit directions ϕ 1 and ϕ 2 . That is, the angle ϕ = arctan(u 3 , u 1 ) should converge to ϕ 1 or ϕ 2 , which can be detected by using an appropriate tolerance. Note that at the discontinuity set, the leading term of the expansion ofu is not linear but constant, and the trajectories tend to the discontinuity in finite time. That is, when the trajectory gets sufficiently close to the limit direction ϕ i , the time for reaching the rolling state can be extrapolated from u andu * (x, ϕ i ). Then, the dynamics is switched back to the rolling differential equations, and the state is set explicitly to u 1 = u 3 = 0. (See the right panel of Fig. 14.)

The results
By implementing the algorithm presented above, the results can be seen in the top-right panels of Figs. 11, 12 and 13. In Fig. 11, we can see that adding the effect of slipping has almost no effect to the diagram. The boundaries of the periodic rolling comes from the separatrices from the rolling phase space, which are not affected by slipping.
In Fig. 12, the effect of slipping is significant. Instead of the deadlock-type boundary, a similar but stricter condition appears from the slipping of the ball. Moreover, the region of the periodic rolling is now fully embedded to the region where the ball falls out from the rim and no score is achieved. In Fig. 13, we can see the similar effect; moreover, the second region of periodic rolling vanishes.

Effect of dissipation
Up to this point, we neglected the dissipation during the rolling motion of the ball. Now, we investigate the effect of small dissipation and its consequence to the long-term rolling solutions. We can identify two main sources of the dissipation: the combined effect of contact deformations and friction, and air resistance. In this paper, we pay attention to the friction effects.

Modelling the dissipation from the finite size contact
If the bodies are not considered to be completely rigid then in (23), the contact force F C complemented with a contact moment M C acting at the contact point C.
In case of finite contact stiffness, a small finite contact area is initiated around the theoretical contact point.
In the case that we consider rolling motion from the point of view of the rigid body model, some regions of the contact area is sticking and some regions are slipping. The resulting distributed slipping friction can be reduced to a force (a small addition to F C ) and a moment (appearing of a new term M C ). We focus on the qualitative effect of the dissipation and assume the simplest relevant models. That is, the contact moment is considered in the linear form where C d and C r are the coefficients of the drilling resistance and the rolling resistance, respectively. These coefficients are assumed to be small but finite. By adding the term (82) to the right-hand side of the Euler Eq. (25), the rolling dynamics (33)-(35) is modified tȯ The simulation can be performed as we did before, and the rolling dynamics is replaced by the new rolling Eq. (83).

The results
It can be seen that the system (83) does not have equilibrium points any more. That is we expect also from the physical point of view: In the presence of dissipation, the energy loss makes the ball fall from the rim sooner or later. However, if the dissipation is small then the solutions remain close to those without energy loss. By including drilling resistance and rolling resistance one-by-one, we get the bottom panels of Figs. 11, 12 and 13. The bottom-left panels contain the results with C d /( jm R 2 ) = 0.02 · 1/s and C r = 0 (only drilling resistance) and the bottom-right panels correspond to C r /( jm R 2 ) = 0.02 · 1/s and C d = 0 (only rolling resistance).
In Fig. 13, we can see that the region of periodic rolling vanishes. In the region of the former periodic rolling (denoted by thin dotted lines in the figure), the ball falls outside for both types of dissipation.
For the smaller and more realistic values of ω 3 depicted in Figs. 11 and 12, we get a more interesting result. The simulation shows that in the presence of a small amount of drilling resistance, all former periodic solution becomes falling out from the rim. However, the presence of a small amount of rolling resistance pushes the periodic motions to fall into the rim. Note that in both cases, the dissipation is small, and the ball can spend a long time on the rim close to the periodic rolling behaviour.
The detailed effect of the dissipation model is not topic of this paper. However, from simulations from a few example models, it seems that the falling in or out is sensitive to the dissipation models and parameters. The combined effect of drilling and rolling resistance, the nonlinear extensions of the model (82) and the air resistance makes the problem even more complicated. Possibly, even the pattern of the basketball with the regions of different surface roughness modifies the dissipation, which introduces the orientation of the ball as new state variables.
We can conclude that from the models without dissipation, we can determine the regions of initial conditions, where the ball is rolling around the rim for a very long time in the presence of dissipation. However, it still seems to be hard to determine whether the ball falls in or out from the rim in real physical circumstances.

Conclusions
From the analysis of the rigid body model of the ball rolling on the rim, we obtain the following results: -We showed that the phase space of the system contains a two-parameter family of steady motions. Physically, this means that the ball is rolling around the rim with a constant orthogonal and circular angular velocity. These results coincide with those in [18]. -We found that this set of steady motions is divided into four branches by singularity lines where no steady motions exist. -The physical restrictions were determined which decide whether the steady motions are realisable. Namely, we analysed the restrictions from linear instability, slipping, falling, and limited kinetic energy. The most significant realisable set of steady motions coincides with the motions experienced in basketball games.
-Based on the symmetries of the system, we analysed the global dynamics of the phase space of the rolling ball. It was possible to categorise the possible types of solutions. It was shown that for realistic parameters, the general behaviour of the rolling ball corresponds to periodic solutions. -We showed that in addition to the total mechanical energy, two other independent conserved quantities can be found, which cannot be expressed algebraically, but are probably related to the linear or angular momentum of the ball. -We explored the structure of realisable periodic rolling solutions around the steady motions. In the presence of falling and slipping, two types of boundaries were detected where periodic solutions cease to exist. By using numerical simulations, it was determined for changing initial conditions whether the periodic solution is realisable or the ball falls inside or outside the basket. -By assuming a small dissipation, the periodic rolling solutions become long-term rolling solutions, and they fall from the rim after going around several times. It was showed that changing the dissipation models modify sensitively whether the ball falls inside or outside the rim. This explains why it is hard to predict the outcome of the long-term rolling motions of the ball.
The presented results could be extended to different directions. First, a detailed exploration of the initial conditions could be carried out by using the numerical simulations. While some of the possible motions presented above are clearly in good agreement with the basketball dynamics experienced in practice, it would be interesting to reproduce these results in experiments. Moreover, the present analysis follows the behaviour of the ball during its continuous contact with the rim. It would be interesting to analyse the effects of collisions of the ball with the rim before getting to the long-term rolling case. That analysis may restrict the set of initial conditions to those regions which are realisable after a real basketball shot.
It seems that the results of this section are closely related to the properties of other similar systems with rolling bodies. In [24] the analysis of the rolling disk relies on classic results on symmetries and integrability of the system, which make it possible to reduce the dynamics to a two-parametric family of second-order differential equations. The application these results to the case of the rolling basketball would be a topic of further research work.