On the Stability of Digital Position Control with Viscous Damping and Coulomb Friction

In this paper, we investigate the combined effect of viscous damping and Coulomb friction on sampled-data mechanical systems. In these systems, instability can occur due the sampling of the applied discrete-time controller which is compensated by the two different physical dissipation effects. In order to investigate the interplay between these, we focus on how the stable domain of operation is extended by the dry friction compared to viscous damping. We also show that dry friction causes concave envelope vibrations in this extended region. The analytical results, presented in the form of stability charts, are verified by a detailed set of simulations at different representative control parameter values.


Introduction
Positioning is a basic task in automation, where the control system aims to drive a device into a desired position. Applications subjected to position control, especially in industrial robotics, often demand for high accuracy and at the same time fast operation. On the other hand, the necessary performance which can fulfill these requirements may be limited by the presence of friction and also by the digital nature 1 of the applied motion controllers.
For example, positioning accuracy of a proportional controller can be improved by increasing the proportional control gain, but at the same time the system becomes less robust to parameter variations and might get unstable for large gain values, because these gains are limited by the digital realization of the controller [1]. Despite of these facts, robotic systems are often modeled by neglecting friction, and many times continuous-time techniques are used alone to design their motion controllers. From the engineering point of view, these simplifications can be good approximations in a lot of cases, but sometimes more details have to be considered to obtain a representative model. As it was demonstrated by simulations in [2], special vibrations occur with concave envelope when Coulomb friction compensates for the possible instability caused by the sampling. In reference [3], a passivity based approximation was used to determine the stability limit of a haptic application, while in reference [4], a closed form result was presented for the same stability limit by neglecting the harmonics due to sampling. In these two papers, the only considered source of dissipation was Coulomb friction. Based on [5], in robotic systems the power losses that play important role in the dynamics consist of a viscous damping and a Coulomb friction terms in general. This is also confirmed by the measurement results and the presented representative set of system parameters in [4]. In this paper, we focus on how the two basic mechanical dissipation effects, dry friction and viscous damping, modify the core dynamics and stability behavior of sampled-data systems. Although there are several friction models exist with different levels of complexity [6,7], we will show that the simple Coulomb model with viscous damping can explain the dynamic behavior observed in [4,8] when the digital effects of the applied controller are also considered. This paper is organized as follows. In Section 2, the mathematical model of a sampled-data system is introduced where the combined viscous damping and dry friction is considered as the sources of dissipation. The stability analysis of the corresponding discrete-time position controller is presented in Section 3 with comparing the effect of dry friction and viscous damping. The simulation results presented in Section 4 confirm the theoretical predictions, and Section 5 concludes the paper.

