Dynamic analysis of an over-constrained parallel mechanism with the principle of virtual work

ABSTRACT This research presents the mathematical modelling of kinematic and complete dynamic analysis of a novel over-constrained parallel mechanism, which consists of two universal-prismatic-revolute joint limbs and one revolute-revolute-universal joint limb. The kinematic model is constructed based on the closed-loop vector method and the velocity Jacobian matrix is deduced, velocity-mapping relationships between all moving components and moving platform are also performed. Afterwards, inertia and applied forces are analysed, the complete dynamic equations with the classical Stribeck friction model of the proposed structure is established based on the principle of virtual work. A theoretical numerical example is given to solve kinematics and dynamics solutions, and theoretical forces from developed dynamic formulation are verified by the physic model simulation in Simscape and the rigid-flexible coupling model simulation in Adams. A good agreement between the theoretical results and multi-body software simulation is found.


Introduction
Due to the excellent mechanical characteristics, such as high stiffness, good dexterity, and large ratio of load capacity to mass [1], the parallel mechanism which is composed of several platforms and several limbs, has been widely applied in many fields, such as ankle rehabilitation manipulator [2], motion simulation platform [3], welding robots [4] and automatic assembly equipment [5]. Among these applications, kinematic and dynamic analysis is the fundamental research of motion control, performance evaluation, and optimal design, so it is quite essential to obtain an accurate dynamic model for parallel mechanism.
However, as a result of the over-constrained parallel mechanism's multiple closedloop configuration and complicated internal forces relationship between the moving platform and active joints, the forward kinematics and dynamic analysis of parallel manipulators are much more difficult than serial robots [6,7]. In general, most of the forward kinematic problems of parallel mechanisms cannot be solved analytically, since the kinematic equations are highly nonlinear and strong coupling [8]. Therefore, it is necessary to study the forward kinematics solution from the numerical point of view. At present, classical methods mainly fall into three categories: polynomial method, numerical method, and intelligent algorithm. Polynomial method is to transfer the kinematic equations to be a univariate polynomial equation by algebraic elimination, such as Sylvester's dialytic elimination [9] and Bezou's elimination method [10], however, the obtained polynomial equation is usually an 8th-, 16th-or higher-order equation and possesses multiple solutions, which might be not efficient enough for the practical applications since a unique solution needs to be determined among them. The numerical method is to solve the forward kinematics by a numerical iterative algorithm, such as Newton-Raphson algorithm and related improved methods [11,12], however, these methods are highly sensitive to the initial guesses and cannot guarantee the desired level of accuracy. Intelligent algorithms [13,14] are able to obtain accurate solutions for the forward kinematics, while these algorithms often take a long time for training and cannot meet the requirements of real-time control. Thus, obtaining an efficient and accurate algorithm for the forward kinematics of the parallel mechanism is still a challenge. For the dynamic analysis of parallel manipulators, Lagrange method and Newton-Euler method are the most popular methodologies in recent years. McGrath et al. [15] adopted Lagrange method to form a generalized formulation for the equations of motion (EOM) of n-link open-loop chains which was used for single support walking models. Dasgupta et al. [16] obtained the dynamics model of 6-SPS (six sphericalprismatic-spherical joint limbs, S and P stand for spherical joint and prismatic joint, respectively) parallel mechanism, in closed form through the Newton-Euler method. Chen et al. [17] established the dynamic models of an over-constrained parallel manipulator by combining the Newton-Euler approach and natural orthogonal complement method. Lagrange method obtains independent equations of motion by mathematical analysis method, which is standardized in the modelling process and without considering internal force. However, when a parallel manipulator possesses intricate structure and coupling degrees of freedom (DoF), the derivation of Lagrange equations is extremely complicated. As an intuitive modelling methodology with clear physical meaning, Newton-Euler method is also widely used in dynamic modelling. However, for an overconstrained parallel mechanism, force analysis is a statically indeterminate problem, Newton-Euler equations and constraint equations are always integrated by complicated internal forces, which makes it rather difficult to solve. Besides, the principle of virtual work is another effective dynamic modelling method for the multibody system [18][19][20], Li et al. [21] analysed equivalent dynamics model of 3-PRC (three prismatic-revolutecylindrical joint limbs, R and C stand for revolute joint and cylindrical joint, respectively) parallel mechanism with this method, however, rotational inertias of legs were neglected which might make it imprecise. Pedrammehr et al. [22] developed the closed-form dynamics of hexarot manipulators by adopting Alembert's form of virtual work principle and performed the validation in SimMechanics. Except for the above methodologies, Kane method [23,24], screw theory [25], Gibbs-Appell [26,27], and Hamilton method [28] have also been commonly adopted in dynamic analysis.
Among these considerable achievements, the principle of virtual work is considered as an adequate methodology for over-constrained mechanism dynamics. There are two main reasons. Firstly, due to the internal forces or moments are eliminated from the dynamic model, complicated internal force or moment analysis could be avoided. Secondly, in the model-based control, inverse dynamics are usually used to compute the actuation forces as feed-forward information to track desired motions, without considering internal forces or moments, the calculation for dynamic formulation could be more efficient, which improves the real-time control process.
In this work, to obtain an accurate and efficient dynamic model for the overconstrained parallel mechanism, which involves two universal-prismatic-revolute joint limbs and one revolute-revolute-universal joint limb, namely 2UPR-RRU parallel mechanism (U stands for the universal joint), the principle of virtual work is applied. In addition, the inclusion of friction forces in the derived rigid body dynamic equations is also discussed. And the correctness of the theoretical dynamic model is validated by comparing it with the results obtained from the virtual prototype simulation. In the remainder of paper, a brief review of the over-constrained 2UPR-RRU parallel mechanism and its kinematic analysis is presented in Section 2, the velocity and acceleration analysis of moving components are discussed in detail. In Section 3, the generalized forces, including the inertial forces, inertial moments, gravity forces, external forces, and friction models are deduced, then the dynamic model is established with the principle of virtual work. Afterwards, the virtual ideal prototype simulation in Simscape and the rigid-flexible coupling model simulation in Adams are conducted to verify the derived model in Section 4. Finally, some concluding remarks are summarized in Section 5.

