Application of a special class of smooth models of the resultant friction force and moment occurring on a circular contact area

Some special cases of a larger class of smooth models of the resultant friction force occurring on a planar contact area are developed and presented. These models are constructed under assumptions of classical Coulomb’s law of friction on each element of contact and instantaneous transition between two modes: fully developed sliding and rigid stick state without local slips and deformations of the contact area. They are able to model cases of different values of static and kinetic friction factors. The considered models are very effective tools for fast and reliable numerical simulations of mechanical systems with frictional contacts. They are applied in modelling and numerical analysis of a special kind of mechanical system, i.e. a kind of a pendulum with special frictional driving. The rotational joint of the pendulum is elastically suspended in the motion plane. The pendulum is driven by circular frictional contact with a rotating disk. Examples of self-excited bifurcation and chaotic dynamics as well as stick–slip behaviour of the pendulum are presented.


Introduction
There are plenty of examples of mechanical systems, in dynamical properties and behaviour of which friction plays a crucial role. Friction can be a desired phenomenon or not, but in both cases, its appropriate and efficient mathematical modelling is a key part of analysis and synthesis of mechanical systems with frictional contacts. The classically understood friction model is a relation between a single component of the friction force and the one-dimensional relative displacement of the contacting bodies. This relation can have different levels of complexity, beginning with the classical Coulomb's law of friction and ending with advanced relations, where additional state variables are often defined. These kinds of models can be directly applied in a mathematical description and analysis of dynamical systems with frictional contacts, where at each element of the contact, the same relative motion of the contacting surfaces occurs. However, in real life, one can encounter many examples of mechanical systems, where the above assumption cannot yield correct results. As examples, dynamics of rolling bearings, billiard balls, different kinds of tops, the wobblestone, polishing machine, disk clutches and many others can be given. Exact and correct results can always be obtained by detailed physical modelling and space discretisation in the vicinity of contact. Nonetheless, this approach leads to a computation cost increase and is inappropriate in fast numerical simulations. This is the reason of interest of many researchers in looking for simple approximate models of contact forces, which would be suitable for fast and realistic simulations of certain classes of mechanical systems with frictional contacts.
Contensou [1] has noticed that if the product of the normal component of relative angular velocity of the contacting bodies and the size of contact is sufficiently large, one should take the coupling between the friction force and moment into account. He has proposed an integral model of the resultant friction force under the assumption of fully developed sliding and Coulomb's law of friction valid on each element of the circular contact. The results obtained by Contensou were then significantly developed by Zhuravlev [2], who presented an exact analytical solution to Contensou's integral model and also proposed the corresponding approximant models based on the Padé expansions. Further developments and generalisations of the approximant models of the contact forces, including rolling resistance and assuming elliptical contact area, are proposed in the work [3]. An integral model of the friction force and moment was constructed under the assumption of classical Coulomb's law and fully developed sliding over the planar contact zone of an arbitrary shape and general contact stress distribution. In the next step, a special model of contact distribution over the elliptic area was constructed, which can be understood as a kind of generalisation of Hertzian normal stress distribution, where a special distortion related to the rolling resistance was modelled. Certain original approximate models of the resultant friction force and moment, which are based on Padé approximants and their generalisations, have been presented in [3]. Special regularisations of these models, which enable to avoid singularities for vanishing relative motion as well as to take into account different values of static and kinetic friction coefficients, can be also found in [4,5].
The topic of the present work is a part of the larger field concerning modelling and simulation of multibody systems with frictional contacts. The modelling methods dealing with those kinds of problems can be divided into two main branches, i.e. non-smooth models and regularised models [6]. In the regularised approaches, known also as penalty or complaint methods, the contact force is assumed as a continuous function of penetration of contacting bodies [7][8][9]. The frictional forces are then often formulated by the use of regularised Coulomb's model of friction [10][11][12][13]. The regularisation can be a result of a need of simplification of the model but also very often corresponds to the physical nature of a real object (compliance of the bodies, microslips, etc.).
Stamm and Fidlin [14] have proposed a regularised two-dimensional model of friction forces appearing on a finite plane area based on elasto-visco-plastic theory, but requiring discretisation of the contact area. The authors have applied their model in modelling and analysis of a disk-on-a-disk system being, to a certain extent, a counterpart of a disk clutch, where alternating sliding and sticking solutions can occur [15]. The work [16] presents results of analytical studies of a similar system, where the approximations based on Taylor's expansion of the friction force and moment for fully developed sliding are used.
In the present work, the authors apply their earlier developed models of the resultant friction force and moment in modelling and numerical simulations of the mechanical system being a certain modification of the disk-on-a-disk system analysed in the works [15,16]. The present work is an extended version of the conference publication [17]. In comparison with the work [17], analysis of solution's accuracy as a function of the "smoothing parameter" ε is included.

