Principal vectors and equivalent mass descriptions for the equations of motion of planar four-bar mechanisms

The use of principal points and principal vectors in the formulation of the equations of motion of a general 4R planar four-bar linkage is shown with two kinds of methods, one that opens kinematic loops and one that does not. The opened kinematic loop approach analyses the moving links as a system with a tree connectivity, introducing reaction forces for closing the loops. Compared with the conventional Newton–Euler method, this approach results in fewer equations and constraint forces, whereas the mass matrix entries remain meaningful, but there is a stronger coupling between the equations. Two equivalent mass formulations for the closed kinematic loop approach are presented, which need not open the loop and introduce loop constraint forces in the equations of motion. With the method of complex joint masses, the mass of the links closing the loops is represented by real and virtual equivalent masses, deﬁning the principal points. The principle of virtual work with the inclusion of inertia terms is used to derive the equations of motion. As an example the dynamic balance conditions of the four-bar linkage are derived. With the method of the equivalent mass matrix it is shown how a constant mass matrix can be used to describe the dynamics of binary links with an arbitrary mass distribution. A four-bar linkage could be modelled with only three truss elements instead of the conventional result with three or more beam elements, which gives a signiﬁcant reduction of the computational complexity.


Introduction
This article presents a sequence of three ways to derive the equations of motion of planar linkages. The planar 4R four-bar mechanism is used as an example. It is intended to find methods which present the equations of motion in a meaningful and insightful way such that the influence of varying the individual design parameters on the dynamics can be grasped. Planar 4R four-bar mechanism with its parameters along A 0 A 3 , and its y-direction is obtained by rotating the x-direction by a right angle in the positive direction.
The relative and absolute positions of points will be represented by two-dimensional vectors, which are indicated by an overbar. They will also be considered complex numbers, where the absolute value is the length of the vector, denoted by the same symbol without an overbar, and the argument is the angle of the positive direction of the vector with the x-axis. For instance, the link lengths and their orientations are denoted byl k = l k exp(iθ k ) = l k cos θ k + il k sin θ k , where θ k is the angle of link k with the global x-direction as shown in Fig. 1. For the fixed link 4, θ 4 = π .
The scalar product is defined as the usual inner product of two vectors, for which a notation with angle brackets is chosen. Ifx andȳ are two complex numbers, the inner product of the corresponding vectors can be expressed as x,ȳ = 1 2 (x * ȳ +ȳ * x ) = Re(x * ȳ ), where a superscript asterisk denotes a complex conjugate and Re indicates the real part.
With the above notation, the loop-closure conditions can be expressed as where0 is the zero vector or the complex zero. The angle θ 1 is chosen as the independent coordinate and the angles θ 2 and θ 3 and all position coordinates of the joints and the centres of mass of the links are dependent coordinates. Several ways to determine the dependent coordinates can be found in the literature; see, for instance, [14]. Here, we use a way, similar to a method presented in [15], that limits the introduction of trigonometric functions and their inverses in the equations and is numerically stable if the kinematic problem itself is stable. Furthermore, it makes the choice between the two possible configurations simple. With the vectorl 23 =l 2 +l 3 = −l 1 −l 4 having the length l 23 = (l 4 − l 1 cos θ 1 ) 2 + (l 1 sin θ 1 ) 2 = l 2 1 + l 2 4 − 2l 1 l 4 cos θ 1 , we introduce the additional lengths d 23 and h 23 shown in Fig. 1 as The triangle inequalities |l 2 − l 3 | ≤ l 23 ≤ l 2 + l 3 must be satisfied to give meaningful results. The combination of signs comes from the two possible configurations: the plus sign corresponds to the case in which A 2 is to the left of the directed line A 1 A 3 , as drawn in Fig. 1, and the minus sign corresponds to the case in which this point is to the right of this line.
If it is assumed that the Grashof condition is satisfied and link 1 is the crank, no branching between these two configurations can occur. The solution becomes inherently sensitive to perturbations in the data if the length l 23 is small or if l 23 is close to one of the bounds of the triangle inequalities given above. The dependent link orientation angles θ 2 and θ 3 are obtained froml 2 = (d 23 + ih 23 )l 23 l 23 ,l 3 =l 23 −l 2 (4) as the arguments ofl 2 andl 3 .

