Dynamic response analysis for multi-degrees-of-freedom parallel mechanisms with various types of three-dimensional clearance joints

Abstract

For spatial multibody systems, the dynamic equations of multibody systems with compound clearance joints have a high level of nonlinearity. The coupling between different types of clearance joints may lead to abundant dynamic behavior. At present, the dynamic response analysis of the spatial parallel mechanism considering the three-dimensional (3D) compound clearance joint has not been reported. This work proposes a modeling method to investigate the influence of the 3D compound clearance joint on the dynamics characteristics of the spatial parallel mechanism. For this purpose, 3D kinematic models of spherical clearance joint and revolute joint with radial and axial clearances are derived. Contact force is described as normal contact and tangential friction and later introduced into the nonlinear dynamics model, which is established by the Lagrange multiplier technique and Jacobian of constraint matrix. The influences of compound clearance joint and initial misalignment of bearing axes on the system are analyzed. Furthermore, validation of dynamics model is evaluated by ADAMS and Newton–Euler method. This work provides an essential theoretical basis for studying the influences of 3D clearance joints on dynamic responses and nonlinear behavior of parallel mechanisms.

Keywords

Spatial parallel mechanism dynamics model spatial clearance joints Lagrange multiplier

Introduction

Robots are becoming an indispensable part of modern industry and are playing an increasingly important role.^1,2 For some reasons, such as manufacturing tolerance, assembling errors, and wear phenomenon of the kinematic joints, clearance always exists in kinematics joints, and the performance of the robot is unsatisfactory.^3
–5 Planarization of spatial clearance joints has limitations in studying the dynamics behavior for mechanical systems because the complex motion mechanism of spatial clearance joints is neglected. Moreover, frequent contact and collision events are the source of nonlinearity. It may exhibit chaotic responses under certain conditions and even affect the overall dynamic behaviors of the system.

Over the past decades, issues related to clearance joint in mechanism dynamics have been a hot topic.⁶ It is noteworthy that most of the research were mainly associated with the issue of plane system. However, dynamics research of multi-degrees-of-freedom (DOF) spatial systems considering 3D compound clearance joint has not been reported as far as we know. Yu et al.⁷ eliminated redundant constraints of an overconstrained mechanism and further investigated the dynamics model and performances of the system considering multiple revolute clearance joints. Zhan et al.⁸ modeled joint clearance parameters as interval variables and provided fresh ideas for studying motion reliability of 3-RRR parallel manipulators considering joint clearances. Li et al.⁹ built a rigid-flexible coupling dynamic model to analyze the dynamics behavior of solar array. A rigid frame and flexible parts are described by nodal coordinate formulation and absolute nodal coordinate formulation frames, respectively. Zheng et al.¹⁰ developed a modeling technique based on nonsmooth dynamic method and have studied the effect of stiction on the dynamic response. Zhao et al.¹¹ presented a modeling approach of a piston skirt–liner system with lubricated pairs and evaluated lubrication and dynamic performances of the mechanism. It is worth mentioning that the above work takes into account different factors and provides ideas for establishing an accurate dynamic model with clearance joints. Different from the previous method, Ordiz et al.¹² used Simscape to model the industrial press. Then, ANSYS was used for stress analysis to explore the impact of clearance evolution on service life of the machines. Wu et al.¹³ studied the nonlinear characteristics of a crank-slider mechanism with multiple clearance joints through correlation dimensions. Ambrósio and Pombo¹⁴ and Marques et al.¹⁵ improved pre-existing dynamic modeling methods of the mechanisms containing clearance joint, respectively. Both of these methods simplify the modeling process and improve modeling efficiency.

The above literature works only for planar systems. However, it is necessary to consider 3D models of the system in some cases. Erkaya¹⁶ investigated the influence of spherical clearance joints and flexible joint connection on the vibration of linkage by experiment and numerical simulation. Chen and Sun¹⁷ established a dynamics model of parallel mechanism containing spherical clearance pair. The clearance model in the above article has been extended to three dimensions, and the established dynamic model is more accurate than the plane model. However, the effect of revolute clearance pair is ignored. Using nonsmooth modeling method, Cavalieri and Cardona¹⁸ developed a 3D model of revolute pair, in which axial clearance is neglected. Yan et al.¹⁹ validated the existence of 3D revolute pair with clearance through experiments. They verified correctness of the proposed model, which can describe the axial and radial clearance. Tian et al.²⁰ investigated the influences of lubrication cylindrical joint and flexibility component on dynamics responses of crank-slider system. Then, Liu et al.²¹ developed a four-point model to detect the contact between cylinder elements. Flores et al.²² deducted a kinematic model for revolute clearance joints. They carried out numerical simulations on the double pendulum mechanism. A new enhanced formulation was recently developed by Marques et al.²³ to simulate axial and radial contact of the revolute clearance joint. The simulation results show that it is necessary to consider radial and axial clearances in dynamic modeling process. Subsequently, Marques et al.²⁴ considered different friction models and explored the dynamic response of the RSSR mechanism. As can be seen from the above literature, many scholars have developed different 3D clearance joint frames for application in different fields and have achieved good results.

For multibody systems, to reveal dynamics behavior of the system more accurately and truthfully, various types of clearance joints should be considered in the dynamic model. In this article, a motion frame of 3D revolute clearance joint is proposed. A dynamics modeling method of parallel mechanism considering 3D compound clearance joint is given. The influence of the compound clearance joint on the dynamic response of the spatial parallel mechanism is discussed. Furthermore, evaluating influences of clearance joint on output responses are all prerequisites for analysis, design, optimization, and control of nonlinear systems. The purposes of the article are to analyze the influences of 3D compound clearance joint on dynamics responses of parallel mechanisms.