Topology description and characteristics
The virtual prototype for an over-constrained 2UPR-RRU parallel mechanism is shown in Figure 1, which consists of one fixed platform, one moving platform, one RRU limb, and two UPR limbs. The RRU limb is made up of an active rod and a driven rod, and the UPR limb is formulated by a cylinder and a piston. The first revolute joint of RRU limb and the prismatic joints of UPR limbs are actuated, which are the inputs of the system, the position and orientation of the moving platform is the output. For ease of description, coordinate frames are set in the structure diagram. As shown in Figure 2, a fixed coordinate frameOXYZis attached to the midpoint of B 1 B 2 with z-axis vertical and x-axis pointing towardsB 3 . The moving coordinate framePxyz locates at the midpoint ofA 1 A 2 , the x-axis and y-axis are set along PA 3 andA 1 A 2 , respectively. The fixed platform and the moving platform are the isosceles right triangles.
Referring to the previous research [29], this novel 3-DoF parallel mechanism possesses some characteristics. 1) This mechanism is a novel full-cycle degree of freedom mechanism, possessing one translational DoF along the normal direction of the moving platform, and two rotational DoFs about x-axis in the moving frame Pxyz and Y-axis in the fixed frame OXYZ respectively. And there is no rotational coupling during the moving process of parallel mechanism owing to the fact that the two rotational axes are perpendicular to each other. 2) This mechanism extends the workspace, especially rotation range, owing to the RRU limb emerges a large range of length than other limbs with the same design size, such as RPU limb. Additionally, there is almost no singularity in the workspace. 3) This mechanism is an over-constrained parallel mechanism, based on screw theory, UPR limb applies a constraint force and a constraint couple on the moving platform, RRU limb also applies a constraint force and a constraint couple on the moving platform. Among these constraints, three constraint couples are perpendicular to the moving platform, and two constraint forces are in the same direction. In this way, internal force analysis of this over-constrained parallel mechanism becomes a complicated statically indeterminate problem [30], which makes the dynamic analysis based on Newton-Euler method or similar approaches inefficient.