Modelling of friction forces
Let us consider a circular contact area of dimensionless form depicted in Fig. 1, where (i) the fully developed sliding; (ii) classical Coulomb's law of friction valid on each element of the contact; (iii) constant friction coefficient; (iv) the fact that the relative motion of the contacting surfaces is a planar motion of rigid bodies is assumed. Aside from the aforementioned, also the following relations between real and dimensionless quantities describing the contact are assumed: y), v s = av s and ω s = ω s , where T s and M s are the real resultant friction force and moment of contact reduced to the centre A of contact F, T s and M s -the corresponding dimensionless resultant friction force and moment, μ-friction coefficient, N -real normal loading of contact, a-characteristic dimension of contact (in this case, the real radius of the contact F), σ and σ -real and dimensionless contact pressures, v s and v s -relative real and dimensionless relative linear velocities at the point A, ω s = ω s -relative angular velocity in the plane of contact. Note that in all the introduced quantities, time is real.
Assuming that T s = −T sx e x − T sy e y , M s = −M s e z , v s = v sx e x + v sy e y and ω s = ω s e z , where e ξ is the unit vector along the axis ξ , one can find the following integral expressions: where Since the model (1) requires integration over the contact area at each time step, it is not a suitable tool for fast and reliable numerical simulations. Based on the results presented in the previous works of the authors [3] and assuming a constant dimensionless contact pressure distribution on the circular contact area σ (x, y) = 1/π, one can derive the following two sets of the integral model (1) approximations and The approximations (2-3), after the replacements v sx = v s cosϕ s and v sy = v s sinϕ s , fulfil the following properties of the integral model (1) where f = T sx , T sy , M s , while m and b are arbitrary constants. Note that the above approximations possess the same denominators, which, in general, is not necessary (see [3]), but allows for application of the later presented special form of regularisation.
In order to make comparison of the functions (1-3) simpler, one introduces the spherical coordinates v sx = λ s cos θ s cos ϕ s , v sy = λ s cos θ s sin ϕ s , ω s = λ s sin θ s .
The  Fig. 2. In the further modelling process, the approximations (3) will be used.

Mathematical model of the pendulum
In Fig. 4, a physical conception of the special mechanical system, being a certain modification of the diskon-a-disk system analysed in the works [7,8], is presented. A physical pendulum of a mass M and a moment of inertia B with respect to the mass centre C is rotationally connected with the use of the joint A to a light The governing equations of the presented mechanical system read where X C and Y C are the coordinates of the mass centre C; ϕ-angular position of the pendulum; e = AC-distance between the point A and the mass centre C; k X and k Y -stiffness coefficients of the elements supporting the rotational joint A; c X and c Y -the corresponding damping coefficients; c ϕ -the coefficient of damping in the rotational joint A; μ-kinetic friction coefficient; N -normal loading of contact; X S and Y S -coordinates of the point S; g-gravitational acceleration. As a model of the resultant friction force and moment, the relations (6) are applied. The kinematic arguments of the functions (6) read In some of the figures, we present motion of the point A, the coordinates of which read

Numerical simulations
In all the numerical simulations presented in this section, the following parameters are fixed: g = 9.81m/s 2 , b = 0.452, m = 0.765 and ε = 10 −3 . In Fig. 5a, a bifurcation diagram of the system exhibited in Fig. 4, with constant angular velocity ω 0 playing a role of a bifurcational parameter, is presented.   In Fig. 9, the results of analysis of solution's accuracy as a function of the parameter ε are depicted. In order to make two solutions comparable in the time domain, one introduced the synchronising element, i.e. periodically varying angular velocity of the driving disk ω 0 (t) = cos 1 rad s t rad s . Other parameters are the same as in the case of analysis presented in Fig. 8. The corresponding stick-slip solution, obtained for ε = 10 −4 , is presented in Fig. 8a. Then, the difference between the same solution obtained with ε = 10 −3 (b), 10 −3 , (c), 10 −1 (d) and the solution for ε = 10 −4 treated as a pattern is shown in the remaining subplots of Fig. 9.
For all the numerical computations presented in this study, Wolfram Mathematica 10.3 package has been used. For numerical solution to the differential equations, the method NDSolve with the default settings has been implemented.

Conclusions
In the present work, examples of models of the resultant friction force and moment based on previous works of the authors have been presented. They are simple functions, which can be an effective substitute for the exact integral model, suitable for fast and realistic computer simulations of a certain class of mechanical systems with frictional contacts.
These models, in their primary form, concern the case of a fully developed sliding on the contact area and possess singularity for the case of lack of the relative motion. The applied regularisation occurred to be an effective method to avoid that problem and take into account different values of static and kinetic friction coefficients. The drawbacks include stiff differential equations and the change in physical properties of the system near the stick mode.
In order to test the developed models of friction forces, a mathematical model of a special mechanical system is built, which is some modification of the disk-on-a-disk system analysed in the works [15,16], being for ω 0 = 61.5 rad/s a strongly simplified disk clutch. It is expected that the proposed model can exhibit much richer bifurcational dynamics, allowing for testing different aspects of friction models. It should be noted that the presented work is in progress, and only preliminary results are reported. Future construction of the corresponding experimental rig is also considered.
One has to pay his attention to the problem of accuracy of approximation (3). Since there is a difference between this approximation and original integral model (1) seen in Fig. 2, one has to expect the corresponding difference between the bifurcation dynamics (observed in Fig. 5) of the pendulum with the friction force and moment described by these two models. Since numerical simulations with the full integral model (1) are very inconvenient and time-consuming, we do not present this comparison on bifurcation diagrams. On the other hand, one has to bear in mind that the main goal of modelling and numerical simulations is to predict and explain dynamics of real objects. Since both the integral friction model and its approximation can describe the real friction phenomenon most approximately, the differences between the models (1) and (3) become less important. Their usefulness can be eventually validated by experiments. Note that a more general class of these models was successfully verified during experimental investigations of the wobblestone dynamics [18].