The structure of the remainder of this article is as follows: The second section provides a mathematical description of the 3D revolute clearance joint and 3D spherical clearance joint, including kinematics model and contact force model. Nonlinear dynamics model of the parallel mechanism with 3D revolute clearance joint and spherical clearance joint is presented in the third section. A computational procedure for the dynamic model is also presented in this section. In the fourth section, influences of compound clearance joint and axis misalignment on system response are explored, and dynamics model is verified by ADAMS and Newton–Euler method, respectively. Salient conclusions are presented in the fifth section.

Model for 3D clearances joints of parallel mechanism

Mathematical description for 3D spherical clearance joint

Introduction of clearance into a spherical pair released three kinematic restrictions, which can be shown in Figure 1. It is generally considered that only one contact form, point contact, exists in dynamic simulation modeling because of its structural reasons. The geometric centers of the two elements (socket and ball) of the spherical joint are expressed as s_i and b_j . Vectors r _si and r _bj are positions of s_i and b_j in global coordinate frame. Two potential contact points are denoted as P_si and P_bj . Radii of ball and socket can be written as R_bj and R_si .

Figure 1.

The 3D spherical clearance joint.

The clearance of the spherical joint can be described as

c = R_{s i} - R_{b j}

The vector of eccentricity e is defined by

e = {[e_{x} e_{y} e_{z}]}^{T} = r_{b j} - r_{s i}

Relative penetration depth is written as

δ = \sqrt{e^{T} e} - c

Mathematical description for 3D revolute clearance joint

For idea revolute pair, the axis of journal and bearing always keep in coincidence during joint movement. However, when the robot is connected by the joint, clearance inevitably exists. In that case, axis misalignment may occur in the clearance joint. This results in different joint element configurations, namely, no contact, point contact, line contact, and plane contact. The above contact conditions can be considered as a combination of independent axial and radial contact conditions. The schematic diagram of four contact configurations can be found in other scholars’ work.^19,22,23 A 3D revolute clearance joint is shown in Figure 2, and it consists of bearing and journal. The radius of the journal is R_jo . R_b ₁ and R_b ₂ are the inner and outer radii of bearing, respectively. L_j and L_b denote the length of journal and bearing, respectively.

Figure 2.

The 3D revolute clearance joint.

The axial and radial clearances can be defined as