Position analysis
Due to the fact that the two rotational DoFs are decoupled and two rotational axes are perpendicular to each other. The orientation of the moving platform can be described by two Euler angles αandβ, which denote the rotation angle about x-axis and Y-axis respectively, and another Euler angle is set as 0 for there is no angle motion about the z-axis. P denotes the position vector of point P in the frameOXYZ, which is also the position description of the moving platform.
Assuming that the moving platform rotates about Y-axis withβ firstly and then rotates about x-axis withα, finally translates along z-axis withρ. Thus the transformation matrix of the moving frame relative to fixed frame can be obtained, yields cos β sin α sin β sin β cos α ρ sin β cos α 0 cos α À sin α 0 À sin β cos β sin α cos β cos α ρ cos β cos α Then we can calculate the position vector of point P in the frameOXYZ, thus Note thatP Ox ¼ P Oz tan β,P Oy ¼ 0 .Obviously, the position vector of point P in the frameOXYZ can be defined as P Oz tan β 0 P Oz ½ � T . The rotation matrix T b a of the framePxyzwith respect to frame OXYZ is defined According to the relation of closed-loop vectors, we have where a i (i = 1, 2, 3) denotesPA i in the framePxyz,b i denotesOB i in the frameOXYZ,l i denotesB i A i in the frameOXYZ, the length of each limb chain can be obtained by calculating the module length ofl i , namely ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffiffi P Oz tan βþa sin α sin β ð ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi For RRU limb, as shown in Figure 3,θ B3 represents the angle output of the rotary actuator.r 1 and r 2 denote B 3 C 3 andC 3 A 3 in the fixed coordinate frame OXYZ, respectively, the module length are equal to each other. According to the relationship between angles, where θ ¼ arccos 2r 2 À l 3 2 2r 2 ,θ 1 ¼ arctan P Oz tan βþa cos βÀ d a sin βÀ P Oz . Thus, l 1 ,l 2 , and θ B3 are the inputs of the system, α,β, and P Oz can be set as the outputs of moving platform.
Conversely, forward kinematics, which determines the positions and orientations of moving platform P Oz α β ½ �when the input of each limb l 1 l 2 θ B3 ½ �is given, is also essential in many applications. From Equation (5) and Equation (6), we have Combining Equation (9) and Equation (10), and eliminating intermediate variableβ, we have Fortunately, P Oz can be eliminated and Equation (11) is simplified to a cubic equation with one unknown quantityα, namely In this way, P Oz and βcan be solved by combining Equation (7), Equation (8), and Equation (10) mathematically afterα is solved with Equation (12). It is noticeable that forward kinematics have multiple solutions, and the specific values can be determined by the continuity of the moving platform's trajectory and the range of motion of all joints. Thus, the forward and inverse kinematics of 2UPR-RRU parallel mechanism are solved.

Velocity and acceleration analysis of moving platform
The position of the moving platform is a function of the generalized coordinates βandP Oz . So the velocity of the moving platform could obtain by taking the time derivative, namely where the velocity of the moving platform is given byv h . The angular velocity of the moving platform with reference to frame OXYZ is given byω h , which is obtained from the derivatives of Euler angles.
Combining Equation (13) and Equation (14), we obtain where Taking the time derivative of Equation (15), we obtain wherea h ¼ a hx a hy a hz � � T is the acceleration vector of the moving platform, ε h ¼ ε hα ε hβ ε hγ � � T is the angular acceleration vector.

Velocity analysis of limbs
Referring to the principle of motion synthesis, the velocity of revolute joints A i on the moving platform can be expressed as wherea oi denotes the denotesPA i in the frameOXYZ, namelya oi ¼T b a a i . Then, the velocity along the limbs is whereu ei denotes the unit direction vector of ith limbs, namely u ei ¼ l i = l i j j. For RRU limb, the velocity of A 3 can be obtained in terms of angular velocity, thus where ω B3 represents active rod angular velocity andω C3 describes driven rod angular velocity. c 3 denotes the unit direction vector of ω B3 and ω C3 . Owing to the structure constraints, the angular velocity direction is always parallel to the y-axis, Dot multiplying Equation (19) byr 2 and substitute Equation (17), we have Equation (18) and Equation (20) can be assembled in a matrix form, namely whereJ d represents the Jacobian matrix of the 2UPR-RRU parallel mechanism and describes the relationship between the input velocities and output velocities. Besides, the velocity of A i can also be expressed in terms of the angular velocity of ith limb, whereω qi denotes the angular velocity vector of ith limbs. Cross multiplying Equation (22) both sides byu ei and we have Substituting Equation (17) into Equation (23), the angular velocity of UPR limb can be derived as whereS u ei ð Þrepresents the skew-symmetric matrix [31] of the vector u ei . Therefore, Equation (24) is developed as follows where Obviously, the velocity for the centre of mass of the cylinder can be expressed as whereL si describes the distance between the centre of mass of cylinder and the universal joint on the fixed platform, it's a constant. The angular velocity of the cylinder isω si , which is consistent with the limb's angular velocityω qi . Consequently, the velocitymapping relationship between the centre of mass of the cylinder and the moving platform is obtained as where Similarly, the velocity for the center of mass of the piston is whereL fi describes the distance between the centre of mass of the piston and the joint on the moving platform, it's also a constant. Hence, the velocity-mapping relationship between the centre of mass of piston and moving platform can be also deduced, yields where For RRU limb, the active rod angular velocity can be expressed as matrix form according to Equation (20), thus where Similarly, the velocity and angular velocity for the centre of mass of the active rod can be written as follows v rc1 ω B3 Meanwhile, the driven rod angular velocity is obtained by dot multiplying Equation (17) and Equation (19) withr 1 , thus where Then, the velocity for the centre of mass of the driven rod can be written as So the velocity-mapping relationship between the centre of mass of driven rod and moving platform can be written, v rc2 ω C3 where