Velocity and acceleration analysis
The rates of the dependent angles are obtained from the time derivative of the loop-closure condition (1), The dependent angular velocities are obtained by taking the inner products of this equation withl 3 andl 2 , which results inθ with the first-order kinematic transfer functions Taking another time derivative of the loop-closure condition yields i(l 1θ1 +l 2θ2 +l 3θ3 ) − l 1θ

Fig. 2
System of principal vectors which describes the motion of the link masses with respect to the system centre of mass in S, with principal points P 1 , P 2 and P 3 The expressions for the dependent angles and angular velocities can be substituted in the right-hand sides to obtain explicit results. In most cases, the dependent angles and their time derivatives will be kept in the equations, where it has to be understood that they are functions of the independent angle and its time derivatives. Another form to write the angular velocity and acceleration of link 2 that only contains the angles θ 1 and θ 3 and their derivatives is Use has been made of the loop closure condition (1) and its time derivatives (5) and (8).

Opened kinematic loop approach
In this section, the dynamics of the four-bar linkage are analysed by the method of principal vectors in the way Fischer used them [2]. While Fischer used them only for deriving Lagrange's equations of motion of open trees of links, here the method of principal vectors is applied in the analysis of a closed-loop linkage, where the principle of virtual work is used to obtain the equations of motion. The advantage of using the principle of virtual work is that the equations of motion are often found in a simpler way and it leads to more insight in how parameters contribute to the equations of motion. In Fig. 2, the four-bar mechanism is shown with a system of principal vectorsā 12 ,ā 21 , a 23 andā 32 . These vectors fix the locations of the principal points P k (k = 1, 2, 3) within the Fig. 3 Moving links as separate system of three links in series with the principal vectors transformed into links which trace the system centre of mass in S links as illustrated. These principal points, also known as barycentres [4], are fixed within the links and are the centres of mass of the reduced systems of the links [1,2,5,6]. In this case, the reduced system of link 1 has the mass of link 1, m 1 , in its centre of mass S 1 , a moment of inertia I 1 about its centre of mass and the mass of link 2, m 2 , together with the mass of link 3, m 3 , placed in A 1 . The reduced system of link 2 has the mass m 2 in its centre of mass S 2 , a moment of inertia I 2 about its centre of mass, the mass m 1 placed in A 1 and the mass m 3 placed in A 2 . The reduced system of link 3 has the mass m 3 in its centre of mass S 3 , a moment of inertia I 3 about its centre of mass and a mass m 1 together with a mass m 2 placed in A 2 . The principal point P 1 is the centre of mass of the reduced system of link 1, the principal point P 2 is the centre of mass of the reduced system of link 2, and the principal point P 3 is the centre of mass of the reduced system of link 3.
The principal vectors form parallelograms with the linkage and their lengths are constant, equal to the distances between the principal points and the link joints. This means that the graphical construction can be regarded as a linkage tracing the system centre of mass S during its motion along with the four-bar linkage. This is shown in Fig. 3, where the vectors have been transformed into links with revolute pairs at the principal points, the centre of mass and the two intermediate points B 1 and B 2 .
In the opened kinematic loop approach, the closed kinematic loop of the mechanism is opened at one or more joints to form a system with a tree connectivity. In particular, the mechanism is removed from its supports to form a free system. The constraint forces are included as additional variables. The number of variables to describe the system and hence the number of resulting equations can be reduced by reimposing one of the supports. The free system for the four-bar mechanism is shown in Fig. 3. At the joints A 0 and A 3 the constraint forcesR A0 andR A3 are included as additional variables. The positions of the link centres of mass S k with respect to the link coordinate system defined by the vectorl k and the first joint as the origin is given by the vectors k = s k exp(iσ k ) (k = 1, 2, 3), as shown in Fig. 1. This means that the absolute vector pointing from the first joint to the centre of mass is expressed as a complex number ass klk /l k . For this free system of planar bodies interconnected by revolute joints in a tree, the generalized coordinates can be chosen as the position coordinates of the point S and the three rotation angles of the links about the principal points. The individual rotations leave the point S invariant and owing to the definition of the principal points, leave the centre of mass of the system invariant, which is the point S. This will be confirmed by explicitly calculating the position of the centre of mass.
In vector notation, the principal vectors can be calculated from with m tot = m 1 + m 2 + m 3 . These equations are comparable to the notation without vectors in [5,6]. The derived vectorsā 14 andā 34 are calculated as Note that the vectorsl k can be expressed as the difference of two principal vectors or derived vectors associated with the same link, so they are fixed in this link. That S is the centre of mass of the mechanism can also be seen from the expressions for the positions of the centres of mass of the links in terms of the position vector of S, which is denoted byr S , and the principal vectors. As the positions of the centres of mass in the links can be expressed in the principal vectors and the derived vectors, we have Indeed, m 1r1 + m 2r2 + m 3r3 = m totrS . If the constraint for joint A 0 , is reimposed and used to eliminate the position of the centre of massr S from the equations, the position vectors of the centres of mass of the links become By taking the time derivatives of these equations two times, the accelerations of the centres of mass of the links, multiplied by their respective masses, are found to be The virtual displacements of the centres of mass are easily found from and the virtual displacements of the joint A 3 , considered a point of link 3, as The corresponding reaction force of the linkage on the ground at A 3 isR A3 . The virtual work of the inertia terms, the applied forces and moments, and the constraint forces are given by where M k is the resultant applied moment with respect to the centre of mass of link k and F k is the resultant applied force, represented as a complex number, as illustrated in Fig. 1. Substituting the quantities in this equation, while keeping the link angles provisionally as independent variables, one obtains three equations of motion for the three orientation angles of the links of the form In this equation, the coefficients of the symmetric mass matrix, M jk , are given by The diagonal terms are the moments of inertia of the reduced systems of the links about their first joint and, hence, they are constant. The off-diagonal terms M ij (i = j) have a simple structure if they are expressed in terms of principal vectors: they are minus the total mass multiplied by the inner product of the link vector l i and the principal vector of link j pointing towards the link. They are zero if the two vectors are perpendicular to each other. The skew-symmetric coefficients of the convective inertia terms, C jk , are given by There are no couplings in the angular velocities in the convective inertia terms because the absolute rotation angles of the links are taken as coordinates, so the acceleration of a vector only contains the angular velocity of the link to which this vector belongs. The skew symmetry comes from the fact that the vectors in the inner products change place, but the factor i remains at the first term, so, because for general vectorsx andȳ we have ix,ȳ = − iȳ,x , the terms change their sign and the matrix becomes skew-symmetric. The coefficients for the reaction forces, D jx and D jy , are The generalized applied forces, Q j , are given by These are the moments on the links about their first joint with the forces on the links connected to the second joint shifted towards this second joint. The constraint equations to enforce the connection of the joint A 3 to the ground arē Differentiating these constraints two times with respect to time yields the constraints on the accelerations,l The five scalar equations (22) and (28) form a system of differential-algebraic equations that can be solved for the three unknown accelerations, that is, the three angular accelerationsθ k (k = 1, 2, 3), and the two components of the reaction forces, R A3x and R A3y , as these equations are linear in these five unknowns. Once the three angles of the links and their velocities and accelerations are known, together with the reaction force at A 3 , the reaction force at A 0 can be found from the equation of motion of the centre of mass, wherer S follows from Eq. (16) andR A0 is the force of the mechanism on the foundation at A 0 as shown in Fig. 1.
Alternatively, a single equation of motion can be obtained by multiplying the second equation of (22) by Θ 2 and the third equation by Θ 3 and adding these to the first equation, which eliminates the constraint forces and yields the equation of motion of the mechanism in θ 1 if the expressions for the dependent quantities are substituted. The reaction forces can then be found by substituting the accelerations in the original differential-algebraic equations.