\{\begin{cases} c_{a} = (L_{j} - L_{b}) / 2 \\ c_{r} = R_{b 1} - R_{j o} \end{cases}

For dealing with point contact, the contact impact detection of 3D revolute clearance joint can be regarded as the contact impact detection problem of two impact rings located at the top and bottom of bearing. Therefore, as shown in Figure 3, four potential contact points are located in the top impact ring, where P_b and P_j are used to detect contact between bearing top surface and journal flange, M_b and M_j are used to detect the impact between bearing inner wall and journal. Similarly, the four potential contact points corresponding to the bottom impact ring are Q_b and Q_j , N_b and N_j , respectively. Centers of journal and bearing are o_j and o_b , respectively. J_T and B_T are the centers of the top surface of journal and bearing, respectively. Unit vectors e _j and e _b in global coordinate frame are aligned with the axis of journal and bearing.

Figure 3.

Radial contact diagram: (a) Top impact ring, (b) bottom impact ring, and (c) top view.

Position vectors of J_T and B_T can be given by

\{\begin{cases} r_{JT} = r_{o j} + (L_{j} / 2) e_{j} \\ r_{BT} = r_{o b} + (L_{b} / 2) e_{b} \end{cases}

where r _oj and r _ob are position vectors of o_j and o_b in global coordinate frame.

The angle between journal flange and bearing is ζ, and the angle between vector e _j and bearing is ξ. Construct a plane consisting of three points J_T , B_T , and J_B , the unit vector of normal direction of the plane is $n_{1} = e_{j} × (r_{JT} - r_{BT}) / e_{j} × (r_{JT} - r_{BT})$ . Normal contact direction is derived as $n_{M} = n_{1} × e_{j}$ . The projection of vector e _b on vector n ₁ is Z ₁ point, and its position vector can be given by

r_{z 1 b} = {(r_{B T} - r_{o b})}^{T} n_{1} n_{1}

The location vector of contact point M_b of bearing can be written as

r_{M b} = r_{BT} + R_{b 1} n_{1} \times \frac{(r_{BT} - r_{z 1 b})}{| | (r_{BT} - r_{z 1} b) | |}

Relative penetration depth $δ_{M}$ caused by radial impact is evaluated as

δ_{M} = R_{j} - n_{M}^{T} (r_{M b} - r_{JT})

The position of potential contact point M_j can be obtained by combining equations (7) and (8)

r_{M j} = r_{M b} + δ_{M} n_{M}

For the opposite sides of contact state, the analysis process is similar.

Figure 4 shows the schematic of axial contact state. It should be noted that potential contact points of axial contact and radial contact are not necessarily in the same plane. However, the position of potential axial contact is always related to e _b . Moving vector e _j to the starting point of vector e _b , construct a plane that contains both vectors e _j and e _b . Normal vector of the plane is expressed by $n_{3} = e_{j} \times e_{b} / e_{j} \times e_{b}$ . Normal contact direction of points P_b and P_j can be written as $n_{p} = - e_{j}$ .

Figure 4.

Axial contact diagram.

The potential contact point at bearing top surface and the corresponding potential contact point on the journal flange are P_b and P_j , respectively. Their corresponding position vectors can be written as

\{\begin{cases} r_{P b} = r_{B T} + R_{b 2} (e_{b} \times n_{3}) \\ r_{P j} = r_{p b} - δ_{p} e_{j} \end{cases}

where $δ_{p}$ , caused by contact between bearing top surface and journal flange, is the relative penetration depth and can be given by

δ_{p} = e_{j}^{T} (r_{P b} - r_{J T})

The analysis of contact states occurring in the bottom surface of bearing and journal flange is similar to the analysis above.

The above analysis is about point contact cases. When bearing and journal are parallel to each other, ζ = ξ = 0, the line contact state may occur in the radial direction, and the plane contact state may occur in the axial direction. This is a particular case of the above deduction.

Calculation of contact force for clearances joints

With the motion of clearance joint, contact-impact events frequently occur at clearance joints. An appropriate contact force model is a vital part of restoring the dynamic behavior of kinematic pair and predicting global dynamics behaviors of the mechanism.

L-N continuous contact force model based on Hertz contact law is widely employed to describe normal contact force of contact-impact events. Dissipation term is introduced to reflect the energy dissipation in the contact process.^9,19,21 In this article, L-N model²⁵ is employed to describe normal contact force of clearance joints in point contact state

\{\begin{matrix} F_{k}^{N} = K δ_{k}^{n} + H δ_{k}^{n} {\dot{δ}}_{k} & δ_{k} \geq 0 & (k = M, N, P, Q) \\ F_{k}^{N} = 0 & δ_{k} < 0 \end{matrix}

where the first term $K δ_{k}^{n}$ is estimated by Hertz contact model and used to describe elastic forces, the second term $H δ_{k}^{n} {\dot{δ}}_{k}$ is used to assess energy dissipation. The value of n is 3/2. $δ_{k}$ represents the relative penetration depth at each contact point. ${\dot{δ}}_{k}$ is relative normal contact velocity of two contact bodies.

For the radial contact between two cylinders in a parallel state, based on Hertz pressure distribution, nonlinear models such as Johnson model, Radzimovsky model, and Goldsmith model have been developed.²⁶ These models are implicit functions that require numerical iteration to evaluate contact forces without explaining energy dissipation. It is noteworthy that point contact and line contact on the revolute joint does not have the same contact stiffness as in a spherical clearance joint in a real scenario, since they have different topologies. To simplify the calculation, this work adopts an approximate strategy. Hence, equation (12) is used to describe the relationship between force and penetration of line contact state.^23,27 The exponent n ranges from 1 to 3/2. A more general relation, equation (13), is used to describe normal contact force $F_{A}^{N}$ under plane contact state

\{\begin{matrix} F_{A}^{N} = \frac{E_{j o} δ_{S} π (R_{b 2}^{2} - R_{b 1}^{2})}{L_{b}} & δ_{s} \geq 0 \\ F_{A}^{N} = 0 & δ_{s} < 0 \end{matrix}

where $E_{j o}$ and L_b are Young’s modulus and length of bearing $δ_{S} = δ_{p} or δ_{S} = δ_{Q}$ .

Friction is a commonplace event in which the direction is relative to the direction of motion of two objects. Many friction models were developed to describe the friction phenomena reasonably. Coulomb friction model, one of the most classical models, is expressed as a function of normal force and relative tangential velocity. However, some drawbacks also limit the application of the model, such as discontinuity around zero velocity, neglecting Stribeck effect and viscous characteristics. After that, some bristle-based friction models have been reported.^28

–31 These models can describe some complex phenomena such as microsliding but usually lead to inefficiency and introduce additional parameters that are difficult to identify. In this work, a continuous velocity-based friction model³² is utilized to assess tangential friction force at the contact surface of clearance joint. This model can describe static, dynamic, and viscous friction phenomena without increasing computational time.^33,34 It has continuity and differentiability and is suitable for control and multibody system simulation, which can be written as

\begin{array}{l} F_{m}^{T} = F_{m}^{N} μ_{k} tanh (4 v_{t} / v_{s}) + \\ F_{m}^{N} (μ_{s} - μ_{k}) v_{t} / v_{s} {(\frac{1}{4} {(\frac{v_{t}}{v_{s}})}^{2} + \frac{3}{4})}^{2} m = (k, R, A) \end{array}

where v_t denotes relative tangential velocity between contact elements, v_s is transition velocity, μ_s and μ_k are the static and dynamic friction coefficients.

Contact force and moment of potential contact point acting on bearing are given by

\begin{array}{l} F_{m} = F_{m}^{N} n_{m} + F_{m}^{T} t_{m} \\ M_{m} = r d_{m} \times F_{m} m = (k, R, A) \end{array}

where t _m is the tangential unit vector of potential contact point m and rd _m is the radius vector of potential contact points m. Spherical joint only appears point contact state, contact force and moment can be calculated through equations (12), (14), and (15).

Detection of the instant of contact is correlated with the penetration depth of two discrete integration steps t_n and t_n+1 . It has also been reported in other literature.^9,17,21 Once the collision is detected, contact force should be added to motion equation of multibody system as an external force. Here, note that the normal force does not produce torque. The friction torque acting on the clearance joint should be calculated according to equation (15).

Dynamics model of parallel mechanism with two type clearance joints

Figure 5 shows a 4-UPS/RPS parallel mechanism sketch with 5DOF, including a fixed platform, a moving platform (end effector), and five variable-length kinematic chains. Five kinematic chains can be divided into two categories: UPS (U is universal joint, P represents prismatic joint, and S represents spherical joint) kinematic chains and RPS (R is the revolute joint) kinematic chain. In each kinematic chain, the prismatic joint situated at $P_{i} = (i = 1, \dots,5)$ can be actuated. Five kinematic chains are linked at moving platform through five spherical joints located at $S_{i} = (i = 1, \dots,5)$ . The centers of revolute joint and universal joints are named R ₁ and $U_{i} = (i = 1, \dots,5)$ , respectively. In this model, the revolute joint R and spherical joint S ₁ distributed on the RPS kinematic chain are regarded as imperfect joints. The geometrical and mass properties of the 4-UPS/RPS parallel mechanism are presented in Table 1.

Figure 5.

Kinematic sketch of 4-UPS/RPS parallel mechanism.

Table 1.

Geometrical and mass properties of 4-UPS/RPS parallel mechanism.

Component	1	2	3	4	5	6	7	8	9	10	11
Length (m)	l_s					l_t					—
Length (m)	0.76					0.84					—
Mass (kg)	10.9	10.54				8.57					18.64
I_xx (kg.m²)	0.009	0.006				0.002					0.458
I_yy (kg.m²)	0.684	0.643				0.549					0.237
I_zz (kg.m²)	0.681	0.643				0.549					0.237
Joint	R ₁	U ₂	U ₃	U ₄	U ₅	S ₁	S ₂	S ₃	S ₄	S ₅	—
Distribution radius (m)	r ₁	r ₂				r ₃					-
Distribution radius (m)	0.717	0.645				0.202					-
Distribution angles (rad)	$ϕ_{1}$	$ϕ_{2}$	$ϕ_{3}$	$ϕ_{4}$	$ϕ_{5}$	$ϕ_{1}$	$ϕ_{2}$	$ϕ_{3}$	$ϕ_{4}$	$ϕ_{5}$	—
Distribution angles (rad)	0	π/4	3π/4	5π/4	7π/4	0	2π/5	4π/5	6π/5	8π/5	—

As illustrated in Figure 5, global coordinate frame $O_{A} - X_{A} Y_{A} Z_{A}$ is attached at mass center of the fixed platform. Furthermore, all the vectors mentioned later are represented in this framework. Similarly, the local coordinate frames $o_{i} - x_{i} y_{i} z_{i} (i = 1, \dots,11)$ are attached at the mass center of active components, respectively.

The position and orientation of active components are parameterized by variable

\begin{array}{l} q = {[q_{1} q_{2} \dots q_{11}]}^{T}, \\ q_{i} = [x_{i} y_{i} z_{i} α_{i} β_{i} γ_{i}] (i = 1, \dots, 11) \end{array}

where ${[\begin{matrix} x_{i} & y_{i} & z_{i} \end{matrix}]}^{T}$ are position vectors of the origins of local coordinate frames in global coordinate frame. $α_{i}$ , $β_{i}$ , and $γ_{i}$ are Euler angles.

Rotation matrix from local coordinate frames $o_{i} - x_{i} y_{i} z_{i} (i = 1, \dots,11)$ to global coordinate frame $O_{A} - X_{A} Y_{A} Z_{A}$ can be defined as

\begin{array}{l} R_{i} (α_{i}, β_{i}, γ_{i}) = R_{i} (X, α_{i}) R_{i} (Z, β_{i}) R_{i} (Y, γ_{i}) \\ = [\begin{matrix} c β_{i} c γ_{i} & - s β_{i} & c β_{i} s γ_{i} \\ c α_{i} s β_{i} c γ_{i} + s α_{i} s γ_{i} & c α_{i} c β_{i} & c α_{i} s β_{i} s γ_{i} - s α_{i} c γ_{i} \\ s α_{i} s β_{i} c γ_{i} - c α_{i} s γ_{i} & s α_{i} c β_{i} & s α_{i} s β_{i} s γ_{i} + c α_{i} c γ_{i} \end{matrix}] \end{array}

Universal joint allows relative rotation along perpendicular axes. The kinematic constraint can be written in the form

\begin{array}{l} Φ_{U k} = [\begin{array}{l} o_{i}^{u} - r_{j} - R_{j} o_{j}^{u} \\ {(R o_{i})}^{T} \cdot R o_{j} \end{array}] = 0_{4 \times 1} \\ (k = i = 2, \dots, 5; j = 7, \dots, 10) \end{array}

where r _j denotes position vectors of origins of local coordinate frames in global coordinate frame. R _j is the rotation matrix from local coordinate frames to global coordinate frame. $o_{i}^{u} = r_{2} {[0 cos (ϕ_{i}) sin (ϕ_{i})]}^{T}$ and $o_{j}^{u} = {[- l_{s} / 200]}^{T}$ are position vectors of U_i in global coordinate frame and local coordinate frame, respectively. Ro _i is a unit vector along z_j direction in global coordinate frame. $R o_{j} = R_{j} {[0 1 0]}^{T}$ is unit vector along y_j direction.

Prismatic joint, as the actuated joint, allows relative translational along the joint axis. In this mechanism, the constraint equation is defined as

\begin{array}{l} Φ_{P k} = [\begin{array}{l} p_{i} \times q_{j} \\ p_{i} \times d_{n} \\ f_{i}^{T} \cdot f_{j} \end{array}] = 0_{5 \times 1} \\ (k = i = 1, \dots, 5; j = 6, \dots, 10) \end{array}

where the unit vectors p _i and q _j are aligned with joint axis of x_i and x_j . d _n is a unit vector from the origin o_i to o_j , f _i and f _j are two unit vectors along axes z_i and y_j , respectively.

The kinematic constraint of the perfect revolute joint can be written as

Φ_{R} = [\begin{array}{l} o^{r} - r_{6} - R_{6} o_{6}^{r} \\ {(R o_{r})}^{T} \cdot R o_{j 1} \\ {(R o_{r})}^{T} \cdot R o_{j 2} \end{array}] = 0_{5 \times 1}

where o ^r denotes position vectors of the revolute joint in global coordinate frame. Ro _r is a unit vector and can be written as $R o_{r} = {[0 0 1]}^{T}$ , $R o_{j 2} = R_{j} {[0 1 0]}^{T}$ , and $R o_{j 1} = R_{j} {[1 0 0]}^{T}$ .

Kinematic constraints of the perfect spherical joint can be expressed as

\begin{array}{l} Φ_{S k} = r_{11} + R_{11} o_{i}^{s} - r_{j} - R_{j} o_{j}^{s} = 0_{3 \times 1} \\ (k = i = j = 2, \dots, 5) \end{array}

where r ₁₁ is the position vector of o ₁₁ in global coordinate frame. R ₁₁ is the rotation matrix from $o_{11} - x_{11} y_{11} z_{11}$ to global coordinate frame $O_{A} - X_{A} Y_{A} Z_{A}$ . $o_{i}^{s} = r_{3} {[0cos (ϕ_{i}) sin (ϕ_{i})]}^{T}$ and $o_{j}^{s} = {[l_{t} / 2 0 0]}^{T}$ are position vectors of S_i in local coordinate frames $o_{11} - x_{11} y_{11} z_{11}$ and $o_{j} - x_{j} y_{j} z_{j}$ , respectively.

Driving constraint equation of five variable length kinematic chains can be given by

\begin{array}{l} Φ_{D k} = [| | r_{i} - r_{j} | | - f_{i} (t)] = 0_{5 \times 1} \\ (k = i = 1, \dots, 5; j = 6, \dots, 10) \end{array}

where r _i denotes the position vectors of origins of local coordinate frames in global coordinate frame and f_i (t) is the preset law of motion for each actuated joint.

In this work, the revolute joint (R ₁) and spherical joint (S ₁)are assumed to be clearance joints, and then, they impose no kinematic constraints. Hence, they do not contribute to the overall constraint equation of the mechanism. Combining equations (18) to (22), one obtains an expression for constraint equations of the system as

Φ = {[\begin{matrix} Φ_{U 1} \dots Φ_{U 4} & Φ_{P 1} \dots Φ_{P 5} \\ Φ_{S 2} \dots Φ_{S 5} & Φ_{D 1} \dots Φ_{D 5} \end{matrix}]}^{T} = 0_{58 \times 1}

As the numerical iteration progresses, the accumulation of numerical errors will lead to violation of constraint and even nonconvergence of results. Based on the principle of feedback control, Baumgarte’s approach³⁵ controls the numerical error by modifying acceleration constraint equations. Dynamics equation for 4-UPS/RPS parallel mechanism with both 3D revolute clearance joint and 3D spherical clearance joint utilizing Baumgarte stabilization method can be written in matrix form as

[\begin{matrix} M & Φ_{q}^{T} \\ Φ_{q} & 0 \end{matrix}] [\begin{array}{l} \ddot{q} \\ λ \end{array}] = [\begin{array}{l} Q \\ γ - 2 α_{1} \dot{Φ} - β_{1}^{2} Φ \end{array}]

where M is a 66 × 66diagonal matrix, representing system mass matrix. Φ _q is Jacobian matrix of constraint equations Φ. q includes generalized acceleration vectors of the mechanism. γ represents right sides of acceleration constraint equations. λ contains Lagrange multipliers associated with constraints. Q represents generalized external forces, such as gravity vector and contact forces. α ₁ and β ₁ are Baumgarte feedback parameters.

The computational methodologies for solving the dynamic equations of the mechanism considering both the 3D revolute clearance joint and 3D spherical clearance joint are shown in Figure 6.

Figure 6.

The computational strategy for dynamic equations.

Dynamics response analysis of parallel mechanisms with two types of clearance joints

In this section, the numerical simulation results of the mechanism considering the 3D revolute clearance joint and 3D spherical clearance joint are analyzed. Simulation parameters of clearance joints are presented in Table 2. Initial positions and velocities of each component are deduced by kinematics analysis of this mechanism. The effect of initial axis misalignment of the 3D revolute clearance joint on dynamic response is assessed by comparing the two cases. In the first case, axis of journal and bearing coincides, and in the second case, the axis misalignment occurs at the initial moment, that is, the bearing rotated 0.035° clockwise around coordinate axis x ₆. The motion trajectory of end actuator is defined as (unit: s, m)

\{\begin{cases} X = 1.21 - 0.03 cos (π t) \\ Y = 0.03 sin (π t) \end{cases}

Table 2.

Simulation parameters of clearance joints.

Parameter	Values	Parameter	Values
Ball radius R_sj (m)	0.0226	Length of bearing L_b (m)	0.12
Socket radius R_bj (m)	0.0228	Length of journal L_j (m)	0.118
Bearing outer radius R_b ₂ (m)	0.0175	Young’s modulus of socket and bearing (GPa)	90
Bearing inner radius R_b ₁ (m)	0.0125	Young’s modulus of ball and journal (GPa)	207
Poisson ratio of socket and bearing	0.32	Poisson ratio of ball and journal	0.27
Integration step(s)	10⁻³	Restitution coefficient c_e	0.9
Integration tolerance	10⁻⁶	Baumgarte-α ₁, β ₁	5
Integrator scheme	Ode15s	Transition velocity v_s (m/s)	0.001
Dynamic friction coefficient	0.1	Static friction coefficient	0.15

As shown in Figures 7 and 8, the dynamic behavior at the 3D spherical and revolute clearance joints is given in preparation for the following analysis of global dynamic behavior of the parallel mechanism. Figure 7(a) and (b) shows center trajectory diagram of spherical joint and the corresponding time history diagram. The difference in the initial conditions (axis alignment and axis misalignment) caused the center trajectories not to coincide before 0.515 s. The system gradually stabilizes as it enters the continuous contact state (the trajectories are both concentrated in the rectangular area in Figure 7(a)). From the second cycle, the time history diagram of the two trajectories basically coincides. Besides, in the case of initial alignment, the spherical joint does not move in the Z direction, as shown in Figure 7(b). Figure 7(c) shows center trajectories of bearing relative to journal with initial axis alignment. Combining with Figure 7(d), it can be observed that center trajectory of top and bottom coincides basically. Figure 7(e) and (f) are center trajectory and time history diagrams in the case of initial axis misalignment. Unlike initial axis alignment, there are apparent differences between the top and bottom center trajectories, which continue until 0.515 s, especially in the Y direction. Contact forces of clearance joints are shown in Figure 8(a) and (b). Considering misalignment of axis, contact force of high amplitude appears in the initial stage. When the system enters the stable stage, gravity plays a dominant role. The contact forces of clearance joints are in a continuous contact state. In the stable stage, the collision time of the revolute is comparatively consistent with that of the spherical joint, the two clearance joints are at the upper and lower ends of the same kinematic branch, and the peak contact force of revolute joint is less than that of spherical joint. Driving force curves for prismatic actuators considering initial axis alignment and misalignment are shown in Figure 8(c) and (d). The maximum values of driving force are related to the transient response. When considering axis misalignment, the peak values of the driving forces are higher than the case of axis alignment. Affected by the contact impact, oscillation positions of driving force curves have been aligned to the contact force curves.

Figure 7.

Dynamic behavior at clearance joint. (a) Center trajectory in spherical joint, (b) time history diagram in spherical joint, (c) center trajectory in revolute joint with initial axis alignment, (d) time history diagram in revolute joint with initial axis alignment, (e) center trajectory in revolute joint with initial axis alignment, and (f) time history diagram in revolute joint with initial axis misalignment.

Figure 8.

Contact force and driving force. (a) Contact force at spherical joint, (b) contact force at revolute joint, (c) driving force with initial axis alignment, and (d) driving force with initial axis misalignment.

The enormous contact force at the initial stage is related to the initial conditions. For multibody dynamics systems with clearance joint, the initial stage response depends on the initial configuration and is transient. It can be observed that the contact force generated by the initial misalignment of the bearing axis is much greater than the initial alignment. The collision phenomenon described based on the penalty method takes extremely short. The forced displacement caused by the driving will cause a larger contact force. After the system enters a stable state, the contact force level returned to normal. Therefore, in dynamic analysis, the response after the initial transient stage is often more concerned. For example, in the literature,^16,36
–38 the initial transient response is filtered.

Next, ADAMS was used to verify the dynamic model. The virtual prototype model is imported into ADAMS and then materials are assigned to each component. Fixed platform is connected to the ground by a fixed joint. Four universal joints are added to the fixed platform according to the joint distribution angle. As shown in Figure 5, the ideal spherical joints and prismatic joints are added to the moving platform and kinematic branch chains, respectively. The joint motions of prismatic joint are selected as translational. Create a cylindrical journal and a cylindrical bearing at the revolute clearance joint and add them to the fixed platform and the RP branch, respectively. The contact type is defined as cylinder to cylinder. Similarly, two balls with different radii are added to the PS branch chain and moving platform. The contact type is defined as sphere to sphere. Normal force and tangential force are defined as impact and Coulomb, and set the same simulation parameters as in the article. From the analysis of Figures 7 and 8, results indicated that initial axis misalignment only greatly impacts the initial stage. Therefore, axis alignment is considered in ADAMS model. From the observations in Figure 9(a) and (b), the compound clearance joint effect on the vibrational levels of displacement curves is inconspicuous. However, as shown in Figure 9(c) to (f), the velocity and acceleration curves of the moving platform fluctuate significantly. These fluctuations occur almost simultaneously with contact at clearance joints.

Figure 9.

Dynamic responses of end effectors: (a) X-displacement, (b) Y-displacement, (c) X-velocity, (d) Y-velocity, (e) X-acceleration, and (f) Y-acceleration.

Compared with numerical simulation results, the curve of ADAMS has little difference in fluctuation degree and amplitude, and the vibration is more severe in Y direction, as displayed in Figure 9(c) and (d). From the acceleration curves (Figure 9(e) and (f)), it can be seen that there are several large fluctuations due to them is an alignment of the axis, which is also supported by contact force curve. The contact impact at the clearance joint is transmitted to the end effector, which can cause a high-frequency response of acceleration, and further to induce velocity fluctuations. It ultimately leads to a decrease in positioning accuracy. In short, clearance will bring a series of adverse effects, such as high-frequency vibration and low positioning accuracy. In this article, the dynamic model of 4-UPS/RPS parallel mechanism is established based on the Lagrange multiplier technology, and the influence of the compound clearance joint on the dynamic response of the system is further analyzed. The constraint equation library of corresponding types of joints can be established through the joint constraint equations given in the article, which can be easily applied to this type of parallel robot. For redundant constraints, LU factorization and SVD factorization techniques can be used to identify and delete redundant constraints in the system.

The dynamic model of the 4-UPS/RPS mechanism established by the constraint joint library using Lagrangian multiplier technology is compared with our previous work (Newton–Euler method).³⁹ Comparison curve of displacement response is shown in Figure 10. The maximum displacement errors of the two methods in the X and Y directions are 2.28 × 10⁻⁴ m and 9.5 × 10⁻⁴ m, respectively, as shown in Figure 10(a) and (b). Coefficient of determination (R ²) is used to evaluate agreement observed in the results and $R^{2} = \sum (y_{i} - \bar{y}) / \sum ({\hat{y}}_{i} - \bar{y})$ . y_i and ${\hat{y}}_{i}$ are the i’th data calculated by the Lagrange multiplier technique and Newton–Euler method, respectively, and $\bar{y}$ is the mean value of the data obtained by the Lagrange multiplier technique. As shown in Figure 10(a) and (b), R ² of the displacement in the X and Y directions of the model established by the Newton–Euler method and Lagrange multiplier technique is 0.999988 and 0.999952, respectively. A similar method was used to obtain R ² of the ADAMS axis alignment result and the Lagrange multiplier technique result. As shown in Figure 9(a) and (b), R ² of the displacement in the X and Y directions are 0.999963 and 0.999927, respectively. The above results show that the agreement between the three methods is excellent. The results of Lagrange multiplier technique and Newton–Euler method are closer. However, the modeling and programming of Newton–Euler method is more complicated. According to the established constraint joint library, it can be directly integrated into the constraint equation of the corresponding mechanism. Although the Lagrangian multiplier introduction increases the computational cost, the modeling has better generalizability and easy to code in a computer program. Compound clearance joint effect on the displacement of end effector is inconspicuous, but it dramatically affects the motion state of the joints. At the same time, the existence of the clearance joint has caused a series of adverse effects such as noise, high-frequency vibration, and wear of the joint. This work has a guiding effect on the structure design and performance optimization of the mechanism.

Figure 10.

Comparison of results between Newton–Euler method and Lagrange multiplier technique. (a) X-displacement and (b) Y-displacement.

Conclusions

As an essential nonlinear source, the contact phenomenon of clearance joint brings abundant dynamic behaviors to the system. The influences of the compound clearance joint on the dynamic responses of the system have received little attention in the previous studies.

The kinematics models of 3D spherical joint and revolute joint with axial and radial clearances are established to this end. According to different contact conditions, the corresponding normal contact force model is provided. The friction phenomena between contact elements are described by a continuous velocity-based friction model.

The dynamics model of the parallel mechanism is derived using Lagrange multiplier technique and Baumgarte stabilization method. The compound clearance joint and the contact forces are embedded in the dynamics model within the framework of multibody dynamics.

The effect of compound clearance on the end effector is mainly reflected in the speed and acceleration response. The initial misalignment of the bearing axis contributes to the large contact force of the transient response and has little effect on the subsequent time history diagram.

From our results, this method can effectively analyze the dynamic response of parallel mechanisms with various types of clearance joints and study the effect of compound clearance on parallel mechanisms. This work is of great significance to the structural design and performance optimization of the parallel mechanism. However, the work should further consider the influence of components’ flexibility and temperature on the system to expand its application to more industrial robots.

Nomenclature

Symbols
B_B, B_T	Bottom and top centers of bearing
c	Clearance of spherical joint
c_a, c_r	Axial and radial clearance of revolute joint
e	Eccentricity vector of spherical joint
e _j, e _b	Unit vectors aligned with the axis of journal and bearing
J_B, J_T	Bottom and top centers of journal
L_j, L_b	Length of journal and bearing
M	System mass matrix
n	Normal contact direction of spherical clearance joint
n ₁, n ₂	Unit vectors of normal direction of auxiliary planes
o_j, o_b	Centers of journal and bearing
q	Generalized coordinate vector
Q	Generalized external forces
r _Mb, r _Mj	Position coordinates of potential radial contact points of revolute clearance joint
r _pb, r _pj	Position coordinates of potential axial contact points of revolute clearance joint
r _psi, r _pbj	Position coordinates of potential contact points of spherical clearance joint
r _si, r _bj	Position coordinates of socket and ball center
R _j	Rotation matrix
R_jo	Radius of journal
R_b1, R_b2	Inner and outer radii of bearing
v _n, v _t	Relative normal and tangential velocities of spherical contact elements
v_s	Transition velocity
α₁,β₁	Baumgarte feedback parameters
Φ _q	Jacobian matrix of constraint equations
Φ _Dk	Driving constraint equations of the system
Φ _Pk	Kinematic constraints of prismatic joint
Φ _Sk	Kinematic constraints of perfect spherical joint
Φ _Uk	Kinematic constraints of universal joint
λ	Lagrange multipliers
μ_s,μ_k	Static and dynamic friction coefficients
ζ	Angle between journal flange and bearing
ξ	Angle between vector e _j and bearing
δ	Relative penetration depth of spherical
δ_M, δ_P	Relative penetration depths of revolute clearance joint

Footnotes

Declaration of conflicting interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This work was supported by Shandong Key Research and Development Public Welfare Program [Grant No. 2019GGX104011] and Natural Science Foundation of Shandong Province [Grant No. ZR2017MEE066].

ORCID iD

Chen Xiulong

References

Kilin

Bozek

Karavaev

, et al. Experimental investigations of a highly maneuverable mobile omniwheel robot. Int J Adv Robot Syst 2017; 14(6): 1–9.

Bozek

Ivandic

Lozhkin

, et al. Solutions to the characteristic equation for industrial robot’s elliptic trajectories. Teh Vjesn 2016; 23: 1017–1023.

Zheng

Wang

Guo

, et al. Dynamic modeling and error analysis of planar flexible multilink mechanism with clearance and spindle-bearing structure. Mech Mach Theory 2019; 131: 234–260.

Bai

Jiang

, et al. Reducing undesirable vibrations of planar linkage mechanism with joint clearance. J Mech Sci Technol 2018; 32: 559–565.

Zhang

Hou

. Accuracy analysis of a spherical 3-DOF parallel underactuated robot wrist. Int J Adv Manuf Tech 2015; 79: 395–404.

Tian

Flores

Lankarani

. A comprehensive survey of the analytical, numerical and experimental methodologies for dynamics of multibody mechanical systems with clearance or imperfect joints. Mech Mach Theory 2018; 122: 1–57.

Zhang

Wang

. Dynamic performance of over-constrained planar mechanisms with multiple revolute clearance joints. Proc Inst Mech Eng C 2018; 232: 3524–3537.

Zhan

Zhang

Jian

, et al. Error modelling and motion reliability analysis of a planar parallel manipulator with multiple uncertainties. Mech Mach Theory 2018; 124: 55–72.

Wang

Huang

. Dynamics analysis of planar rigid-flexible coupling deployable solar array system with multiple revolute clearance joints. Mech Syst Signal Process 2019; 117: 188–209.

10.

Zheng

Zhang

Wang

. Modeling and simulation of planar multibody systems with revolute clearance joints considering stiction based on an LCP method. Mech Mach Theory 2018; 130: 184–202.

11.

Zhao

Zhang

Fang

, et al. Modeling and analysis of planar multibody system with mixed lubricated revolute joint. Tribol Int 2016; 98: 229–241.

12.

Ordiz

Cuadrado

Cabello

, et al. Prediction of fatigue life in multibody systems considering the increase of dynamic loads due to wear in clearances. Mech Mach Theory 2021; 160: 104293.

13.

Sun

Wang

, et al. Correlation dimension and bifurcation analysis for the planar slider-crank mechanism with multiple clearance joints. Multibody Syst Dyn 2021; 52: 95–116. DOI: 10.1007/s11044-020-09769-3.

14.

Ambrósio

Pombo

. A unified formulation for mechanical joints with and without clearances/bushings and/or stops in the framework of multibody systems. Multibody Syst Dyn 2018; 42: 317–345.

15.

Marques

Roupa

Silva

, et al. Examination and comparison of different methods to model closed loop kinematic chains using Lagrangian formulation with cut joint, clearance joint constraint and elastic joint approaches. Mech Mach Theory 2021; 160: 104294.

16.

Erkaya

. Experimental investigation of flexible connection and clearance joint effects on the vibration responses of mechanisms. Mech Mach Theory 2018; 121: 515–529.

17.

Chen

Sun

. Dynamic modeling of spatial parallel mechanism with multi-spherical joint clearances. Int J Adv Robot Syst 2019; 16(5): 255793401.

18.

Cavalieri

Cardona

. Non-smooth model of a frictionless and dry three-dimensional revolute joint with clearance for multibody system dynamics. Mech Mach Theory 2018; 121: 335–354.

19.

Yan

Xiang

WWK

Zhang

. A comprehensive model for 3D revolute joints with clearances in mechanical systems. Nonlinear Dynam 2015; 80: 309–328.

20.

Tian

Liu

Machado

, et al. A new model for dry and lubricated cylindrical joints with clearance in spatial flexible multibody systems. Nonlinear Dynam 2011; 64: 25–47.

21.

Liu

Tian

. Dynamics and control of a spatial rigid-flexible multibody system with multiple cylindrical clearance joints. Mech Mach Theory 2012; 52: 106–129.

22.

Flores

Ambrosio

Claro

J C P

, et al. Spatial revolute joints with clearances for dynamic analysis of multibody systems. Proc Inst Mech Eng K 2006; 220(4): 257–271.

23.

Marques

Isaac

Dourado

, et al. An enhanced formulation to model spatial revolute joints with radial and axial clearances. Mech Mach Theory 2017; 116: 123–144.

24.

Marques

Isaac

Dourado

, et al. A Study on the dynamics of spatial mechanisms with frictional spherical clearance joints. J Comput Nonlinear Dyn 2021; 12(5): 051013.

25.

Lankarani

Nikravesh

. A Contact force model with hysteresis damping for impact analysis of multibody systems. J Mech Design 1990; 112(3): 369–376.

26.

Pereira

Ramalho

Ambrosio

. A critical overview of internal and external cylinder contact force models. Nonlinear Dynam 2011; 63: 681–697.

27.

Flores

Ambrosio

Claro

JCP

, et al. Kinematics and dynamics of multibody systems with imperfect joints. Berlin: Springer, 2008.

28.

Xiang

WWK

Yan

, Dynamic analysis of planar mechanical systems considering stick-slip and Stribeck effect in revolute clearance joints. Nonlinear Dyn 2019; 95(1): 321–341.

29.

Dahl

. Solid friction damping of mechanical vibrations. AIAA J 1976; 14(12): 1675–1682.

30.

Muvengei

Kihiu

Ikua

. Dynamic analysis of planar multibody systems with LuGre friction at differently located revolute clearance joints. Multibody Syst Dyn 2012; 28: 369–393.

31.

Dupont

Hayward

Armstrong

, et al. Single state elastoplastic friction model. IEEE Trans Autom Contr 2002; 47(5): 787–792.

32.

Brown

McPhee

. A continuous velocity-based friction model for dynamics and control with physically meaningful parameters. J Comput Nonlinear Dyn 2016; 11: 054502.

33.

Marques

Flores

Claro

JCP

, et al. Modeling and analysis of friction including rolling effects in multibody dynamics: a review. Multibody Syst Dyn 2019; 45: 223–244.

34.

Haug

. Simulation of spatial multibody systems with friction. Mech Based Des Struc 2018; 46: 347–375.

35.

Baumgarte

. Stabilization of constraints and integrals of motion in dynamical systems. Comput Method Appl Mech Eng 1972; 1(1): 1–16.

36.

Sun

Wang

, et al. Dynamic analysis of the double crank mechanism with a 3D translational clearance joint employing a variable stiffness contact force model. Nonlinear Dynam 2020; 99: 1937–1958.

37.

Zhang

. Minimizing the influence of revolute joint clearance using the planar redundantly actuated mechanism. Robot Computer Integr Manuf 2017; 46: 104–113.

38.

Flores

Ambrósio

Claro

JCP

, et al. Dynamic behaviour of planar rigid multi-body systems including revolute joints with clearance. Proc Inst Mech Eng K 2007; 221: 161–174.

39.

Chen

Jia

. Dynamic response and nonlinear characteristics of spatial parallel mechanism with spherical clearance joint. J Comput Nonlinear Dyn 2019; 14(4): 041010.