Acceleration analysis of limbs
For UPR limb, the angular acceleration is derived by calculating the derivative of Equation (24), where _ u ei the derivative of unit direction and _ v i is the velocity derivative of point A i . Thus,_ u ei The accelerations for the centre of mass of cylinder and piston are obtained by taking the time derivative of Equation (26) and Equation (28) respectively, namely wherea qi is the linear acceleration along with the limbs, which can be obtained by calculating the derivative of Equation (18), where For RRU limb, taking the derivative of Equation (17) and Equation (19) respectively, we obtain whereε B3 represents the angular acceleration of the active rod, ε C3 represents the angular acceleration of the driven rod,a 3 represents the acceleration of point A 3. Dot multiplying Equation (39) and Equation (40) byr 1 , we have Similarly, dot multiplying Equation (39) and Equation (40) byr 2 , we have As soon as angular accelerations are known, the accelerations for the centre of mass of the active rod and driven rod can be obtained by differentiating Eq. (31) and Eq. (33), namely where a rc1 represents the acceleration of active rod, a rc2 represents the acceleration of driven rod.
As a result, the velocities, angular velocities, accelerations, and angular accelerations of all the moving components are deuced.

Inertia and applied forces analysis
The gravity of the moving platform is given by where m h denotes the mass of the moving platform. The inertia force of the moving platform in the frame OXYZ is given by The inertia moment of the moving platform in the frame OXYZ is given by where J m represents the 3 × 3 inertia tensor of moving platform with respect to point P.
Assuming that external forcef h0 and external momentn h0 are applied on point P respectively, which needs to be emphasized is that they are described relative to the fixed coordinate frame OXYZ and when the centre of mass of moving platform does not locate on point P, the external moment introduced by gravity should be considered. Thus, the resultant of the inertia and applied forces acting on the moving platform can be expressed as Similarly, the gravity, inertia force, and inertia moment of the cylinders are where m si denotes the mass of the cylinder. J ssi represents the 3 × 3 inertia tensor of the cylinder with respect to its centre of mass.T b sfi is the transformation matrix of the cylinder relative to the fixed frame OXYZ, which can be defined as where k ei denotes the unit direction of the axis of revolute joints A i (i = 1,2). Thus, the resultant of the inertia forces can be expressed as For the pistons in UPR limbs, we have where m fi denotes the mass of the piston. J ffi represents the inertia tensor of piston with respect to its centre of mass. Similarly, the resultant of the inertia forces acting on the pistons can be expressed as For active rod and driven rod in RRU limb, we have In these equations, m r1 and m r2 denotes the mass of active rod and driven rod respectively. J rr1 and J rr2 represents inertia tensor of rods. Similarly, T b r1 is the transformation matrix of the rod relative to the fixed frameOXYZ, which can be implemented with where r ei is the unit direction of the vector r i (i = 1,2).

Friction model
For actual industrial robots, the steady-state error or oscillation in the system caused by the nonlinear joint friction such as Coulomb friction, viscous friction play a nonnegligible aspect in the dynamic characteristics of mechanisms. Therefore, it is necessary to incorporate the friction contribution into the rigid body dynamic equations for a parallel robotic system. In this section, the classical Stribeck friction model [32] is introduced in actuator joints. For each prismatic joint, the friction model can be expressed as where F c denotes Coulomb friction force, F s denotes static friction force, v s is stribeck velocity, E denotes the viscous friction coefficient, the empirical exponent δdepends on the material surfaces and ranges. Similarly, the friction model for revolute joint can also be established. These friction model parameters can be identified with some experiments [33].