Mechanical model
The equation of motion of a dynamic system subjected to discrete-time proportional controller with zero desired position, x d = 0 , where a Coulomb friction model is also considered can be written in the general form where x(t) represents the generalized coordinate as a function of time t corresponding to the modeled degree-of-freedom, and m is the generalized mass that takes its meaning based on the definition of x . In addition, b denotes the coefficient of the generalized viscous damping and f C is the magnitude of the generalized dry friction force. Based on the applied Coulomb model, the friction force can take any values from the interval [ − f C , f C ] at rest. In addition, t j denotes the j th sampling instant, τ is the sampling time and x( t j ) denotes the sampled generalized coordinate at the beginning of the j th time interval according to the considered zero-order-hold digital to analog signal reconstruction.
In order to obtain a compact discrete-time model with reduced number of free parameters, the dimensionless time T = βt is introduced with β = b / m . With this, the dimensionless sampling instant is T j = jβτ = jθ , where θ denotes the dimensionless sampling time and the sampled control input is denoted by x( T j ) = x j . Based on these, the equation of motion in Eq. (1) can be rewritten as where σ = m f C / b 2 , p = m k p / b 2 and the prime denotes differentiation with respect to the dimensionless time.
By assuming that the direction of motion does not change in a given sampling period, the piecewise linear system of Eq. (2) can be solved in closed form. Then, by defining the discrete state vector x j = [ x j x j ' ] T and considering also the velocity reversals, the following switched discrete mapping can be derived with ε = e −θ − 1 .

Stability analysis 3.1 Neglecting dry friction
First, the case is investigated when the effect of dry friction is neglected. The corresponding results can serve as reference for analyzing the stabilization effect of dry friction. It follows, that Eq. (3) simplifies to a linear map x j+1 = A x j .
Using this linear map, the evolution of the discrete state vector can be determined for any initial conditions which results in a multi-dimensional geometric series [1,9]. Therefore, the system is asymptotically stable if, and onlyif, the magnitude of the complex eigenvalues z n = exp ( s k θ) of the state transition matrix A are less than one, i.e., z n s k n k < , = , ⇔ ( )< , = , , ...
It has to be noted, that this condition can be analytically further analyzed by applying the Mobius transformation z = (w + 1 ) / (w − 1) and using the Routh-Hurwitz stability criteria [1,9].
The corresponding stable domain of control parameters is illustrated in the plane of dimensionless sampling time θ and dimensionless proportional gain p in Fig. 1. It is noted, when the effect of viscous damping is also neglected, the desired position x d is always unstable [4]. (3) In Fig. 1, the yellow region with the solid black contour represents the stable domain of the control parameters. The hatched regions correspond to different cases when z 1,2 ∈R . The gray hatching lines indicate when both eigenvalues are positive, i.e, 0 < Re( z 1,2 ) < 1 . In this case the system shows firstorder dynamics; no vibrations will develop. In the narrow upper region hatched in red −1 < Re(z 1 ) < 0 or −1 < Re(z 2 ) < 0. Here, the control force alternates, and at the boundary of stability period doubling vibrations may occur. In the solid yellow region, the characteristic multipliers z 1,2 ∈C with Im(z 1,2 ) ≠ 0 and the system oscillates similar to an underdamped secondorder system.

Approximation of dry friction
Using the method of describing function (DF) analysis, the effect of dry friction can be approximated. The application of the method results in an effective, frequency and amplitude dependent transfer function in general, called the describing function [10,11].
Focusing only on the approximation of the effect of dry friction, the viscous damping term is neglected in Eq. (1) which gives where u denotes the control force. This system is represented in the form of a block-diagram at the top of Fig. 2. This 2π periodic square-wave function can be expanded into a Fourier series as π ω ω π ω ω π π π π π π C C C d d cos c os sin ( ( ) .
n Neglecting the harmonics, i.e., the cases when n ≥ 2 , the fundamental friction force component is The corresponding block-diagram is shown in the lower part of As it can be seen, the describing function N(V) operates as a velocity dependent effective viscous damping coefficient b C . When the motion is initiated with the following conditions then the amplitude of the velocity function becomes V = v 0 . Therefore, for the stability analysis, the same method can be used which was applied in Section 1, but with higher damping is the inherent viscous damping and b C is an effective quantity, originating from Coulomb friction.

Combined friction and viscous damping
The dynamic behavior of the system presented in Eq.
Note that the coefficients of the characteristic polynomial, here, were directly determined by using the matrix sum determinant identity published in [12].
Based on Eq. (14), the characteristic multipliers z 1,2 can be determined as When ρ > 1 , there would be exponentially unstable oscillations around the reference position without the damping effect of dry friction. Therefore, in this case, the motion can be characterized as a damped oscillator with negative viscous damping term which models the unstable dynamic behavior. By also neglecting the harmonics due to sampling, a continuoustime oscillator model can be derived in the form where parameter f 0 = f C / k p describes the same sticking region as that of the original sampled-data system ∆ = f C / k p . For the solution of Eq. (18), ζ is considered to be in the range 0 < ζ < 1 , and the initial conditions are selected as x 0 > 0 and v 0 = 0 . With these, the solution is until the first velocity reversal, which happens at t = π / ω d with

Verification of results
In order to verify the results presented above, a rotating disk is considered with zero desired angular position x d driven by a discrete-time proportional controller where the continuoustime driving torque is realized by a zero-order-hold. The system parameters collected in Table 1 are the effective model parameters of the experimental setup presented in [4]. By using the parameters collected in Table 1, the corresponding stable domains of control parameters are illustrated in the plane of sampling time τ and proportional gain k p in Fig. 4.  In this figure, based on the methods presented in Section 3.1 and Section 3.3, the yellow and the green regions represent the stable domain of the control parameters, neglecting and considering the effect of dry friction. The gray hatched region represents the stability region using the approximation presented in Section 3.2. In our simulation study, we are focusing on systems with sufficiently small sampling times providing that the frictionless system can only experience a Neimark-Sacker type bifurcation at the corresponding limit of stability. This is the domain which is practically the most important. Fig. 4 shows the stable domains calculated by using describing function analysis in hatched gray. The solid green area represent the stable domain obtained with the analytical method presented in Section 3.3. Clearly, the describing function based analysis predicts the stability limit reasonably well, but it does not give a conservative estimate for the practically more important small sampling times.
The difference between the two results are verified by simulations presented in Fig. 5, where the blue and green time histories correspond to the cases when the nonlinear dry friction and the effective viscous damping were considered in addition to the inherent viscous damping, respectively.
The presented simulation results correspond to the control parameters pairs at point A, B, and C in Fig. 4. Among these, point A is located near to the intersecting stability boundaries determined by the exact formula in Eq. (23) and by the stability condition in Eq. (24) that uses the effective viscous damping. From Fig. 5 it is clear that the simulated system is stable in both cases, but the time histories qualitatively do not match. The other two simulation results at point B and C show that close to the stability limit the effective viscous damping based method gives false predictions. The simulation was implemented in Mathematica by using the command NDSolveValue with event handlers for the detection of velocity reversals and for the realization of a uniform sampling. The representative τ − k p simulation parameter pairs are collected in Table 2. Typical dynamic behaviors when the analyzed system is stable and damped by both viscous damping and Coulomb friction are presented in Fig. 6. The top chart correspond to point D in Fig. 4 which point is selected from the boundary of stability without friction. In this case viscous damping compensates for the destabilizing effect of sampling alone and the developing vibrations have linearly decaying amplitudes. The other two charts labeled by E and F show vibrations where the Coulomb friction stabilizes the otherwise unstable motion and where the two different damping effects strengthen each other. Similar to charts A and C in Fig. 5, chart E in Fig. 6 also presents vibrations with a concave envelope. Here, this is less visible due to the weaker instability caused by sampling (see Fig. 4).

Conclusion
In this study, the stabilization effect of combined viscous damping and Coulomb friction was investigated on the dynamics of discrete-time controlled mechanical systems. It was shown, that using the method of describing functions, the stability of sampled-data system can successfully predicted in a wide range of the control parameters, but close to the stability limit this prediction is not reliable. The simulation results also show that time histories obtained with the effective viscous damping do not qualitatively match the time histories of the original non-linear model. When Coulomb friction stabilizes the otherwise unstable motion, the damped vibrations have concave envelope, and decay much faster than the vibrations that correspond to the effective viscous damping model of Coulomb friction.