Closed kinematic loop approach with equivalent mass formulations
In this section, a different way of using the principal vectors is shown and two methods to represent the mass distribution are presented. In the first method, the linear momentum of the centre of mass of a link is replaced by two linear momenta at the joints by means of a complex mass representation, which do not contribute to the moment of inertia. In the second method, the mass is represented by a mass matrix, which is only related to the translations of the joints. The methods have in common that the kinematic loop is not opened and no explicit constraint forces appear in the equations.

Method of complex joint masses
Let us consider a general binary planar link k, which can be any link of the four-bar linkage in the present context, with joints p and q having the respective position vectorsr p andr q . The vector connecting the joints isl k =r q −r p with a length l k and the position of the centre of mass relative tol k andr p is given bys k , as in Fig. 1, so the absolute position of the centre of mass can be written asr k =r p +s klk /l k . If we multiply the position of the centre of mass with the mass of the link, m k , we obtain an expression for the centre of mass of a link in terms of the positions of the two joints at its ends, For convenience, the position of the centre of mass of the link is written in terms of two distances ass k = e k + if k , where e k is the real part, or the distance along the line connecting the two joints, and f k from the imaginary part is the distance with the appropriate sign to the line connecting the two joints of the link, as is shown for link 2 in Fig. 3. The result in Eq. (30) can be interpreted as the expression for determining the centre of mass from two − im k f k /l k andm kq = m ksk /l k = m k e k /l k + im k f k /l k , whose sum is m k . This representation by equivalent complex masses was first proposed by Freudenstein [8]. The application to shaking force balancing was later called the complex mass method [9]. In the case that the centre of mass of the link is on the line connecting the two revolute joints, so f k = 0, these are the two usual real equivalent masses, one in each joint. In the general case we have a real mass m a k = m k (1 − e k /l k ) and an imaginary mass −im c k = −im k f k /l k at node p and a real mass m b k = m k e k /l k and an imaginary mass im c k at node q, as was introduced in [5,6]. There, m a k and m b k were called real equivalent masses and m c k a virtual equivalent mass. With these three equivalent masses, a graphical construction of modelling the masses of open-and closed-loop kinematic chains was presented and applied for the determination of the principal points. The linear momentum can be found by differentiating Eq. (30) with respect to time and can be interpreted as consisting of two contributions of the equivalent complex masses at the joints. The real part of the complex mass gives a linear momentum in the same direction as the velocity vector at the joint, but the imaginary part gives a linear momentum that is perpendicular to the direction of the velocity vector.
It should be stressed that no equivalence for the angular momentum can exist if f k is not zero: a velocity v of the centre of mass in the direction of the link vector gives rise to an angular momentum for the equivalent masses, 2m k f k v, which cannot be compensated for by adjusting the moment of inertia of the link, so there is no dynamical equivalence. In the special case that f k = 0, a reduced moment of inertia I k,red = I k − m k e k (l k − e k ) gives a dynamically equivalent model.
In Fig. 4, a new way to apply the principal vectors for the analysis of the four-bar linkage is presented. Here the principal vector system is located about the frame, including links 1, 3 and 4, and the graphical construction traces the centre of mass, S, of the combination of links 1, 2 and 3 for all motions of the four-bar linkage. The principal points in this case are P 1 , P 3 and P 4 , of which the first two are not the same as their namesakes in Figs. 2 and 3. Also the principal vectorsā 14 ,ā 41 ,ā 43 andā 34 are a new set. The locations of the principal points and the principal vectors can be calculated by using the complex equivalent mass description. As explained, the coupler mass m 2 can be redistributed by an equivalent complex mass m 2 (1 −s 2 /l 2 ) = m a 2 − im c 2 in A 1 and another equivalent complex mass m 2s2 /l 2 = m b 2 + im c 2 in A 2 . The projection of these masses on the links 1, 3 and 4 was presented graphically in [5,6] and is shown in Fig. 5. Here a real equivalent mass m a 2 is located in A 1 , a real equivalent mass m b 2 is located in A 2 and a virtual equivalent mass m c 2 is located twice about each principal point as illustrated. With the methods fully explained in these references, the principal points are found as the centres of mass of the three resulting reduced mass models in a way that is similar in form to the method without equivalent masses. The results can be rewritten in terms of complex vectors as Clearly,l 4 =ā 41 −ā 43 . The derived vectorsā 12 andā 32 are obtained from The centre of mass of the system, S, can be found from the parallelogram construction shown in Fig. 4. The correctness of this construction can easily be seen by writing the position of the centre of mass of the system starting from one of the joints or principal points. The position of the centre of mass of the linkage S in terms of the principal vectors isr S = −ā 14 −ā 41 −ā 34 .
This can be compared with Eq. (16), where the same vector is expressed in terms of a different set of vectors. The vectors locating the centres of mass multiplied by their respective masses of the movable links expressed in the principal vectors and the derived vectors become From Eq. (34), the virtual displacements are found as Two differentiations with respect to time give the accelerations, The equations of motion are again obtained from the virtual work expression (21), without the explicit contribution of the constraint forces, which leads to the equations of motion having the form where Q c j are the generalized constraint forces, which will be eventually eliminated. The non-zero entries of the mass matrix now become (s 3 The only two non-zero terms in C jk are C 13 and C 31 = −C 13 , Note that there are off-diagonal terms in the mass matrix and the convective inertia terms, so the replacement of the mass of link 2 by complex masses at the nodes does not necessarily lead to a full dynamic equivalence. The generalized applied forces are given by It is seen that Q 1 is the moment of the forces on link 1 with respect to the joint A 0 , including the force on link 2 that has been distributed over the joints A 1 and A 2 ; likewise, Q 3 is the moment of the forces on link 3 about the joint A 3 , including the contribution of the distributed force of link 2. The reduced equation of motion, in which the constraint forces have been eliminated, are obtained by multiplying the second equation of (37) by Θ 2 , the third equation of (37) by Θ 3 and adding the resulting equations. Also, the relations (11) can be used to remove the dependency on θ 2 in a first step of the reduction and then to combine the two remaining equations to a single equation of motion. Finally, the expressions of the dependent angles, angular velocities and angular accelerations can be substituted to obtain a single secondorder differential equation in θ 1 , the degree of freedom.
The virtual displacements used for the calculation of the reactions are chosen as a virtual displacement of the centre of mass, δr S , and a rotation of the whole system, without changing its relative configuration, about the point S, denoted by δθ S . The corresponding virtual displacements of the centres of mass of the bodies are The resultant reaction force and moment are now obtained from the virtual work equation with these variations and the calculated accelerations.

Method of the equivalent mass matrix
Instead of expressing the positions of the centres of mass of the links in terms of the principal vectors, in this section equivalent mass matrices for the links are derived. For this purpose, the constraint that the distance between the two joints of a link is constant is released. This means that the links are considered to be pseudorigid bodies [16], which are allowed to undergo a uniform dilatation in the plane. A body of this kind can be approximately realized in so-called auxetic metamaterials with a negative value of Poisson's ratio close to −1 [17]. The joints of the binary link as described above p and q have coordinatesr p andr q , respectively, and define local coordinate axes directed from p to q and in the perpendicular direction with their origin at p. The position of a material point of the link is described by the dimensionless coordinates ξ and η, scaled with the distance l between the points p and q. In a way similar to Eq. (30), the position of a point in the plane of the link can be expressed as (the index k for a specific link is omitted here) The position of the centre of mass of the link is given by ξ = e/ l, η = f/ l, its mass is m, and its moment of inertia with respect to its centre of mass is I . The mass is independent of the deformation, but the moment of inertia scales with the square of l, so I/l 2 is constant. The velocity and acceleration of the point, as well as its virtual displacement, are given bẏ The virtual work due to the inertia terms can be determined by evaluating the integral − δr,r dm, which yields − Re[δr px , δr py , δr qx , δr qy ] 6 Single truss element (a) and the example four-bar mechanism (b) and the components of the force vector F r are given by By taking appropriate combinations of variations, the reaction forces, or their resultant force and moment, can be obtained.

Truss element
The mass matrix in Eq. (46) can be used to specify the mass properties of a planar truss element. This element as illustrated in Fig. 6(a) has two nodal points, p and q, with their Cartesian positions as the nodal coordinates, x = [r px , r py , r qx , r qy ] T , and its elongation as its only generalized strain, ε 1 = l − l 0 = (r qx − r px ) 2 + (r qy − r py ) 2 − l 0 , where l 0 is the reference length. This truss element was applied in the context of multibody system dynamics by van der Werff [18]. This element with the present extended mass matrix allows us to model the four-bar system with only three truss elements instead of three (or more) beam elements, which improves the efficiency of the calculations. However, the applied forces and moments have to be replaced by nodal forces as described in Eq. (47). The equations of motion are fully equivalent to Eq. (48). The rotation angles θ are not directly available as results of calculations, but they can be extracted from the nodal coordinates of the element by noting cos θ = (r qx − r px )/ l, sin θ = (r qy − r py )/ l and using the two-argument arctangent function. However, in the modelling, a massless rigid beam element is laid over the truss element representing the crank to be able to choose the Fig. 7 Angles of the links of the example four-bar mechanism crank angle as the degree of freedom. This adds one coordinate and two generalized strains to the model. In total, there are nine coordinates (the eight Cartesian coordinates of the joints and the rotation angle of the beam), of which four are fixed to the ground, five generalized strains (one for each truss element and two for the beam element of which one is dependent and can be ignored), and one degree of freedom. For a model with rigid beams only, six rigid beam elements are needed, one for each movable link and one for connecting the centre of mass to each link with mass. This leads to a system with 17 coordinates (14 Cartesian coordinates and three angles), of which four are fixed, twelve generalized strains (two for each beam element) and one degree of freedom. It is seen that the model with the truss elements is smaller and simpler.
As an example to verify the equations, dimensions and masses of a four-bar mechanism are chosen as (see Fig. 6 There is no gravity.
The system was modelled in the multibody dynamics program SPACAR [11] and the motion was simulated, starting from a rest position with θ 1 = 1 rad over 5 s. The results for the link angles are shown in Fig. 7 and the reaction forces at the support points are given in Fig. 8. For 501 values of the crank angle and the crank angular velocity obtained by the simulation, the accelerations and reaction forces as obtained by the different methods were checked. All results were in agreement up to truncation errors in the numerics. The initial and final values of the link angles, angular velocities and angular accelerations are given in Table 1.

Shaking force balance and shaking moment balance conditions
The conditions for dynamic force balance can be found from the open-loop approach and the principal vectors obtained in it [5], but they become almost trivial if the parallelogram construction for the system with equivalent complex masses in Fig. 4 is considered. The centre of mass S of the linkage remains stationary with respect to the ground link 4 if the principal vectorsā 14 andā 34 are zero. This gives the conditions or written out in scalar form, If the masses of all links are replaced by equivalent masses, these conditions mean that the total equivalent mass at the moving joints must be zero. On the other hand, the resultant reaction forces or shaking forces can be found from the equations with equivalent mass matrices by taking the same virtual displacements for all Fig. 9 Components of the reaction forces at the support points for the example balanced four-bar mechanism The centre of gravity can easily be found, and hence the linear momentum of the linkage as a whole and of the individual links. Also the sum of the reaction forces on the ground can be easily obtained, which is advantageous for shaking force balancing. A possible disadvantage is that there are still more equations than degrees of freedom, which leads to a system of differential-algebraic equations.
If the mass of one or more links is represented by equivalent complex masses at the joints as explained in Sect. 4.1, the centre of mass and hence the linear momentum of the system can easily be found; the conditions of shaking force balancing (58) for the four-bar linkage become almost trivial. The mass description remains relatively simple. A disadvantage is that the equivalent masses are not fully dynamically equivalent, so the equations of motion cannot be derived by applying standard techniques to the equivalent masses only.
Principal points and principal vectors are useful for the opened kinematic loop approach, as the off-diagonal terms of the mass matrix and the convective inertia terms are easily obtained. To a lesser extent, this applies also for the method with equivalent complex masses, although this is not so obvious in the considered simple example. In both approaches, the motion of the centre of mass of the system is easily found with the principal vectors, which yields conditions for shaking force balance.
The description with equivalent mass matrices in Sect. 4.2 includes the complete dynamics, although the mass matrix is filled and does not have a block diagonal structure. The equations of motion can initially be described in terms of the joint coordinates, without making use of the angles of the links. The independent coordinate can be chosen as one of the angles or one of the position coordinates of the joints. This mass description can be used for the finite truss element, leading to simpler models in terms of the kind of element used and the number of coordinates and constraints.
The methods presented can be applied to other linkages, with similar advantages and disadvantages, although the equations can become more complicated. General planar mechanisms with links interconnected by pin joints can be directly treated by the methods shown here. An example is the seven-body mechanism form [20], which was modelled by truss elements and simulated in another publication [21].
Principal points and principal vectors can also be defined for systems of spatial bodies interconnected with spherical joints [2,4] and they can be used to obtain the equations of motion. Some other types of joints, such as a revolute joint, a universal joint or a homokinetic coupling, can be described by adding constraints on the relative motion in a spherical joint. The method of equivalent masses can be extended to some special spatial systems. In general, a binary link between spherical joints cannot be modelled in this way, as the orientation of this link is not determined by the positions of two its joints and even the linear momentum is generally not determined by the velocities of these two joints. Six classical spatial truss elements which interconnect four points can be used to describe the dynamic properties of a rigid body [21].

Conclusions
This article has led us from a familiar approach with opened kinematic loops via a method with equivalent complex masses to a method using a constant mass matrix, which were all applied to formulate the equations of motion of a general 4R planar four-bar mechanism. In the first two methods, the use of principal vectors has been shown. The principle of virtual work with the inclusion of inertia terms was used to derive the equations of motion.
The opened kinematic loop approach, compared with the conventional Newton-Euler method, results in fewer equations and constraint forces, while the mass matrix entries remain meaningful, but there is a stronger coupling between the equations. For the closed kinematic loop approach with equivalent complex masses, no explicit loop constraint forces are introduced in the equations. The mass of the coupler link could be modelled onto the other links of the four-bar linkage by using real and virtual equivalent masses, defining the principal points. The complex masses give a correct representation of the linear momentum, but they do not give a full dynamic equivalence if the centre of mass of the link is not on the line connecting the two joint positions. This makes this approach useful in problems mainly involving linear momentum, such as shaking force analysis.
With the method of the equivalent mass matrix, it was shown how a constant mass matrix can be used to describe the dynamics of binary links with an arbitrary mass distribution. This seems to lead to the simplest form of the equations of motion, but has as a disadvantage that the link angles appear as derived quantities and are no longer directly present in the equations of motion.
The constant mass matrix, which makes use of the joint coordinates only, can describe the mass properties of a truss element in a finite element formulation, which is fully dynamically equivalent to that of a rigid link. In addition, a body which is allowed to undergo a uniform dilatation, as is approximately the case in some auxetic metamaterials, can be described by this element. As an example, the general 4R four-bar linkage was modelled with only three truss elements instead of three or more beam elements, which is a significant reduction in the complexity of the model.
Although all three presented methods can be used to obtain the equations of motion and the reaction forces of the four-bar linkage, the equations remain complicated for all methods if an explicit form is aimed at. For different purposes, different methods to derive the equation can be useful. In particular, the methods which use principal points and principal vectors have advantages if the main interest is in the linear momentum and the resultant reaction forces, which can be used to obtain conditions for dynamic force balance.
The described methods are a step towards obtaining insight into the dynamic equations such that they can direct the synthesis process towards desired dynamic conditions.