Dynamic equation
Base on the principle of virtual work, the work produced by the active force/torque on the driving limbs should be equivalent to the work produced by inertia force/moment on the mechanism, and the friction forces can be regarded as the active force, considering D'Alembert's principle, it can be stated as denotes the active forces for all the active actuators, F qi (i = 1,2) represents the active forces on prismatic joints of limb UPRs.T q3 represents the active torque on the first revolute joint of RRU limb. τ T fi ¼ τ q1 τ q2 τ q3 � � denotes the friction forces/torque for these active actuators, δq denotes the instant and small virtual motions for all the active actuators,F T j represents the inertia and applied forces acting on the moving components of the mechanical system, δp j represents the small virtual motions for these moving components. Dividing both sides of Equation (57) by timeδtand substituting Equation (48), Equation (51), Equation (53), and Equation (54) into Equation (57), we have Substituting Equation (15), Equation (21), Equation (27), Equation (29), Equation (31), and Equation (34) into Equation (58), the dynamic equation with friction model of 2UPR-RRU parallel mechanism can be obtained after simplification, yields where It is noticed that J dh describes the relationship between the input velocities and output velocities directly, combining Equation (15) and Equation (21), we obtain v q1 v q2 ω B3 Obviously, J dh is a 3 × 3 Jacobian matrix, when the eigenvalue of J dh equals 0, there exists no inverse matrix of J dh , active forces may become numerically unstable in numerical calculation. In physic, it means the mechanism is in a singular configuration, the values of the active forces may become very large or discontinuous, which may cause a problem for the robot. Fortunately, it has been proved that there is almost no singularity in the workspace [29].

Validation with Simscape
In this section, to validate the developed mathematical model for the 2UPR-RRU mechanism above, the results of the developed dynamic formulation based on the virtual work principle are taken into comparison with the results obtained by Simscape. Although Simscape is one of the physical modelling methods, which provides an effective way to quantify the actuator forces with respect to their physical nature [34], the model in Simscape cannot be used for real-time control and it is difficult for engineers to develop control laws or conduct dynamic performance evaluation without dynamic equations.   To validate the ideal dynamic equation, the nonlinear joint friction model is ignored in this section. The verification procedure is illustrated in Figure 4. Inverse dynamics is the analysis of the forces and torques required on the actuators to produce a given motion trajectory. For developed dynamic formulation, the desired trajectory path, velocities, and accelerations of the moving platform are given to the analytical model as an input. Meanwhile, the physical properties of the parallel mechanism are also required which are listed in Table 1. Then the actuated joint motions and torques are obtained directly with codes in Matlab.
The inertia tensor of the moving platform (unit: kg·m 2 ) respect to point P is For the physic model in Simscape, it is complicated to derive the active forces directly from the moving platform trajectory. So with the assistance of inverse kinematics, we obtain the three actuator motion trajectories which are set as the inputs for the model, and subsequently, we can calculate the forces and torque with the physic model. To do this, the 3D model of the 2UPR-RRU parallel mechanism has been established in SolidWorks, and then the 3D model is converted to the XML format which can be imported into Simscape conveniently. As a result, the model of the 2UPR-RRU parallel mechanism including  a detailed description of blocks is illustrated in Figure 5. It is worth noticing explicitly that this model involves all the physical properties of the proposed mechanism.
As shown in Figure 5 (a), the model consists of a fixed platform, three limbs, a moving platform, converter blocks, inputs, and outputs. As it is mentioned above, three input blocks get the actuator motion data from the inverse kinematics. The Simulink-PS Converter block and PS-Simulink Converter block transforms the units we need. The actuator forces and torque are collected and shown in output blocks. Blocks UPR Limb1, UPR Limb2, and RRU Limb represent three limbs respectively. Figure 5 (b) shows the detailed information of model blocks of UPR limb, which is composed of universal joint, cylinder, prismatic joint, piston, and revolute joint. The first block and last block connect to the fixed platform and moving platform respectively. INPUT block and OUTPUT block represent the displacement and driving force of prismatic joint respectively.
In this validation procedure, the trajectory of the moving platform with respect to time is defined by three generalized coordinates, namely P Oz t ð Þ ¼ À 0:42 À 0:03 sinðπtÞ The units of the above expressions are metre, rad and rad respectively. The mentioned trajectory is a complex trajectory since orientation and position of the moving platform change all the time.The external force is set asf h0 ¼ ½0; 0; 100� T Nand the external moment is set asn h0 ¼ ½0; 0; 0� T þ T b a r h À � � m h gðNmÞ, r h is the vector from point P to the centre of mass of the moving platform. The displacements of limbs are deuced with inverse kinematics, which are depicted in Figure 6. Substituting this trajectory and external force into the developed dynamics model and physic model in Simscape respectively, we obtain the force curves and torque curve of active prismatic actuators and the active rotary actuator, which are depicted in Figure 7.
As shown in Figure 7, the magnitude and changing trend of inputs are fully consistent with each other. Namely, the results of proposed dynamic equations based on the principle of virtual work and physic model in Simscape are almost equivalent, which indicates the correctness of the proposed dynamic model. Meanwhile, there are two speculated reasons for the differences in result. First, the mass of joints are not considered in the developed dynamic formulation, however, the model in Simscape exists these components. Second, the inertia tensor of the moving components in the dynamic formulation may not be accurately consistent with the physic model in Simscape, due to the complex shapes. Additionally, to compare this method with existing methods, we also used the Newton-Euler method to obtain the dynamic equations of this manipulator, with constrained forces analysis, the computational time with dynamic formulation has been reported as 0.572s, however the results with virtual work only cost 0.223s, which shows the excellent computational efficiency of the developed dynamic model.

