Issue 
Mechanics & Industry
Volume 24, 2023



Article Number  28  
Number of page(s)  20  
DOI  https://doi.org/10.1051/meca/2023017  
Published online  18 August 2023 
Regular Article
Rigidbody inverse dynamics of a spatial redundantly actuated parallel mechanism constrained by two point contact higher kinematic pairs
^{1}
Lab of Mechatronics, Xi'an College of Technology, 710025 Xi'an, China
^{2}
School of Mechanical Engineering and Automation, Fuzhou University, 350116 Fuzhou, China
^{3}
Institute of Robotics and Automatic Information System, College of Artificial Intelligence, Nankai University, Tianjin 300350, China
^{*} email: cche943@sina.com
Received:
4
August
2022
Accepted:
30
May
2023
This paper presents a comparative study of the rigidbody inverse dynamics of a spatial redundantly actuated parallel mechanism constrained by two point contact higher kinematic pairs (HKPs). Firstly, its constrained motions are analysed comprehensively, then four different models are formulated by the generalized momentum approach and the LagrangeD'Alembert formulation to explore its inverse dynamics. In each method, the first model is built by employing the method directly to the mechanism. In the second model, the dynamic model of its nonredundantly actuated counterpart free of HKPs is built by this approach first, then the constraints from HKPs are modelled, to finally reach the model of the redundantly actuated parallel mechanism (RAPM) where that of its counterpart is utilised as the core. The four models give rise to equivalent numerical results, and the second model in both methods of the RAPM can alleviate the strong coupling between the parasitic motion variables and degrees of freedom (DOFs), boosting the computational speed as fast as that of its nonredundantly actuated counterpart without simplification or loss of accuracy. The comparisons between the mechanism and its counterpart validate that the HKP constraints greatly increase the computational complexity, and the torques required by the parasitic motions of the end effector are significantly smaller than those by the corresponding DOFs.
Key words: Inverse dynamics / higher kinematic pair / spatial parallel mechanism / redundant actuation
© C. Cheng et al., Published by EDP Sciences 2023
This is an Open Access article distributed under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
1 Introduction
In the food industry, there is a great interest in evaluating the timevarying dynamics of newly developed food textures during the chewing process. To this end, a machine that can accurately replicate the humanlike chewing behaviours in a biomimetic fashion can come into play a significant role in terms of the threedimensional (3D) chewing motions and bite forces. Inspired by the masticatory system of human beings, a spatial redundantly actuated parallel mechanism (RAPM) constrained by two point contact higher kinematic pairs (HKPs) has been developed [1,2]: the base is the skull, the six RSS (revolutesphericalspherical) kinematic chains are the primary chewing muscles, the end effector is the mandible, and the two point contact HKPs are the left and right temporomandibular joints (TMJs), respectively. The underlined letter means the joint is active in the chain. At the early stage of our study in the bioinspired masticatory robot, the six RSS kinematic chains were designed to play the role of primary chewing muscles in the human's masticatory system, as in Chapter 4 of [3]. However, during the motions of the lower mandible, it is not only driven by chewing muscles, but also guided at TMJs at its left and right condyles. In this regard, to better replicate the chewing system and produce humanlike chewing behaviours, TMJs were modelled using two HKPs as in [1,2]. Readers interested in the human masticatory system and chewing robotics can refer to [3] for a comprehensive description.
Before the applications of the RAPM in practice, developing a precise and computationally efficient inverse dynamic model is a prime requirement. In fact, in our previous work [2], an initial attempt at the rigidbody inverse dynamics of the RAPM of interest has been made via the hybrid of the NewtonEuler's law and the Lagrangian equations. However, due to the complicated constraints from HKPs at the end effector in the RAPM, the model was quite sophisticated and errorprone. For instance, the forces at HKPs without friction effects are ideal constraints on the end effector which cannot be included in the Lagrangian formulation, whilst they must be considered in the NewtonEuler's law. Meanwhile, the Lagrange formulation introduces many unknown Lagrange multipliers, greatly increasing the number of variables in the equations of motion (EOMs). In these regards, finding efficient solutions to its inverse dynamics is a strong motivation for the study in this paper.
From the literature review, the conventional methodologies in dynamics of serial mechanisms and PMs without redundant actuation have been successfully applied in a wide range of RAPMs. For instance, the NewtonEuler's law was used in a planar 3DOF RAPM to find its maximum dynamic loadcarrying capacity and the impact of actuation redundancy [4]. The Euler's equation was applied in a 4RRR planar RAPM to achieve its shaking force/moment balancing [5]. In [6], this law was employed to solve the inverse dynamics of the spatial RAPM with parallelograms inside a serialparallel hybrid machine. One can find that under this law, the free body diagram must be drawn to each body, and there are a number of intermediate variables associated with the constraint wrenches in the joints. Therefore, its computational complexity and efficiency would be negatively influenced by the increased number of actuators and moving links in RAPMs. In [7], the Lagrange formulation with unknown Lagrangian multipliers was used to address the dynamics of a planar RAPM with two translational DOFs. The multipliers representing the magnitude of constraint forces were eliminated using the null space of the differential matrix of the closedloop constrained equation. The identical procedure was also employed in [8,9]. The utilisation of the Lagrange equations in the modelling process is accompanied by the introduction and elimination of Lagrangian multipliers, and it involves a number of partial derivatives. As a result, the computational demands are not light. The principle of virtual work was used in evaluating the superiority of redundant actuations in the dynamic performance of an 8PSS spatial RAPM in [10]. In [11], this principle was applied to compare the dynamic performance of planar PMs with both actuation and kinematic redundancies. It was employed in [12] to calculate the joint torque to compensate for disturbances to the mobile platform of a dualarm underwater robot. It also has been used in [13–15] where the EOM was used for modelbased control. In this principle, the ideal constraint wrenches that do not generate virtual work are not taken into account, and it does not require as many partial deviations as in the Lagrange formulation. As such, compared to the foregoing two methods, its computational cost is lower.
Apart from them, another two methods in the literature attract our attention. The first is the newly developed generalized momentum approach in [16,17]: it has been successfully applied in both the conventional GoughStewart platform and a 6PUS (prismaticuniversalspherical) PM. In this approach, the kinetic and the gravitational components of generalized forces acting on a moving body are computed by the differentiation of its generalized momentum and potential energy, respectively, and they can be expressed by the generalized coordinates. In addition, the contributions to the total actuating torques by each moving body can be found separately. As such, compared to the classical NewtonEuler' law, it is characterised by being wellstructured, straightforward, and systematic. The second method is the classical LagrangeD'Alembert formulation, by which the dynamic models of the planar 2DOF RAPM from [8] and 3DOF purely translational PM were built in [18,19], respectively. In this approach, the dynamic model of each rigid body is firstly built by the Lagrange formulation, then following the D'Alembert principle, the EOM of the mechanism is reached. Compared to the procedure by the Lagrange's equations, it is not indispensable to introduce and then eliminate the unknown Lagrangian multipliers. In addition to the two methods' respective advantages, they have two common features: the constraints imposed by the closed loops in the mechanism are embedded via the velocity relationships between the end effector and the chains, obviating the use of explicit constraint equations; they do not take constraint wrenches in joints that do not generate virtual work into consideration. As such, a minimal set of variables and equations is generated and it is straightforwardly beneficial to the modelbased motion control.
Even though the successful applications of these two methods in RAPMs are not quite extensive, due to these foregoing merits, an attempt is going to be made to extend them to the RAPM with HKPs which has a sophisticated topology, and this is the first reason to use them in this paper. The second is, the computational cost of these two methods in RAPMs has not been compared to date, despite their common advantages. Thus, it is necessary to find which is faster. The third reason is, the software at hand i.e., Matlab does not support redundant actuation mode in simulations of inverse dynamics. Thereby, it is required to use two methods simultaneously to crosscheck the correctness of the numerical results, and the foregoing two methods are adopted in this paper.
In the meantime, from the view point of mechanism, there are primarily two methods to implement actuation redundancy in PMs as stated in [20]: the first is adding identical kinematic chains as the existing ones in a PM, and the second is actuating some of the passive joints of PMs. By a sharp contrast, from the analysis in Section 3 in this paper, redundant actuations are realised by imposing two point contact HKPs at the end effector of a 6RSS PM. Evidently, by this process, both the RAPM and its nonredundantly actuated counterpart, namely, the 6RSS PM have identical numbers of moving bodies and actuators, respectively. Nonetheless, intuitively, the kinematic variables in the RAPM are more strongly coupled, leading to more complex kinematics and dynamics, and the computational demands would be higher than that of its counterpart. Furthermore, studies on this sort of RAPMs with HKP constraints are quite rare: from the literature, the publications are almost about RAPMs with lower kinematic pairs and their end effectors are not constrained by the base directly. Thus, under the adopted methods, finding a shortcut that can alleviate the coupling among variables to increase the computational speed of the dynamic model is the second strong motivation for conducting the study in the paper. In fact, the discovery of this shortcut would indicate the nature of HKP constraints in the inverse dynamics of the RAPM is explored deeply.
In this paper, the rigidbody inverse dynamics is studied. It is assumed that all bodies and joints are rigid, frictionless, and free of clearances in the mechanism. The inertia of the spherical joints is rather small and then ignored. The remainder of the paper is organized as follows: a detailed description of the mechanism is given in Section 2. In Section 3, the kinematics of the mechanism including the constrained end effector and the kinematic chains is analysed comprehensively. The inverse dynamics of the RAPM is analysed by the generalized momentum approach in Section 4, in Model 1, the approach is directly employed to the RAPM, while in Model 2, firstly the dynamic model of the 6RSS PM is built by this approach, then the two constraints from HKPs are modelled, finally, the model of the RAPM is built by incorporating the model of the 6RSS PM as its core. In Section 5, an identical procedure as in Section 4 is implemented by using the LagrangeD'Alembert formulation. Section 6 is devoted to the numerical computations and discussions in the four models by the two methods, and the influences of the HKP constraints. Finally, some conclusions are made in Section 7.
The three contributions of this paper are:

By using the inverse dynamic model of the 6RSS PM as the core of that of the RAPM can alleviate the coupling in its parasitic motion variables and increase its computational speed as fast as that of the 6RSS PM.

The insight into the nature of HKP constraints on the inverse dynamics has been provided clearly in terms of the modelling process, the computational cost, and the numerical results.

The applications of the generalized momentum approach and the LagrangeD'Alembert formulation have been extended to the RAPM constrained by HKPs.
2 Description of the mechanism
The scheme of the RAPM constrained by two HKPs is illustrated in Figure 1. The maxilla (i.e., the base) is fixed on the ground and the movable mandible (i.e., the end effector) is connected to the base by six independent kinematic chains. The maxilla, to which the inertia frame {S} is assigned, is not shown in the figure for a clear exhibition of moving bodies. This frame consists of a horizontal X_{S}Y_{S} plane perpendicular to the vertical Z_{S} axis. A frame {M} is established at the mass centre O_{M} of the end effector. The origins and orientations of {S} and {M} overlap when the mechanism is at the home position, that is, the maxilla and the mandible are in the occlusal state. The origin O_{M} is used as the reference point to describe the mandibular translations, and its orientations with respect to {S} are expressed by Bryant angles. In each chain, the crank G_{i}S_{i} (i =1,…,6) is driven by a basemounted rotary actuator, whose shaft connects it with a rotational joint at G_{i}, and a coupler S_{i}M_{i} that joins the crank and the end effector via two spherical joints at its two ends S_{i} and M_{i}, respectively. The layout of six chains is in accordance with that of the six primary chewing muscles. The rotation of the ith actuator with respect to {S} is described by the actuator frame{C_{i}}attached at G_{i}. In it, the. .axis is directed from G_{i} to S_{i}, the axis runs through the driving shaft of the actuator, and the axis completes the frame, obeying the righthand rule. A frame {N_{i}} is attached at the mass centre E_{i} of S_{i}M_{i} to describe its motions with respect to {S}. The axis points from S_{i} to M_{i}, the axis is parallel to the cross product of two unit vectors defined along the and X_{S} axes, and the axis is defined by the righthand rule.
Two HKPs modelling the left and right TMJs are formed by the two condyle balls being in contact with the articular surfaces. The two contacts between the left and right condyles of the mandible and the maxilla at TMJs, are realised by the two condyle balls being always constrained within a condylar socket in the mechanism prototype, as shown in Figure 2. The width of the socket equals the diameter of the condyle ball; thus, it can always guarantee the point contact during the movement of the mandible. By this design, the motion of the condyle ball centre is always constrained onto a surface, which is offset from the upper and lower surfaces of the socket by the ball radius. Thereupon, it is clear that the end effector is actuated by six chains and constrained by the base at the two HKPs simultaneously. Besides, the end effector is shown in Figure 2d. the green component is a manufactured mechanical part, while the red lower mandible and golden teeth are scanned from a human cadaver by computer tomography in the School of Dentistry at the University of Otago, Dunedin, New Zealand [3]. The prototype of the mechanism is provided in Figure 3.
Fig. 1 A schematic view of the RAPM constrained by two point contact HKPs, where ① right condyle ball, ② left condyle ball, ③ articular surface of right TMJ, ④articular surface of left TMJ [21]. 
Fig. 2 The RAPM: (a) CAD model, (b) magnification of the right HKP, (c) prototype of the HKP [2], (d) end effector. 
3 Kinematics of the mechanism
A priori development of accurate kinematic analysis is fundamental to the inverse dynamics. In this section, the kinematics of the end effector and the chains in terms of their displacements, velocities, and accelerations is analysed comprehensively.
3.1 Constrained motions of the end effector
A secondorder surface was used as the workspace of the centre of the condylar ball in [1,2,21,22]; its shape and size was designed by referring to [23–25]. However, on the one hand, when the mechanism tracks the real chewing trajectories of healthy human subjects, for instance, the one in Section 6 which serves as an illustrative example in this paper, the surface where the path of the condylar ball centre is situated can be approximated as a flat one. On the other hand, it is quite difficult to derive explicit analytical expressions of the parasitic motions in the RAPM under the secondorder surface. In these regards, in this paper where the chewing system is explored from the viewpoint of mechanical dynamics and with an emphasis on the constrained dynamics of the end effector, the surfaces in {S} where the left and right condyle ball centres T_{L} and T_{R} slide on are designed as flat (unit: mm):
From the KutzbachGrübler criterion, the mechanism now has four DOFs, but the information on which four to choose is not given. They are to be derived from a rigorous computation below: the coordinates of T_{i} (i = L,R) in {S} can be expressed as
where denotes the 3◊1 position vector of O_{M} in {S}, is the rotation matrix from {S} to {M}, R_{X} (α) , R_{Y} (β), andR_{Z} (γ)are three rotation matrices about the X_{M}, Y_{M}, and Z_{M} axes by three Bryant anglesα, β, andγ, respectively. It is worth noting that in this paper, a matrix/vector/scalar in local frames owns a leading superscript on its left to denote the specific frame it refers to, but those in {S} omit their superscripts for the sake of convenience and clarity. The six motion variables of the end effector are grouped and expressed as
From equation (2), one can obtain
where is the ith(i = 1,3) row of . Putting equation (4) into equation (1) produces
Because of the leftright symmetry ofMO_{M}T_{L}andMO_{M}T_{R}in {M}, a summation and a subtraction of the two equations in equation (5) sidewise yield
where MO_{M}T_{L (1)}and MO_{M}T_{L (3)}are the first and third terms ofMO_{M}T_{L}, respectively. From these computations, it is found Z and γ are transferred from DOFs to parasitic motions and they are functions of q_{EE}, which is a 4 × 1 vector by grouping four DOFs as
and it constitutes the task space of the mechanism. To characterise the instantaneous configuration of the mechanism, equation (3), or both equations (6) and (7) ad hoc are needed. In other words, the RAPM can still perform motions in six directions with four DOFs and two parasitic motion variables. Regarding this, redundant actuations in the mechanism are essentially caused by constraints from the base directly onto the end effector, which is completely different from the two methodologies mentioned in [20]. It is also worth noting that though the workspace of the centre of the condylar ball is simplified as a flat surface as in equation (1), a strongly nonlinear and sophisticated relationship between Z/γand q_{EE} in equation (6) can be observed.
After figuring out the DOFs of the end effector, its motions can be defined below. The angular velocity is
where
0_{3} is 3 × 3 the zero matrix, and M_{J} denotes the 6 × 4 Jacobian matrix between X_{EE} and q_{EE}, namely
From it, the second timederivative of X_{EE} can be derived as
The translational velocity of the end effector is expressed by
where E_{3} is the 3 × 3 identity matrix.
3.2 Kinematics of the ith chain
The inverse kinematics of the mechanism, i.e., that consists of a system of six decoupled equations expressed byq_{EE}, has already been derived in Section 3.2 of [22]. Nonetheless, the motions of the coupler S_{i}M_{i}(i = 1,…,6) are still needed for its rigidbody dynamics. Due to the two spherical joints at S_{i} and M_{i}, the coupler can rotate around the three orthogonal axes of {N_{i}}. The rotation around axis is a passive DOF for it is not controllable; this rotational range is quite small thanks to the physical restrictions from the used spherical joints in the mechanical design, however. In these regards, it is assumed that there is no axial rotation in the coupler. As such, two Euler angles β_{i} and γ_{i} around the and axes, respectively, are used to express the rotation of S_{i}M_{i} in {S} in terms of the two rotational matrices R_{Y} (β_{i}) and R_{Z} (γ_{i}). Thereafter, the coordinate vector of the coupler can be expressed as
where is the orientation of S_{i}M_{i} in {S} at the initial configuration of the mechanism at hand, and∥S_{i}M_{i}∥is the length of S_{i}M_{i}. From the geometry of the mechanism in {S}, the vector of the coupler can also be found from the difference in the position vector of M_{i} and S_{i}:
where and is the coordinate vector of M_{i} in {M}; O_{S}G_{i}is the constant position vector of G_{i} in {S}; and , in which is the orientation of G_{i}S_{i} in {S} at the initial configuration of the mechanism,R_{Z} (θ_{i}) is the rotation matrix about the axis by θ_{i}, and ∥G_{i}S_{i}∥is the length of G_{i}S_{i}.
Substituting equation (12) and θ_{i} = θ_{i} (q_{EE}) into equation (13) produces
where . For the sake of convenience, a3 × 1 generalized vector is defined as , which consists of the joint space of the ith chain and completely specifies its configuration. It is highlighted that is the function of q_{EE}, i.e., .
To derive the relationship between the first timederivatives of andq_{EE}, equation (13) can be rewritten as
The two sides can be expressed by and q_{EE}, respectively. The first timederivative of equation (15) yields
where
Moreover, one can find that
where , and the second timederivative of is given as
So far, ˙q_{ri} ˙riand ˙ ˙q_{ri}¨ rihave been derived as quantities intimately associated with the motions of the ith chain.
The angular velocity of S_{i}M_{i} is
where . Substituting equation (17) into equation(19) renders
The coordinate vector of the mass centre E_{i} in {S} is
The translational velocity of E_{i} can then be found by the differentiation of equation (21) with respect to time as
Upon substitution of equation (17) into equation (22), it yields
As far as the crank G_{i}S_{i} is concerned, its mass centre G_{i} is passing through the rotational shaft of the actuator, as such . Its rotational velocity in {S} is
where is the third column of .
Heretofore, the motions of all moving bodies in the mechanism have been explicitly expressed. These are all the ingredients needed for the construction of the EOM based on the two dynamic methods in the following two sections. In each of them, two models are built by the generalizedmomentum approach and the LagrangeD'Alembert formulation, respectively. In Section 4, the generalizedmomentum approach is directly used to the RAPM in Model 1; in Model 2, by following an intuitive procedure, the EOM of the 6RSS PM free of HKPs is firstly written, then the two constraints from HKPs are modelled, to finally reach the EOM of the RAPM where the model of the 6RSS PM is used as its core. An identical procedure is implemented in Section 5 by the LagrangeD'Alembert formulation. The structure of dynamic modelling in the following two sections are given in Figure 4.
Fig. 4 The structure of dynamic modelling by the two methods in Sections 4 and 5. 
3.3 Singularity analysis
To analyse the singular configurations of the mechanism, a 6◊1 constraint vector is defined as
Its first timederivative is
where
and their sizes are 6◊6 and 6◊4, respectively. It is also noted that the ith(i = 1,…,6) row of L (θ, q_{EE}) is the function of θ_{i} and q_{EE}, then the nondiagonal entries of J_{A} must be zero, i.e., J_{A} is a diagonal matrix.
From Section 5.3 of [26], there are three sorts of singularity conditions:
1: An inverse kinematic singularity occurs when J_{A} is rankdeficient, i.e.,
Thus, at least one diagonal entry of J_{A} is zero. In this case, a motion of the actuators does not lead to the displacement of the end effector.
2. A direct kinematic singularity occurs when J_{W} is rankdeficient, i.e.,
In this case, the end effector can possess uncontrollable motions in some directions while all actuators are locked.
3. Combined singularities occurs when both equations (27) and (28) are valid. In this case, the end effector can remain stationary while the actuators undergo some infinitesimal motions. On the other hand, it also can undergo infinitesimal displacements in some directions while all actuators are locked.
4 Generalized momentum approach
In writing the EOM by this approach, it must be kept in mind that the forces acting at the end effector via the two HKPs are ideal constraint forces rather than active forces, as such, they do not appear in the formulation. In this section, Model 1 is built by directly applying this approach to the RAPM, and in Model 2, the dynamic model of the nonredundantly actuated counterpart, namely, the 6RSS PM is firstly built, then the HKP constraints are modelled, to arrive at the dynamic model of the RAPM finally.
4.1 Model 1
Under this approach, the kinetic and gravitational components of the generalized forces of each moving body is obtained by differentiating its generalized momentum and potential energy with respect to time and q_{EE}, respectively. The EOM of the RAPM is acquired by the law of power equivalence: the power generated by the generalized forces of all rigid bodies in the mechanism is equivalent to that by actuating torques from actuators and external wrenches from the environment.
4.1.1 The end effector
On the one hand, the generalized momentum of the end effector is computed as
where is the 6 × 6 inertia dyad, m_{EE} and I_{EE} are the mass and the inertia matrix with respect to O_{M} and expressed in {S}, respectively. By equations (8) and (11), it can be found that , and the size of q_{EE}_{1} is 4 × 1.
The kinetic component of the generalized force acting on the end effector can be computed as
On the other hand, the potential energy of the end effector is
where g = 9800 m/s^{2} is the gravitational acceleration. The gravitational component of the generalized force is
4.1.2 The ith chain
As far as the coupler S_{i}M_{i} of the ith kinematic chain is concerned, by virtue of equations (20) and (23), its generalized momentum is computed as
where is the 6 × 6 inertia dyad of S_{i}M_{i}, and are the mass and the inertia matrix with respect to O_{M} and expressed in {S}, respectively, and . The kinetic component of the generalized force acting on S_{i}M_{i} is computed as
The potential energy of S_{i}M_{i} is
See Formula in PDF
whereO_{S}E_{i (3)} is the third term of O_{S}E_{i}. The gravitational component of the generalized force is
For the crank G_{i}S_{i}, its generalized momentum is
where is a constant6 × 1vector, and the scalar is the rotational inertia of G_{i}S_{i} around the axis of frame{C_{i}}. As such, its kinetic component of the generalized force is
Because G_{i} is a fixed point in {S}, its potential energy is constant and computed as
in whichO_{S}G_{i (3)} is the third term of O_{S}G_{i}. Thereby, its gravitational component of the generalized force is
4.1.3 The overall mechanism
In the entire mechanism, the instantaneous power generated from the sum of the generalized forces on all the moving bodies is balanced by that from the input torques and bite force. As such, one can write
whereτis the6 × 1actuating torque vector by the six actuators, and is the6 × 1wrench vector by the bite force F_{B}. Substituting equations (8), (11), and (17) into equation (41), and rearranging the term by give rise to
4526+++++++where is the generalized force vector by F_{B}, is the 6 × 4 Jacobian matrix mapping q˙_{EE } into θ˙, i.e., θ˙=J_{θ1}⋅q˙_{EE}, and is the first row of. Because q_{EE} is the independent generalized coordinate vector of the entire mechanism under study, q˙_{EE} is free to vary, thus one can obtain
which is the inverse dynamic model of the RAPM and it contains four equations and six unknowns.
4.2 Model 2
In this model, the EOM of the 6RSS PM without HKP constraints is firstly written, then the constraints from two HKPs are modelled, to finally reach the model of the RAPM, where the model of the 6RSS PM is utilised as its core. However, it must be kept in mind that now the kinetic and gravitational components of generalized forces are functions of rather than . In the 6RSS PM, the step about the instantaneous power equivalence among the sum of the generalized forces on all the moving bodies, the input torques, and the bite force is as follows:
where , , , , and are the kinetic and gravitational components of the generalized force acting on the end effector, the ith coupler and the ith crank, respectively. They can be computed in an identical manner as their counterparts in the last section; thus, the process is omitted to save pages. After the two constraints from HKPs are exerted on the end effector, its DOFs are changed from X_{EE} to q_{EE}. Then using equations (8), (9), and (11), equation (44) can be rewritten as
where is the generalized force vector by the bite force F_{B}, J_{θ2} is the 6 × 6 Jacobian matrix mapping X˙_{EE} to θ˙, namely, θ˙=J_{θ2}⋅X˙_{EE}. Because q˙T_{EE} is free to vary, it can be rewritten as
which is the model of the RAPM with four equations and six unknowns. It must be emphasised that even though the generalized forces in the bracket are ultimately functions of , they are directly expressed and computed by . In this manner, the strong coupling in the parasitic motion variables Z and by q_{EE} is completely alleviated, reducing the computational complexity to a huge degree, which will be exhibited in Section 6.2.1.
One can easily find that, by substituting only equations (8) and (11) into equation (44), and “deleting” the independent vector X˙_{EE} produces
which is the EOM of the 6RSS PM without HKPs. That is to say, from the comparison between equations (46) and (47), the EOM of the RAPM can be easily generated by leftmultiplying at the two sides of that of the 6RSS PM, which is used as the core of the dynamic model of the RAPM.
5 LagrangeD'Alembert formulation
In writing the EOM of the RAPM by this formulation, there are two models to achieve the EOM in this section: in Model 1, the EOM is directly formatted to the RAPM; while in Model 2, firstly, the dynamic model of the 6RSS PM is built, and next, the two HKP constraints are modelled, to finally obtain the model of the RAPM of interest.
5.1 Model 1
In this model, this formulation is directly applied to the RAPM: at first, the mechanism is virtually cut at M_{i}, then the dynamics of the end effector and the chains is formulated in the task space and the joint space independently by the Lagrangian formulation; next, the complete dynamic model of the mechanism without constraints is built. Finally, the model of the RAPM is acquired via the D'Alembert formulation.
5.1.1 The end effector
The kinetic energy of the end effector is computed as
whereM_{EE}is the inertia dyad that has been given in equation (29), and M_{EE1} is its4 × 4 mass matrix that can be expressed as
Its potential energy has already been derived in equation (31). Regarding these, its Lagrange function is
One can find that
by virtue of which the Lagrangian formulation of the end effector is
where
See Formula in PDF
are the 4 × 4 Coriolis/centrifugal force matrix and the4 × 1gravitational force vector of the end effector, respectively.F_{EE1} (q_{EE}) is the4 × 1generalized force vector corresponding to the four DOFs of the end effector and has been computed in equation (42); evidently, it is a function of q_{EE}. Now the dynamic model of the end effector has been formulated straightforwardly in the task space.
5.1.2 The ith chain
The EOM of a serial chain by the Lagrangian formulation is well developed and can be found in many textbooks on robotics. As such, the details are not given and only the final EOM is written as
where and are themass matrix, theCoriolis/centrifugal force matrix, thegravitational force vector, and thegeneralized force vector of the ith chain, respectively, and is the actuating torque provided by the ith actuator. Note that the model is built in the joint space, and .
5.1.3 The overall mechanism
For the complete RAPM, the model without constraints can be written as
where
are the generalized coordinate vector, the mass matrix, the Coriolis/centrifugal force matrix, the gravitational force vector, and the generalized force vector of the six chains without constraints, respectively.
Following the procedure of the D'Alembert formulation, one can obtain that
where
and and are the virtual displacement vectors formed by q_{r} and q_{EE}, respectively. Because δq_{EE}is an independent vector, it can be found that
Substituting equations (17) and (18) into equation (56) and transposing the result, it can be equivalently written as
where
In equation (57), there are four equations and six unknowns, indicating the mechanism is under actuation redundancy.
5.2 Model 2
In this model, the EOM of the 6RSS PM without HKP constraints is written, as the function of X_{EE},X˙_{EE},X¨_{EE}. Then the constraints from HKPs are modelled, to finally reach the model of the RAPM, where the model of the 6RSS PM is incorporated as its core.
Likewise, the Lagrange equation of the 6RSS PM free of constraints among rigid bodies is written as
where and have identical structures as their counterparts in equation (54). However, it must be remembered that they are functions of X_{EE},X˙_{EE},X¨_{EE }rather than q_{EE},q˙_{EE},q¨_{EE}. Specifically, differentiating equation (15) with respect to time once and twice, and rearranging the results for the six chains prodocue
By the D'Alembert formulation, one can obtain
where
Substituting equation (59) into equation (60) and transposing the result give rise to
where
and J_{θ2} can be found in equation (46).
After the two HKP constraints are imposed to the end effector, its DOFs have transferred from X_{EE} to q_{EE}. Thus, from equation (9), one can obtain
Putting it into equation (61) and “deleting” the free term produce
This set of equations is the model of the RAPM with four equations and six unknowns.
From this process, there are two things should be emphasised: the first is, though in equation (61), the terms inside the bracket are functions of q_{EE},q˙_{EE},q¨_{EE }ultimately, the numerical values of X_{EE},X˙_{EE},X¨_{EE }computed from equations (6), (9), and (10) are directly used in them. That is to say, the model of the 6RSS PM is used as the core of the model of the RAPM. It can alleviate the coupling in the parasitic motion variables and boost the computational efficiency significantly, and this will be illustrated in Section 6.2.1. The second is, by “deleting” the term in ((61), one can find the following equation
is the EOM of the 6RSS PM. In this regard, the model of the RAPM can be formatted by leftmultiplying to the two sides of the EOM of the 6RSS PM, which is an identical manipulation as in Section 4.2.
6 Numerical computations and discussions
To verify the methods and models in this paper, the mechanism is commanded to follow a real incisor trajectory in free chewing by a healthy human subject which lasts 5 seconds. The trajectory from different perspectives is given in the first four subplots in Figure 5. The experimental setup and procedure to obtain the mastication movements from human subjects can be found in Chapter 6 of [3]. The corresponding mandibular motions in terms of the four elements of q_{EE} are computed, using the same method in equation (38) of [1]. Next, via equation (6), the numerical values of parasitic motion variables Z andγin terms of q_{EE} can be obtained. Finally, the first and second timederivatives of these six motion variables are also computed via equations (9) and (10). The values of X_{EE},X˙_{EE},X¨_{EE }versus time are provided in the following six subplots.
Four interesting features in this chewing trajectory in Figure 5 can be observed. Firstly, the sagittal view in the X_{S}Z_{S} plane is like a straight line, indicating the trajectory is nearly concentrated on one plane. Secondly, from the last six subplots, during the first one second, the motions of the mandible are quite small. After that period, it performs a rhythmic chewing motion with larger amplitudes of Z andβthan those of X, Y, α and γ, and the first and second timerates of Z andβare also larger than their counterparts, indicating the chewing motions are dominated by the mouth opening and closing movements. In the third place, the lateral movement indicated by the values of Y,Y˙,Y¨ are quite small. Finally, the trajectories of α and γ are nearly coincident in the time history, so are their first and second timerates. In following this predefined trajectory, numerical values of six diagonal entries of J_{A} in equation (26) are exhibited in Figure 6. It clearly shows that none of them is zero and as a result, equation (27) is not satisfied and an inverse singularity cannot occur. Meanwhile, det(J_{W}) always equals four in the time history and then equation (28) is invalid. In these regards, the planned trajectory is not passing the three sorts of singular configurations as mentioned in Section 3.3.
Concerning the 6RSS PM, its six DOFs and their first and second timederivatives as a function of time are identical as those in the final six subplots in Figure 5. As a consequence, in the two mechanisms, the end effectors undergo identical motions in six directions, and each of the six chains in them also undergoes identical motions in terms of the numerical values of . In tracking this prescribed chewing trajectory, the 6RSS PM is at or near singular configurations using the analysis of the RAPM as in Section 3.3. Correspondingly, an experimentally measured 3D bite force in {S} on peanuts by a healthy human subject on the molars in Figure 7, is exerted onto the right molar of the end effector. In this figure, there are in total five bursts in the three components of the bite force in the timeline, and each of the burst corresponds to the jaw closing portion for each chew. The size of the z component is larger than those in the other two directions, because the force in this direction plays a more important role in the food chewing process.
Fig. 6 Numerical values of the diagonal elements of J_{A}. 
6.1 Solutions of two methods
6.1.1 Generalized momentum approach
In both Models 1 and 2 by this approach, there are four equations and six unknowns in equations (43) and (46), indicating there are infinite solutions for the torque distribution theoretically. A closedform solution for the actuating torques can be obtained simply by the pseudoinverse of and , respectively, as
which corresponds to the minimum of actuating torques. From model 1, the input torquesτversus time along the chewing trajectory are exhibited in the first column of Figure 8, and the torque differences between Model 1 and Model 2 in this approach are given in the second column. From the first column, one can see that the five bursts produced by each actuator are tightly following those in the bite force in the time history in Figure 7, rather than those in the free chewing motion in Figure 5. It means the torque cost required by the bite force is larger than that by the inertia and gravitational forces of the mechanism itself. More importantly, in the second column the differences are quite minor, thus, the two models are equivalent in terms of the numerical results.
Fig. 8 Torques from Model 1 and differences between the two models by the generalized momentum approach. 
6.1.2 LagrangeD'Alembert formulation
Via the two models by the LagrangeD'Alembert principle in equations (57) and (64), the actuating torques can also be minimised using the pseudoinverse solution, respectively, as
The actuating torques from Model 1 and the differences between the actuating torques from the two models under this formulation are equivalent to those in Figure 8, respectively. As a consequence, they are not exhibited for the sake of brevity. Furthermore, the equivalent results demonstrate the correctness and accuracy of the developed four models by the two methods.
6.2 Discussion
6.2.1 Computational demands
The computational efficiency of inverse dynamics is crucial to the realtime application of the RAPM, and it is an important factor in judging the performance of the model. To assess the suitability of the built models for realtime control, reliable quantitative measurement of the computational load is critical. The time consumption is summarised in Table 1, where GM and LD are short for the generalized momentum approach and the LagrangeD'Alembert formulation, respectively, and the number 1 and 2 are meaning the first and the second model in the corresponding method, respectively. The procedures under each model are all divided into symbolic and numeric computations, and have been implemented in programs written in Matlab, using an Intel(R) Core(TM) i78700K CPU@3.70GHz and 64GB of RAM. By the symbolic EOMs, it is rather convenient to obtain the mathematical functions of the Jacobian matrix, the kinetic and gravitational components of the generalized forces, etc. in these methods, and they can be called in the numerical computations.
In the dynamic models where the two methods are directly utilised to the RAPM, i.e., in the first two columns in Table 1, the model under the newly developed generalized momentum approach is more computationally efficient than that under the classical LagrangeD'Alembert formulation by about 40% and 43.6% in the symbolic and numeric computations, respectively. The reason for the superiority of the generalized momentum approach roots in the high efficiency in deriving the Coriolis/centrifugal force matrix by the differentiation of the generalized momentum, as stated in [16,17]. In the two models where the EOM of the 6RSS PM is incorporated as the core of the model of the RAPM, namely, from the third and fourth columns, one can also reach the conclusion that the generalized momentum approach is faster than the LagrangeD'Alembert formulation.
From the comparisons between the first four columns in Table 1, it can be seen using the dynamics of the 6RSS PM can boost the efficiency to a huge degree. The total time consumption by the two methods in Model 2 is only about 22.5% and 16.2% of that in Model 1, respectively. The main reason for this considerable improvement is that the resulting complexity of Model 1 by each method is a consequence of using the complex symbolic expressions of Z and γ from equation (6), and the Jacobian matrix M_{J} from equation (9), for q_{EE},q˙_{EE},q¨_{EE} are used directly. By contrast, in Model 2, even though the matrices and vectors are ultimately functions of q_{EE},q˙_{EE},q¨_{EE}, they can be computed by the numerical values of X_{EE},X˙_{EE},X¨_{EE} directly, which can be available from q_{EE},q˙_{EE},q¨_{EE} at the beginning of the computation, and then are fed into equations (46) and (63). In this manner, the strong coupling between the two parasitic motion variables Z and γ and the four DOFs q_{EE} is completely alleviated. In other words, Model 2 is free of the complex symbolic computations brought by Z,γ, and M_{J}. As such, the symbolic mathematical functions, which are called and used in numerical computations, are less complex than those in Model 1, offering simpler possible computational algorithms. In this regard, a shortcut that can greatly boost the computational efficiency is found, by incorporating the dynamics model of the 6RSS PM as the core of that of the RAPM. On this basis, a leftmultiplication ofis simply performed to the two sides of EOMs of the 6RSS PM as in equations (47) and (64). Indeed, this shortcut is inspired intuitively by the comparison between the RAPM and the 6RSS PM: the former is generated by the two HKP constraints onto the end effector of the latter.
Computational time of the RAPM and the 6RSS PM (unit: s).
6.2.2 Influence of HKP constraints
It is remembered that the RAPM is built by imposing two HKP constraints onto the end effector of the 6RSS PM. Thereby, the dynamic model of the 6RSS PM free of HKPs is also built via the foregoing two methodologies, as far as the role of the HKP constraints in the numerical results and the computational cost is concerned. equations (47) and (64) are the EOMs of the 6RSS PM by these two methods, and the variables inside them are functions of X_{EE},X˙_{EE},X¨_{EE}, which are identical as in the final six subplots of Figure 5. By these two approaches, its inverse dynamics problem is reduced to solving a system of six linear algebraic equations in six unknowns. The input torques can be uniquely determined and there is no optimisation as in the RAPM. Their profiles in the time history are similar to those in the first column of Figure 6, and the only difference is the magnitudes of the torques in the 6RSS PM are larger than those in the RAPM. Thus, they are not exhibited to save pages.
Specifically, the contributions the total torques by the end effector, the six couplers, and the six cranks, respectively in the RAPM and the 6RSS PM are to be discovered and compared. in the RAPM, according to the EOM by the generalized momentum approach, the torques required by the abovementioned three parts are given in the three columns of Figure 9, respectively. Clearly, the bursts in each actuator are tightly tracking those in the freechewing trajectory as shown in the final six subplots in Figure 5. In addition, in one subplot, the torques offered by each pair of actuators are more or less equivalent. A very large portion of the total torques is required by the end effector, but those by the couplers and the cranks are quite minor. From the comparison between Figure 9 and the first column of Figure 8, it is obvious that the torque cost by the motion of the RAPM itself is quite small, but the torque required to sustain the bite force occupies a large portion. In the 6RSS PM, the torque distributions by the three parts are experiencing a similar profile as those in Figure 9 but with larger magnitudes, thus they are not exhibited to save pages.
To quantitatively justify and compare the influence by HKPs to the actuating torques, in the RAPM and the 6RSS PM, four indices are set about the numerical values of total actuating torques, and torques required by the aforementioned three parts as
where
and N is the number of sampling points in the time history. Obviously, the total actuating torques τ can also be obtained from Model 2 by the generalized momentum approach, or from the two models by the LagrangeD'Alembert formulation. can also be computed from Model 2 by the generalized momentum approach; however, by the LagrangeD'Alembert formulation, the contributions by each of the three parts to the total actuating torques are not as easily available as from the generalized momentum approach.
The numerical values of these indices of the RAPM and the 6RSS PM are given in the first and the second line of Table 2, respectively. The value of F is always 0.0207N.m from the RAPM under the four models; while in the 6RSS PM, F equals 0.0723N.m from two methods and is larger than that by the RAPM over three times. It means that redundant actuation can minimise the input torques, which is a wellestablished and wellknown opinion and has been proved in a number of publications, even though the way to generate redundant actuations in the mechanism is rather different from the two primary methods in [20]: they are caused by transferring two DOFs into parasitic motions under the HKP constraints. In other words, smaller torques are devoted to the remaining four DOFs, even though they can give rise to two parasitic motion variables. Besides, the torques required to maintain the motion of the aforementioned three parts increase by 72.32%, 52.94%, and 9.09% from the RAPM to the 6RSS PM, respectively. Thereby, both the proportion in the end effector and the torques required by it are the largest. In these regards, an important conclusion can be reached: even though the end effector undergoes identical motions in six directions in terms of displacements, velocities, and accelerations in both the RAPM and the 6RSS PM, the required torques are evidently smaller in the parasitic motion variables Z and γin the RAPM than those by the corresponding two DOFs in the 6RSS PM.
The computational time for the 6RSS PM under the two methods is given in the last two columns of Table 1. It is noted that the symbolic computational time is equivalent to that in Model 2 by each method in the RAPM, respectively, because one common file in Matlab for the symbolic computation is shared by the RAPM in Model 2 and the 6RSS PM. The numerical computational burden of the generalized momentum approach is also lower than that of the LagrangeD'Alembert formulation, which is an identical conclusion found in the RAPM. Meanwhile, as expected, the time used in the two methods in the 6RSS PM is comparable to that in Model 2 by the two methods, respectively, because the core of the model of the 6RSS PM is used in Model 2, and the difference is only a numerical leftmultiplication of .
By comparing the time consumption between the 6RSS PM and Model 1 from each method in the RAPM, where both the two methods are directly used to the two mechanisms, it is found that the two HKP constraints greatly increase the computational complexity. The reason why the RAPM consumes considerably more time is that the complex expressions of parasitic motion variables Z andγ, and the Jacobian matrix M_{J}, etc. caused by the HKP constraints greatly increase the demands in the symbolic computations. On this basis, more complicated numerical computations are processed. By contrast, in the 6RSS PM, Z andγare directly used as DOFs which are much more straightforward, and there is no Jacobian matrix M_{J}.
Fig. 9 Torques required by the end effector, the six cranks, and the six couplers in the RAPM. 
Four indices in the RAPM and the 6RSS PM (unit: N.m).
7 Conclusion
The inverse dynamics of a spatial RAPM constrained by two HKPs was solved systematically via the generalized momentum approach and the LagrangeD'Alembert formulation. The methods and models were thoroughly studied and compared. Generally, studies of inverse dynamics concentrate on RAPMs with only lower kinematic pairs. Thereby, the scientific contribution of this paper is the deep study of inverse dynamics in a spatial RAPM with both lower kinematic pairs and HKPs in terms of revealing the nature of HKP constraints. Specifically, the following conclusions can be drawn:
The two models in each method give rise to identical numerical results; however, by incorporating the dynamic model of the 6RSS PM as the core of that of the RAPM in Model 2 in each method, the computational efficiency has been enhanced remarkably as fast as that of the 6RSS PM without any simplification or loss of accuracy.
The two methods are equivalent in terms of the numerical results; nonetheless, the computational demands are greatly lower in the newly developed generalized momentum approach than in the classical LagrangeD'Alembert formulation.
By comparing the computational requirements between the RAPM and the 6RSS PM, it is discovered that the HKP constraints greatly raise the computational complexity.
Compared to the 6RSS PM with six DOFs, in the RAPM smaller torques are devoted to its four DOFs and two parasitic motion variables. In other words, the torques required by the parasitic motion variables Z and γin the RAPM are smaller than those by the corresponding DOFs in the 6RSS PM.
The models in this paper can provide a great potential to facilitate the realtime motion and/or force control of the RAPM. It is also expected to straightforwardly extend the presented process to a class of RAPMs whose end effectors are constrained by the base directly.
In future work, the new form of LagrangeD'Alembert formulation in [27] for dynamic modelling of PMs will be employed to this RAPM.
Nomenclature
To comprehend the paper conveniently, a list of symbols together with their definitions is provided in the following table at their appearance order.
{S} Inertia frame: {M} Frame established at the mass centre O_{M} of the end effector
G_{i}S_{i}(i=1,…,6): The ith crank
S_{i}M_{i}(i=1,…,6): The ith coupler
{C_{i}}: Frame established at the fixed point G_{i} of the ith crank
{N_{i}}: Frame attached at the mass centre E_{i} of S_{i}M_{i}
T_{j} (j=L, R): Left and right condyle ball centres
O_{S}O_{M}: 3◊1 position vector of O_{M} in {S}
: Rotation matrix from {S} to {M}
R_{X} (α),R_{Y} (β), andR_{Z} (γ): Rotation matrices about X_{M}, Y_{M}, and Z_{M} axes by Bryant anglesα, β, andγ, respectively
X_{EE}: 6◊1 vector grouping six motion variables of the end effector
X˙_{EE},X¨_{EE}: First and second timederivatives of, respectively
p_{i}(i=1,…,6): Parameters describing the shape of surfaces where condyle ball centres T_{L} and T_{R} slide on
X,Y, and Z : Parameters describing translations of the end effector
α,β, and γ: Bryant angles describing rotations of the end effector
q_{EE}: 4◊1 vector grouping four DOFs of the end effector
q˙_{EE},q¨_{EE}: First and second timederivatives of q_{EE}, respectively
ω_{EE}: Angular velocity of the end effector
M_{J}: 6◊4 Jacobian matrix between X_{EE} and q_{EE}
θ: 6◊1 vector containing the displacements of six active joints
β_{i}andγ_{i}: Euler angles expressing the rotation of the ith coupler
∥S_{i}M_{i}∥: Length of S_{i}M_{i}
O_{M}M_{i}andMO_{M}M_{i}: Vectors ofO_{M}M_{i}in {S} and {M}, respectively
O_{S}G_{i}: 3◊1 constant position vector of G_{i} in {S}
: Rotation matrix from {S} to{C_{i}}at the initial configuration of the mechanism
R_{Z} (θ_{i}): Rotation matrix about the axis byθ_{i}
∥G_{i}S_{i}∥: Length ofG_{i}S_{i}
: 3◊1 generalized vector specifying the configuration of the ith kinematic chain
M_{1i}: Jacobian matrix of O_{S}G_{i} + G_{i}S_{i} + S_{i}M_{i} with respect to
M_{2i}: Jacobian matrix ofO_{S}O_{M} + O_{M}M_{i} with respect toq_{EE}
: Angular velocity of S_{i}M_{i}
: Translational velocity of mass centre E_{i} of S_{i}M_{i}
: Jacobian matrix ofO_{S}E_{i}with respect to
: Angular velocity of G_{i}S_{i}
L (θ, q_{EE}): 6◊1 constraint vector to analyse singularities
J_{A}: Jacobian matrix betweenL (θ, q_{EE})and θ
J_{W}: Jacobian matrix betweenL (θ, q_{EE})andq_{EE}
q_{EE1}: Generalized momentum of the end effector in Model 1 under the generalized momentum approach
M_{EE}: 6◊6 inertia dyad of the end effector
m_{EE} and I_{EE}: Mass and inertia matrix of the end effector with respect to O_{M} and expressed in {S}, respectively
: Kinetic component of the generalized force acting on the end effector in Model 1 under the generalized momentum approach
P_{EE1}: Potential energy of the end effector in Model 1 under the generalized momentum approach
: Gravitational component of the generalized force acting on the end effector in Model 1 under the generalized momentum approach
: Generalized momentum of S_{i}M_{i} in Model 1 under the generalized momentum approach
: 6◊6 inertia dyad of S_{i}M_{i}
and : Mass and inertia matrix of S_{i}M_{i} with respect to O_{M} and expressed in {S}, respectively
: Kinetic component of the generalized force acting on S_{i}M_{i} in Model 1 under the generalized momentum approach
: Potential energy of S_{i}M_{i} in Model 1 under the generalized momentum approach
: Gravitational component of the generalized force in Model 1 under the generalized momentum approach
and : Mass and rotational inertia of G_{i}S_{i} around the axis of frame{C_{i}}, respectively
: Gravitational component of the generalized force of G_{i}S_{i} in Model 1 under the generalized momentum approach
τ: 6◊1 actuating torque vector by the six actuators
: 6◊1 wrench vector by the bite force F_{B}
J_{θ1}: 6◊4 Jacobian matrix mapping q˙_{EE} into θ˙
F_{EE1}: 4◊1 generalized force vector by F_{B} in Model 1 under the generalized momentum approach
and : Kinetic and gravitational components of the generalized force acting on the end effector, respectively, in Model 2 under the generalized momentum approach
and : Kinetic and gravitational components of the generalized force acting on the ith coupler, respectively, in Model 2 under the generalized momentum approach
and : Kinetic and gravitational components of the generalized force acting on the ith crank, respectively, in Model 2 under the generalized momentum approach
J_{θ2}: 6◊6 Jacobian matrix mapping X˙_{EE} into θ˙
T_{EE1}: Kinetic energy of the end effector in Model 1 under the LagrangeD'Alembert formulation
L_{EE1}: Lagrange function of the end effector in Model 1 under the LagrangeD'Alembert formulation
C_{EE1}(q_{EE},q˙_{EE}): 4◊4 Coriolis/centrifugal force matrix, gravitational force vector, and generalized force vector of the end effector, respectively, in Model 1 under the LagrangeD'Alembert formulation
: 3◊3 mass matrix, 3◊3 Coriolis/centrifugal force matrix, 3◊1gravitational force vector, and 3◊1 generalized force vector of the ith chain, respectively, in Model 1 under the LagrangeD'Alembert formulation
τ_{i}: Actuating torque provided by the ith actuator
δq_{r}, δq_{EE}, δX_{EE}: Virtual displacement vectors of q_{r}, q_{EE}, and X_{EE}, respectively
F, F_{EE}, , : Indices about the numerical values of total actuating torques, and torques required by the end effector, six couplers, and six cranks, respectively
RAPM: Redundantly actuated parallel mechanism
RSS: Revolutesphericalspherical
PUS: Prismaticuniversalspherical
Conflict of interest
The authors declare that they have no conflict of interest.
Data availability
The datasets generated during and/or analysed during the current study are available from the corresponding author on reasonable request.
Data availability
This work was supported by the Natural Science Foundation of Shaanxi Province (Grant No. 2019JQ467).
Acknowledgements
The authors also sincerely appreciate the valuable time, suggestions and comments from the editors and the anonymous reviewers for improving the quality of the paper.
References
 C. Cheng, W. Xu, J. Shang, Kinematics, stiffness and natural frequency of a redundantly actuated masticatory robot constrained by two pointcontact higher kinematic pairs, in IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS), Hamgburg, Germany (2015), pp. 963–970 [Google Scholar]
 C. Cheng, W. Xu, J. Shang, Optimal distribution of the actuating torques for a redundantly actuated masticatory robot with two higher kinematic pairs, Nonlinear Dyn. 79, 1235–1255 (2014) [Google Scholar]
 W. Xu, J. Bronlund, Mastication Robots: Biological Inspiration to Implementation (SpringerVerlag, Berlin, Germany, 2010) [CrossRef] [Google Scholar]
 J. Wu, X. Chen, L. Wang, X. Liu, Dynamic loadcarrying capacity of a novel redundantly actuated parallel conveyor, Nonlinear Dyn. 78, 241–250 (2014) [CrossRef] [Google Scholar]
 V. van der Wijk, S. Krut, F. Pierrot, J. Herder, Design and experimental evaluation of a dynamically balanced redundant planar 4RRR parallel manipulator, Int. J. Robotics Res. 32, 744–759 (2013) [CrossRef] [Google Scholar]
 J. Hernández, A. Chemori, H. Sierra, J. Anieva, A new solution for machining with RAPKMs: modelling, control and experiments, Mech. Mach. Theory 150, 103864 (2020) [CrossRef] [Google Scholar]
 A. Müller, Problems in the control of redundantly actuated parallel manipulators caused by geometric imperfections, Meccanica 46, 41–49 (2011) [CrossRef] [MathSciNet] [Google Scholar]
 W. Shang, S. Cong, Robust nonlinear control of a planar 2DOF parallel manipulator with redundant actuation, Robot. Comput. Integr. Manufactur. 30, 597–604 (2014) [CrossRef] [Google Scholar]
 D. Liang, Y. Song, T. Sun, G. Dong, Optimum design of a novel redundantly actuated parallel manipulator with multiple actuation modes for high kinematic and dynamic performance, Nonlinear Dyn. 83, 631–658 (2015) [Google Scholar]
 Y. Zhao, F. Gao, Dynamic performance comparison of the 8PSS redundant parallel manipulator and its nonredundant counterpartthe 6PSS parallel manipulator, Mech. Mach. Theory 44, 991–1008 (2009) [CrossRef] [MathSciNet] [Google Scholar]
 J.V. Fontes, M.M. da Silva, On the dynamic performance of parallel kinematic manipulators with actuation and kinematic redundancies, Mech. Mach. Theory 103, 148–166 (2016) [CrossRef] [Google Scholar]
 S. Jin, J. Bae, J. Kim, T. Seo, Disturbance compensation of a dualarm underwater robot via redundant parallel mechanism theory, Meccanica 52, 1711–1719 (2017) [CrossRef] [MathSciNet] [Google Scholar]
 Z. Liu, R. Tao, J. Fan, Z. Wang, F. Jing, M. Tan, Kinematics, dynamics, and load distribution analysis of a 4PPPS redundantly actuated parallel manipulator, Mech. Mach. Theory 167, 104494 (2022) [CrossRef] [Google Scholar]
 X. Liu, J. Yao, Q. Li, Y. Zhao, Coordination dynamics and modelbased neural network synchronous controls for redundantly fullactuated parallel manipulator, Mech. Mach. Theory 160, 104284 (2022) [Google Scholar]
 J. Wu, J. Wang, L. Wang, T. Li, Dynamics and control of a planar 3DOF parallel manipulator with actuation redundancy, Mech. Mach. Theory 44, 835–849 (2009) [CrossRef] [MathSciNet] [Google Scholar]
 A. Lopes, Dynamic modeling of a Stewart platform using the generalized momentum approach, Commun. Nonlinear Sci. Numer. Simul. 14, 3389–3401 (2009) [CrossRef] [Google Scholar]
 A. Lopes, F. Almeida, The generalized momentum approach to the dynamic modeling of a 6dof parallel manipulator, Multibody Syst. Dyn. 21, 123–146 (2009) [CrossRef] [MathSciNet] [Google Scholar]
 H. Cheng, Y. Yiu, Z. Li, Dynamics and control of redundantly actuated parallel manipulators, IEEE/ASME Trans. Mechatr. 8, 483–491 (2003) [CrossRef] [Google Scholar]
 Y. Lou, Z. Li, Y. Zhong, J. Li, Z. Li, Dynamics and contouring control of a 3DoF parallel kinematics machine, Mechatronics 21, 215–226 (2011) [CrossRef] [Google Scholar]
 C. Gosselin, T. Schreiber, Redundancy in parallel mechanisms: a review, Appl. Mech. Rev. 70, 010802 (2018) [CrossRef] [Google Scholar]
 C. Cheng, H. Liao, Elastodynamics of a spatial redundantly actuated parallel mechanism constrained by two higher kinematic pairs, Meccanica 56, 515–533 (2021) [CrossRef] [MathSciNet] [Google Scholar]
 C. Cheng, B. Liu, Y. Li, Y. Wang, Elastodynamic performance of a spatial redundantly actuated parallel mechanism constrained by two pointcontact higher kinematic pairs via a model reduction technique, Mech. Mach. Theory 167, 104570 (2022) [CrossRef] [Google Scholar]
 M. Mesnard, J. Coutant, M. Aoun, J. Morlier, M. Cid, P. Caix, Relationships between geometry and kinematic characteristics in the temporomandibular joint, Comput. Methods Biomech. Biomed. Eng. 15, 393–400 (2012) [CrossRef] [PubMed] [Google Scholar]
 J.H. Koolstra, Dynamics of the human masticatory system, Crit. Rev. Oral Biol. Med. 13, 366–376 (2002) [CrossRef] [PubMed] [Google Scholar]
 S. Siegler, R. Hayes, D. Nicolella, A. Fielding, A technique to investigate the threedimensional kinesiology of the human temporomandibular joint, J. Prosthetic Dent. 65, 833–839 (1991) [CrossRef] [Google Scholar]
 L. Tsai, Mechanics of Serial and Parallel manipulators (John Wiley & Sons Inc, Canada, 1999) [Google Scholar]
 R. AboShanab, Dynamic modeling of parallel manipulators based on LagrangeD'Alembert formulation and Jacobian/Hessian matrices, Multibody Syst. Dyn. 48, 403–426 (2020) [CrossRef] [MathSciNet] [Google Scholar]
Cite this article as: C. Cheng. X. Yuan. Y. Chen. Y. Dang, Rigidbody inverse dynamics of a spatial redundantly actuated parallel mechanism constrained by two point contact higher kinematic pairs, Mechanics & Industry 24, 28 (2023)
All Tables
All Figures
Fig. 1 A schematic view of the RAPM constrained by two point contact HKPs, where ① right condyle ball, ② left condyle ball, ③ articular surface of right TMJ, ④articular surface of left TMJ [21]. 

In the text 
Fig. 2 The RAPM: (a) CAD model, (b) magnification of the right HKP, (c) prototype of the HKP [2], (d) end effector. 

In the text 
Fig. 3 Prototype of the mechanism [2]. 

In the text 
Fig. 4 The structure of dynamic modelling by the two methods in Sections 4 and 5. 

In the text 
Fig. 5 Motions of the end effector in 3D space and time history [22]. 

In the text 
Fig. 6 Numerical values of the diagonal elements of J_{A}. 

In the text 
Fig. 7 3D bite force profiles on peanuts [3]. 

In the text 
Fig. 8 Torques from Model 1 and differences between the two models by the generalized momentum approach. 

In the text 
Fig. 9 Torques required by the end effector, the six cranks, and the six couplers in the RAPM. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.