Validation with rigid-flexible coupling model
The developed dynamic formulation based on virtual work principle and physic model in Simscape are based on the assumption of the rigid body and ideal parallel mechanism, which may not fully represent the real system since there exists flexible bodies, clearance joints, and other non-ideal factors in the real world mechanism. Thus, a comparison to a rigid-flexible coupling model is conducted in this section.
With the aid of two kinds of commercial software, namely Adams and Ansys, a rigidflexible coupling dynamic model for the mechanism is established and the co-simulation analysis is investigated. The flowchart of co-simulation is shown in Figure 8.
Pistons and driven rod of the parallel mechanism are processed to flexible bodies in Ansys, other components are still set as rigid bodies. It is noteworthy that, rigid regions of flexible bodies need to be established to connect with rigid bodies in Adams. The finite element model and rigid region of the piston and the driven rod is shown in Figure 9. After introducing these flexible bodies into Adams and setting up constraints and motion curves, a rigid-flexible coupling dynamic model for 2UPR-RRU parallel mechanism is established, which is depicted in Figure 10.
Similarly, three actuator motion trajectories are set as the inputs for the model, with the assistance of inverse kinematics. The driving forces are depicted in Figure 11. As shown in Figure 11, the magnitude of driving forces have a slight difference from the ideal model. This is most likely due to the fact that the deformation caused by these flexible links are very small, which has tiny effect on the dynamic response of this parallel mechanism. The results also indicate the correctness of the proposed ideal dynamic Figure 11. Input forces at active actuators for limbs with two models. model. However, in the work condition of high frequency, reciprocating motion, and heavy load, the flexibility effect of parallel robot might be obvious, the influence of these nonlinear factors needs to be further investigated [35].
The accuracy of the developed dynamic formulation based on virtual work principle was also performed for several different trajectories, the results show the accuracy and efficiency of the analytical methodology. Ideal rigid dynamic analysis with this methodology only requires the kinematics information of moving components, namely position, velocity, and acceleration. With constrained forces/moments eliminated, it is extremely advantageous for the parallel mechanism with overconstrained wrenches. However, when the flexibility effect of the mechanism cannot be neglected under some circumstances, the methodology for the nonlinear dynamic model should be considered.

Conclusion
To obtain an accurate dynamic formulation, this paper presents a general procedure to conduct the complete dynamics analysis of the over-constrained 2UPR-RRU parallel mechanism. The kinematic analysis was performed which included the position, velocity, and acceleration problems and Jacobian matrices, the complete dynamic model with friction model was developed by the principle of virtual work. Theoretical required active forces by actuators are finally obtained to follow the proposed trajectory of the moving platform. Simulations and comparisons with a Simscape virtual prototype and a rigid-flexible coupling dynamic model showed the correctness of the developed theoretical dynamic model. Without considering the analysis of complex internal forces, the virtual work method is more computationally efficient. This work paves the way for the control strategy research of 2UPR-RRU parallel mechanism in the future. The methodology can also be used for dynamic analysis of similar parallel mechanisms.