Alessandro Cammarata

*University of Catania, Department of Mechanical and Industrial Engineering Italy*

#### **1. Introduction**

84 Serial and Parallel Robot Manipulators – Kinematics, Dynamics, Control and Optimization

Funda, J. (1988). *A Computational Analysis of Line-Oriented Screw Transformations in Robotics*.

Gu, Y.L. & Luh, J. Y. S. (1987). Dual-Number Transformation and Its Applications To

Keler, M. L. (2000). *On the theory of screws and the dual method.* Proceedings of A Symposium

Kiat Teu, K.; Kim, W.; Fuss, F. K. & Tan, J.(2006). The analysis of golf swing as a kinematic

Moon, Y-M. & Kota, S.(2002). Automated synthesis of mechanisms using dual-vector

Page, A.; Mata, V.; Hoyos, J. V. & Porcar, Rosa. (2007) Experimental depermination of

Pennestrí, E. &Stefanelli, R. (2007). Linear Algebra and numerical algorithms using dual numbers. *Journal of Multi-body System Dynamics*, Vol. 18, 3, pp. (323-344) Pennock, G. R. & Mattson, K. G. (1996). Forward position problem of two PUMA-type

Pennock, G. R. & Meehan P. J. (2000). Geometric insight into dynamics of a rigid body using

Sai-Kai, C. (2000). Simbolic computation of Jacobian of manipulators using dual number

Seilig, J.M. (1999). *Geometrical Methods in Robotics*. (1st Edition), Springer, ISBN: 0387947280,

Wang, J.; Liang, H-Z. & Sun Z. (2010). Relative Coupled Dynamics and Control using Dual

Yang, J. & Wang, X. (2010). The application of the dual number methods to scara kinematics. Ying, N.; Kim, W.; Wong, Y. & Kan, H. K. (2004). Analysis of passive motion characteristics

Number. *Systems and Control in Aeronautics and Astronautics (ISSCAA), 2010 3rd* 

of the ankle joint complex using dual Euler angle parameters. Clinical Biomechanics, Vol. 19, 2, (February 2004), pp. (153-160), doi: 10.1016/j.clinbiomech.

instantaneous screw axis in human motions. Error analysis

robots manipulating a planar four-bar linkage payload.

Commemorating the Legacy, Words and Life of Sir Robert Stawell Ball Upon the 100th Anniversary of a Treatise on the Theory of Screws, University of Cambridge,

algebra. *Mechanisms and Machine Theory,* (February 2002), pp. (143-166), doi:

Technical Report, University of Pennsylvania, U. S. A.

chain using dual Euler angle algorithm.

10.1016/S0094-114X(01)00073-8

Robotics.

Trinnity College.

the theory of screws.

New Jersey, U. S. A.

*International Symposium on*

transformations.

2003.10.005

Accurate models to describe the elasticity of robots are becoming essential for those applications involving high accelerations or high precision to improve quality in positioning and tracking of trajectories. Stiffness analysis not only involves the mechanical structure of a robot but even the control system necessary to drive actuators. Strategies aimed to reduce noise and dangerous bouncing effects could be implemented to make control systems more robust to flexibility disturbances, foreseing mechanical interaction with the control system because of regenerative and modal chattering (1). The most used approaches to study elasticity in the literature encompass: the *Finite Element Analysis* (FEA), the *Matrix Structural Analysis* (MSA), the *Virtual Joint Method* (VJM), the *Floating Frame of Reference Formulation* (FFRF) and the *Absolute Nodal Coordinates Formulation* (ANCF).

FEA is largely used to analyze the structural behavior of a mechanical system. The reliability and precision of the method allow to describe each part of a mechanical system with great detail (2). Applying FEA to a robotic system implies a time-consuming process of re-meshing in the pre-processing phase every time that the robot posture has changed. Optimization all over the workspace of a robot would require very long computational time, thus FEA models are often employed to verify components or subparts of a complex robotic system.

The MSA includes some simplifications to FEA using complex elements like beams, arcs, cables or superelements (3–5). This choice reduces the computational time and makes this method more efficient for optimization tasks. Some authors recurred to the superposition principle along with the virtual work principle to achieve the global stiffness model (6–8). Others considered the minimization of the potential energy of a PKM to find the global stiffness matrix (9), while some approaches used the total potential energy augmented adding the kinematic constraints by means of the Lagrange multipliers (10; 11).

The first papers on VJM are based on pseudo-rigid body models with "virtual joints" (12–14). More recent papers include link flexibility and linear/torsional springs to take into account bending contributes (15–19). These approaches recur to the Jacobian matrix to map the stiffness of the actuators of a PKM inside its workspace; especially for PKMs with reduced mobility, it implies that the stiffness is limited to a subspace defined by the dofs of its end-effector. Pashkevich *et al*. tried to overcome this issue by introducing a multidimensional lumped-parameter model with localized 6-dof virtual springs (20).

Finally, the FFRF and ANCF are powerful and accurate formulations, based on FEM and continuum mechanics, to study any flexible mechanical system (21). The FFRF is suitable

optimization all over the workspace of a robot too cumbersome to be faced with conventional formulations. Further, considering that the elastodynamics behavior of a robot changes according to its pose, an accurate analysis should be carried out several times to capture the right response inside the robot's workspace, drastically affecting the computational cost. The focused issue becomes even more complicated when constrained optimization routines, based on indices of elastostatics/dynamics performances, are implemented to work in all the workspace. In the latter case the search for local minimum configurations needs for a simple and robust mathematical framework adapt to iterative procedures. Working with ordinary resources, in terms of CPU speed and memory, can make optimization prohibitive unless the complex flexible multibody formulations would be simplified to meet requirements. Therefore, the second question might be: *Why simplify a complex formulation and do not create a simpler one instead?* The sought formulation is what the author is going to explain in this

On the Stiffness Analysis and Elastodynamics of Parallel Kinematic Machines 87

Let us start considering a generic PKM. It is essentially a complex robotic system with parallel kinematics in which one or more limbs connect a base platform to a moving platform, as shown in Fig. 1. The latter contains a tool, often referred to as *end-effector*, necessary to perform a certain task or, sometimes, as in the case of flight simulators or assembly stations, the end effector is the moving platform itself. The limbs connecting the two platforms are composed of links constrained by joints. A limb can be a serial kinematic chain or an articulated linkage

In performing our analysis we have to choose what is flexible and what is rigid as well. Generally, each body is flexible and the notion of rigid body is an abstraction that becomes a good approximation if strains of the structure are small when compared to displacements. Thus, considering a link flexible or rigid depends on many aspects as: material, geometry, wrenches involved in the process. Besides, some tasks might need for high precision to be accomplished, then an accurate analysis should take care of deformations to fulfill the requirements without gross errors. Here, we model the MP and BP as rigid bodies because, for the most part of industrial PKMs, these are usually one order of magnitude stiffer than the remaining links. The latter will be either flexible or rigid depending on the assumptions made by the designer. As already pointed out in the Introduction, we use 3D Euler beams to represent links, even if the treatment can be extended to superelements, as reported in (22). The formulation is linear and only small deformations are considered. Given a starting pose of the MP, the IKP allows us to find the robot's configuration; then, the elastodynamics—statics—is performed around this undeformed configuration. It might be useful to lock the actuated joints in order to avoid rigid movements and to isolate only the

Below, we summarize all steps necessary to find the elastodynamics equations of a PKM:

1. Solve the inverse kinematic problem (IKP) in order to know the starting pose, i.e. the undeformed configuration, of the PKM. Here, we make the assumption that the configuration is *frozen*, meaning that no coupling of rigid body and elastic motions is

2. Distinguish rigid and flexible links, and discretize the latter into the desired number of

3. Introduce nodes and then nodal-arrays for each node. Each flexible part has two end-nodes

flexible elements. The base and moving platforms are modeled as rigid.

4. Introduce joint-matrices and arrays for each joint.

at its end-sections; the rigid links have a single node at their center of mass.

with one or more closed kinematic chains or loops.

chapter.

flexible modes.

considered.

Fig. 1. Schematic drawing of a PKM: P prismatic joint, R revolute joint, S spherical joint, U universal joint

for large rotations and small deformations, while the ANCF is preferred for large rotations and large deformations. Notwithstanding, both the two formulations usually require great computational efforts to be extended to complex robots or to optimization techniques involving the whole robot's workspace.

In this chapter we propose a formulation to study the elastostatics and elastodynamics of PKMs. The method is linear and tries to combine some feature existing in the literature to build a solid framework, the outlines of which are described in Section 2. We start from a MSA approach based on the minimization of the strain energy in which all joints are introduced by means of constraints between nodal displacements. In this way we avoid Lagrange multipliers and the introduction of joints becomes straightforward. Unlike VJM methods based on lumped stiffness, we use 3D Euler beams to simulate links and to distribute stiffness to the flexible structure of a robot, as described in Section 3. The same set of nodal displacement arrays is used to obtain the generalized stiffness and inertia matrix, the latter being lumped or distributed, as discussed in Sections 4 and 5. Besides, the flexibility of the proposed method allows us to provide useful extension to compliant PKMs with joint flexures. The ease in setting up, the direct control of joints and the speed of execution make the procedure adapt for optimization routines, as described in Section 6. Finally, some feasible applications are described in Section 7 studying an articulated PKMs with four dofs.

### **2. Outlines of the algorithm**

Before introducing the reader to the proposed methodology, we have to clarify the reasons of this work. The first question that might arise is: *Why use a new method without recurring to well confirmed formulations existing in the literature?* The right answer is essentially tied to the simplicity of the proposed formulation. The method is addressed to designers that want to implement software to study the elastodynamics of PKMs, as well as to students and researchers that wish to create their own customized algorithm. The author experienced that the elastodynamics study of a complex robotic system is not a trivial issue, mainly when dealing with some features concerning the optimization of performances in terms of elasticity or admissible range of eigenfrequencies. Performing these analyzes often needs for a global 2 Will-be-set-by-IN-TECH

*MP*

<sup>P</sup> <sup>R</sup>

Fig. 1. Schematic drawing of a PKM: P prismatic joint, R revolute joint, S spherical joint, U

for large rotations and small deformations, while the ANCF is preferred for large rotations and large deformations. Notwithstanding, both the two formulations usually require great computational efforts to be extended to complex robots or to optimization techniques

In this chapter we propose a formulation to study the elastostatics and elastodynamics of PKMs. The method is linear and tries to combine some feature existing in the literature to build a solid framework, the outlines of which are described in Section 2. We start from a MSA approach based on the minimization of the strain energy in which all joints are introduced by means of constraints between nodal displacements. In this way we avoid Lagrange multipliers and the introduction of joints becomes straightforward. Unlike VJM methods based on lumped stiffness, we use 3D Euler beams to simulate links and to distribute stiffness to the flexible structure of a robot, as described in Section 3. The same set of nodal displacement arrays is used to obtain the generalized stiffness and inertia matrix, the latter being lumped or distributed, as discussed in Sections 4 and 5. Besides, the flexibility of the proposed method allows us to provide useful extension to compliant PKMs with joint flexures. The ease in setting up, the direct control of joints and the speed of execution make the procedure adapt for optimization routines, as described in Section 6. Finally, some feasible

applications are described in Section 7 studying an articulated PKMs with four dofs.

Before introducing the reader to the proposed methodology, we have to clarify the reasons of this work. The first question that might arise is: *Why use a new method without recurring to well confirmed formulations existing in the literature?* The right answer is essentially tied to the simplicity of the proposed formulation. The method is addressed to designers that want to implement software to study the elastodynamics of PKMs, as well as to students and researchers that wish to create their own customized algorithm. The author experienced that the elastodynamics study of a complex robotic system is not a trivial issue, mainly when dealing with some features concerning the optimization of performances in terms of elasticity or admissible range of eigenfrequencies. Performing these analyzes often needs for a global

U

universal joint

involving the whole robot's workspace.

**2. Outlines of the algorithm**

*BP*

S

optimization all over the workspace of a robot too cumbersome to be faced with conventional formulations. Further, considering that the elastodynamics behavior of a robot changes according to its pose, an accurate analysis should be carried out several times to capture the right response inside the robot's workspace, drastically affecting the computational cost. The focused issue becomes even more complicated when constrained optimization routines, based on indices of elastostatics/dynamics performances, are implemented to work in all the workspace. In the latter case the search for local minimum configurations needs for a simple and robust mathematical framework adapt to iterative procedures. Working with ordinary resources, in terms of CPU speed and memory, can make optimization prohibitive unless the complex flexible multibody formulations would be simplified to meet requirements. Therefore, the second question might be: *Why simplify a complex formulation and do not create a simpler one instead?* The sought formulation is what the author is going to explain in this chapter.

Let us start considering a generic PKM. It is essentially a complex robotic system with parallel kinematics in which one or more limbs connect a base platform to a moving platform, as shown in Fig. 1. The latter contains a tool, often referred to as *end-effector*, necessary to perform a certain task or, sometimes, as in the case of flight simulators or assembly stations, the end effector is the moving platform itself. The limbs connecting the two platforms are composed of links constrained by joints. A limb can be a serial kinematic chain or an articulated linkage with one or more closed kinematic chains or loops.

In performing our analysis we have to choose what is flexible and what is rigid as well. Generally, each body is flexible and the notion of rigid body is an abstraction that becomes a good approximation if strains of the structure are small when compared to displacements. Thus, considering a link flexible or rigid depends on many aspects as: material, geometry, wrenches involved in the process. Besides, some tasks might need for high precision to be accomplished, then an accurate analysis should take care of deformations to fulfill the requirements without gross errors. Here, we model the MP and BP as rigid bodies because, for the most part of industrial PKMs, these are usually one order of magnitude stiffer than the remaining links. The latter will be either flexible or rigid depending on the assumptions made by the designer. As already pointed out in the Introduction, we use 3D Euler beams to represent links, even if the treatment can be extended to superelements, as reported in (22). The formulation is linear and only small deformations are considered. Given a starting pose of the MP, the IKP allows us to find the robot's configuration; then, the elastodynamics—statics—is performed around this undeformed configuration. It might be useful to lock the actuated joints in order to avoid rigid movements and to isolate only the flexible modes.

Below, we summarize all steps necessary to find the elastodynamics equations of a PKM:


#### **3. Mathematical background and key concepts**

In a mechanical system flexible bodies storage and exchange energy like a tank is able to storage and supply water. The energy associated to deformation is termed *strain energy* while the aptitude of a body for deformation is tied to the property of *stiffness*. For a continuum body the stiffness is distributed and the strain energy changes according to the variation of the displacement field of its points. In a discrete flexible element inner points' displacements depend on displacements of some points or *nodes*.

#### **3.1 Nodal and joint-array**

Robotic links can be well described by means of 3D-beams, that is, flexible elements with two end-nodes, as shown in Fig. 2. The expression of the strain energy depends on the kind of displacements chosen for rotations: slopes, Euler angles, quaternions, and so on, and on the entity of displacement: large or small. Here, we choose a linear formulation based on small displacements and Euler angles to describe rotations. Based on these assumptions, the strain energy is a positive-definite quadratic form in the nodal displacement coordinates of the end-nodes arrays of the beam, where a nodal displacement array **u** = *ux uy uz u<sup>ϕ</sup> u<sup>θ</sup> u<sup>ψ</sup> T* has six scalar *displacements*, three translational and three rotational. Hereafter, subscripts and superscripts will be, respectively, referred to the beam and to one of the two end-nodes of a beam.

Beams belonging to the same link are contiguous, while joints couple beams belonging to two adjoining links. To express the kinematic bond existing between the two nodes, or sections, coupled by the joint, a constraint equation is introduced in which the nodal displacements of the two nodes are tied through the means of an array of joint displacements. Figure 2 describes two beams linked by a prismatic joint of axis parallel to the unit vector **e**. The two nodal arrays **u**<sup>2</sup> <sup>1</sup> and **<sup>u</sup>**<sup>1</sup> <sup>2</sup> are tied together by means of the translational displacement *s* of the prismatic joint *P* and by the 6-dimensional joint-array **h***<sup>P</sup>* = **<sup>e</sup>***<sup>T</sup>* **<sup>0</sup>***<sup>T</sup> <sup>T</sup>* , i.e.,

$$\mathbf{u}\_1^2 = \mathbf{u}\_2^1 + \mathbf{h}^P \mathbf{s} \tag{1}$$

**u**1 1

*<sup>X</sup> <sup>Y</sup>*

Fig. 2. Notation of two 3D Euler beams coupled by a joint.

, for translations, or angular *θ<sup>j</sup>*

*Z*

displacement about the said axis.

**3.2 Strain energy and stiffness matrix**

The strain energy *Vi*(**u**<sup>1</sup>

the end-nodes arrays **u**<sup>1</sup>

either linear *s<sup>j</sup>*

Revolute joint:

Universal joint:

joints are provided:

1

**h***<sup>R</sup>* =

 **0 0 e**<sup>1</sup> **e**<sup>2</sup>

**H***<sup>U</sup>* =

used to obtain joint displacements in terms of nodal displacements.

are the angular displacements about the axes of U.

*<sup>i</sup>* , **<sup>u</sup>**<sup>2</sup>

*<sup>i</sup>* and **<sup>u</sup>**<sup>2</sup>

the strain energy assumes the following expression

containing unit vectors indicating geometric axes, the latter containing joint displacements,

On the Stiffness Analysis and Elastodynamics of Parallel Kinematic Machines 89

where **e** is the unit vector along the axis of the revolute joint R and *θ* is the angular

where **e**<sup>1</sup> and **e**<sup>2</sup> are the unit vectors along the axes of the universal joint U and *θ*<sup>1</sup> and *θ*<sup>2</sup>

Other joints, i.e. cylindrical, spherical and so on, can be created combining together elementary prismatic and/or revolute joints. The described constraint equations are used to consider joints contribute to elastodynamics in a direct way without recurring to Lagrangian multipliers to introduce joint constraints, (10). In the next section the above equation will be

deformations, here expressed via nodal displacements. For the case of a 3D Euler beam, considered in this analysis, the strain energy *Vi* of the *i*th-beam turns into a positive-definite quadratic form of the stiffness matrix **K***<sup>i</sup>* in twelve variables, i.e. the nodal displacements of

The 12 × 12 stiffness matrix **K***<sup>i</sup>* of a 3D Euler beam with circular cross-section, and for the case of homogeneous and isotropic material, depends on geometrical and stiffness parameters as:

*Vi* <sup>=</sup> <sup>1</sup> 2 **u**˜ *T*

, *θ<sup>U</sup>* =

*<sup>i</sup>* ) of a flexible body is a nonlinear scalar function of nodal

*<sup>i</sup>* . By introducing the 12-dimensional array ˜**u***<sup>i</sup>* =

**u**2 1

**e**

2

**u**1 2 P

, for rotations. Below, two more examples of

*<sup>θ</sup>*<sup>1</sup> *<sup>θ</sup>*<sup>2</sup> *<sup>T</sup>* (3)

 **u**1 *i T* **u**2 *i <sup>T</sup> <sup>T</sup>* ,

*<sup>i</sup>* **K***i***u**˜*<sup>i</sup>* (4)

**<sup>0</sup>***<sup>T</sup>* **<sup>e</sup>***<sup>T</sup> <sup>T</sup>* , *<sup>θ</sup><sup>R</sup>* = *<sup>θ</sup>* (2)

**u**2 2

We stress that eq.(1) describes a constraint among displacements, then deformations, and not nodal coordinates. In general, **<sup>H</sup>***<sup>j</sup>* is a 6 <sup>×</sup> *<sup>m</sup>*(*j*) joint-matrix where *<sup>m</sup>*(*j*) is the dimension of the joint-array *θ<sup>j</sup>* . The joint-matrix **H***<sup>j</sup>* and the joint-array *θ<sup>j</sup>* depend on the type of joint: the former

Fig. 2. Notation of two 3D Euler beams coupled by a joint.

containing unit vectors indicating geometric axes, the latter containing joint displacements, either linear *s<sup>j</sup>* , for translations, or angular *θ<sup>j</sup>* , for rotations. Below, two more examples of joints are provided:

Revolute joint:

4 Will-be-set-by-IN-TECH

5. Individuate the couples of bodies linked by a joint and distinguish among three cases:

6. Find the equations expressing dependent nodal-array in terms of independent nodal-arrays, then find the generalized stiffness and inertia matrices. The latter change

7. Introduce the global array **q** of independent nodal coordinates. This array contains all the

8. Expand all matrices expressing each generalized stiffness and inertia matrix in terms of **q** by means of the Boolean matrices **B**<sup>1</sup> and **B**2. Then, sum all contributes to find the

9. Introduce the array **f** of generalized nodal wrenches and, finally, write the elastodynamics

In a mechanical system flexible bodies storage and exchange energy like a tank is able to storage and supply water. The energy associated to deformation is termed *strain energy* while the aptitude of a body for deformation is tied to the property of *stiffness*. For a continuum body the stiffness is distributed and the strain energy changes according to the variation of the displacement field of its points. In a discrete flexible element inner points' displacements

Robotic links can be well described by means of 3D-beams, that is, flexible elements with two end-nodes, as shown in Fig. 2. The expression of the strain energy depends on the kind of displacements chosen for rotations: slopes, Euler angles, quaternions, and so on, and on the entity of displacement: large or small. Here, we choose a linear formulation based on small displacements and Euler angles to describe rotations. Based on these assumptions, the strain energy is a positive-definite quadratic form in the nodal displacement coordinates of the

has six scalar *displacements*, three translational and three rotational. Hereafter, subscripts and superscripts will be, respectively, referred to the beam and to one of the two end-nodes of a

Beams belonging to the same link are contiguous, while joints couple beams belonging to two adjoining links. To express the kinematic bond existing between the two nodes, or sections, coupled by the joint, a constraint equation is introduced in which the nodal displacements of the two nodes are tied through the means of an array of joint displacements. Figure 2 describes two beams linked by a prismatic joint of axis parallel to the unit vector **e**. The two

We stress that eq.(1) describes a constraint among displacements, then deformations, and not nodal coordinates. In general, **<sup>H</sup>***<sup>j</sup>* is a 6 <sup>×</sup> *<sup>m</sup>*(*j*) joint-matrix where *<sup>m</sup>*(*j*) is the dimension of the

. The joint-matrix **H***<sup>j</sup>* and the joint-array *θ<sup>j</sup>* depend on the type of joint: the former

**u**2 <sup>1</sup> <sup>=</sup> **<sup>u</sup>**<sup>1</sup>

<sup>2</sup> are tied together by means of the translational displacement *s* of the

**<sup>e</sup>***<sup>T</sup>* **<sup>0</sup>***<sup>T</sup> <sup>T</sup>*

, i.e.,

<sup>2</sup> <sup>+</sup> **<sup>h</sup>***Ps* (1)

*ux uy uz u<sup>ϕ</sup> u<sup>θ</sup> u<sup>ψ</sup>*

*T*

end-nodes arrays of the beam, where a nodal displacement array **u** =

prismatic joint *P* and by the 6-dimensional joint-array **h***<sup>P</sup>* =

rigid-flexible, flexible-flexible and flexible-rigid.

whether lumped or distributed method is used.

**3. Mathematical background and key concepts**

depend on displacements of some points or *nodes*.

equations.

**3.1 Nodal and joint-array**

beam.

nodal arrays **u**<sup>2</sup>

joint-array *θ<sup>j</sup>*

<sup>1</sup> and **<sup>u</sup>**<sup>1</sup>

independent nodal arrays in the order defined by the reader.

generalized stiffness matrix **K***PKM* and inertia matrix **M***PKM* of the PKM.

$$\mathbf{h}^{\boldsymbol{R}} = \begin{bmatrix} \mathbf{0}^{\boldsymbol{T}} \ \mathbf{e}^{\boldsymbol{T}} \end{bmatrix}^{\boldsymbol{T}}, \quad \boldsymbol{\theta}^{\boldsymbol{R}} = \boldsymbol{\theta} \tag{2}$$

where **e** is the unit vector along the axis of the revolute joint R and *θ* is the angular displacement about the said axis.

Universal joint:

$$\mathbf{H}^{II} = \begin{bmatrix} \mathbf{0} & \mathbf{0} \\ \mathbf{e}^1 \ \mathbf{e}^2 \end{bmatrix}, \quad \mathbf{\theta}^{II} = \begin{bmatrix} \theta^1 \ \theta^2 \end{bmatrix}^T \tag{3}$$

where **e**<sup>1</sup> and **e**<sup>2</sup> are the unit vectors along the axes of the universal joint U and *θ*<sup>1</sup> and *θ*<sup>2</sup> are the angular displacements about the axes of U.

Other joints, i.e. cylindrical, spherical and so on, can be created combining together elementary prismatic and/or revolute joints. The described constraint equations are used to consider joints contribute to elastodynamics in a direct way without recurring to Lagrangian multipliers to introduce joint constraints, (10). In the next section the above equation will be used to obtain joint displacements in terms of nodal displacements.

#### **3.2 Strain energy and stiffness matrix**

The strain energy *Vi*(**u**<sup>1</sup> *<sup>i</sup>* , **<sup>u</sup>**<sup>2</sup> *<sup>i</sup>* ) of a flexible body is a nonlinear scalar function of nodal deformations, here expressed via nodal displacements. For the case of a 3D Euler beam, considered in this analysis, the strain energy *Vi* of the *i*th-beam turns into a positive-definite quadratic form of the stiffness matrix **K***<sup>i</sup>* in twelve variables, i.e. the nodal displacements of the end-nodes arrays **u**<sup>1</sup> *<sup>i</sup>* and **<sup>u</sup>**<sup>2</sup> *<sup>i</sup>* . By introducing the 12-dimensional array ˜**u***<sup>i</sup>* = **u**1 *i T* **u**2 *i <sup>T</sup> <sup>T</sup>* , the strain energy assumes the following expression

$$V\_i = \frac{1}{2} \tilde{\mathbf{u}}\_i^T \mathbf{K}\_i \tilde{\mathbf{u}}\_i \tag{4}$$

The 12 × 12 stiffness matrix **K***<sup>i</sup>* of a 3D Euler beam with circular cross-section, and for the case of homogeneous and isotropic material, depends on geometrical and stiffness parameters as: cross section area *A*, length of the beam *L*, torsional constant *J*, mass moment of inertia *I*, the Young module *E* and shear module *G*. For our purposes, it is convenient to divide **K***<sup>i</sup>* into blocks, i.e.

$$\mathbf{K}\_{i} = \begin{bmatrix} \mathbf{K}\_{i}^{1,1} \mathbf{ K}\_{i}^{1,2} \\ \mathbf{K}\_{i}^{2,1} \mathbf{ K}\_{i}^{2,2} \end{bmatrix} \tag{5}$$

**u**1 2

*<sup>θ</sup>* **<sup>d</sup>**

1

F

<sup>2</sup>,..., **<sup>u</sup>***G*1, **<sup>u</sup>***G*2,..., *<sup>θ</sup>*1, *<sup>θ</sup>*2,...) of this system depends on nodal

<sup>1</sup> **<sup>u</sup>**<sup>2</sup>

On the Stiffness Analysis and Elastodynamics of Parallel Kinematic Machines 91

Let us consider a linkage of flexible beams and rigid bodies connected by means of joints.

displacement arrays of both rigid and flexible parts and on joint displacement arrays. In the following, we will show how to express *V* in terms of a reduced set of independent nodal coordinates. In this process, we will start from three elementary blocks composed of rigid and flexible bodies that will be combined to build a generic linkage, for instance, the limb of a PKM. The first and the third case pertain the coupling between a rigid body and a flexible beam by means of a joint; the second case describes the coupling of two beams belonging to two different links. The choice to use two cases to describe the rigid body-flexible body connection is only due to convenience to follow the order of bodies from the base to the moving platform of a PKM. The reader might recur to a unique case to simplify the treatment. The last part of the section is devoted to some insight on joints' stiffness and

Figure 3 describes a rigid body R coupled to a flexible beam F by means of a joint. The strain

**K**2,1 <sup>2</sup> **<sup>K</sup>**2,2 2

where **G** depends on **d**. By substituting the previous equation into eq.(8) we find that *Va* =

*dVa*/*dθ* = **0***<sup>T</sup>*

where **0**<sup>6</sup> is the six-dimensional zero array. After rearrangements and simplifications, we find

*θ* = **Y**1**u**<sup>1</sup>

<sup>2</sup>, *θ*). If the joint is passive, its displacement *θ* depends on the elastic properties of the

<sup>1</sup> <sup>+</sup> **<sup>Y</sup>**2**u**<sup>2</sup>

<sup>1</sup> and **<sup>u</sup>**<sup>2</sup>

<sup>1</sup>, **<sup>u</sup>**<sup>2</sup>

 **u**1 2 **u**2 2 

<sup>1</sup> and *<sup>θ</sup>* recalling eq.(1), while **<sup>u</sup>**<sup>2</sup>

*<sup>T</sup>* **K**1,1 <sup>2</sup> **<sup>K</sup>**1,2 2

**u**1 <sup>2</sup> <sup>=</sup> **Gu**<sup>1</sup>

energy *Va* of the beam is function of the nodal displacement arrays **u**<sup>1</sup>

*Va* <sup>=</sup> <sup>1</sup> 2 **u**1 2 **u**2 2

<sup>1</sup> from eq.(7): upon combining both the equations, we obtain

<sup>2</sup> can be expressed in terms of **<sup>u</sup>**<sup>2</sup>

system and, therefore, on the two displacements **u**<sup>1</sup>

function of the two mentioned array, i.e. *Va* = *V*(**u**<sup>1</sup>

minimize the strain energy *Va* w.r.t. *θ*:

R

**u**1

1

Fig. 3. Case a): Rigid body-flexible beam.

<sup>1</sup>, **<sup>u</sup>**<sup>2</sup> <sup>1</sup>, **<sup>u</sup>**<sup>1</sup>

**4. Stiffness matrix determination**

The strain energy *V*(**u**<sup>1</sup>

feasible application to flexures.

The array **u**<sup>1</sup>

**u**1

*V*(**u**<sup>1</sup> <sup>1</sup>, **<sup>u</sup>**<sup>2</sup>

**4.1 Case a) Rigid body-flexible beam**

**u**2 2

<sup>2</sup> and **<sup>u</sup>**<sup>2</sup>

<sup>1</sup> + **H***θ* (9)

<sup>2</sup>: thus, it implies that *Va* is only

<sup>2</sup>). In order to obtain the law for *θ* we

<sup>6</sup> (10)

<sup>2</sup> (11)

<sup>2</sup>, i.e.

<sup>1</sup>, in turn, is tied to

(8)

2

**e**

in which **K**2,1 *<sup>i</sup>* = (**K**1,2 *<sup>i</sup>* )*<sup>T</sup>* and the other blocks are defined as

**K**1,1 *<sup>i</sup>* <sup>=</sup> *<sup>E</sup> L* ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ *A* 0 0 0 00 0 <sup>12</sup>*<sup>I</sup> <sup>L</sup>*<sup>2</sup> <sup>000</sup> <sup>6</sup>*<sup>I</sup> L* 0 0 <sup>12</sup>*<sup>I</sup> <sup>L</sup>*<sup>2</sup> <sup>0</sup> <sup>−</sup>6*<sup>I</sup> <sup>L</sup>* 0 00 0 *G J <sup>E</sup>* 0 0 0 0 <sup>−</sup>6*<sup>I</sup> <sup>L</sup>* 0 4*I* 0 0 <sup>6</sup>*<sup>I</sup> <sup>L</sup>* 0 0 04*I* ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ (6a) **K**1,2 *<sup>i</sup>* <sup>=</sup> *<sup>E</sup> L* ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ −*A* 0 0 0 00 <sup>0</sup> <sup>−</sup>12*<sup>I</sup> <sup>L</sup>*<sup>2</sup> 0 00 <sup>6</sup>*<sup>I</sup> L* 0 0 <sup>−</sup>12*<sup>I</sup> <sup>L</sup>*<sup>2</sup> <sup>0</sup> <sup>−</sup>6*<sup>I</sup> <sup>L</sup>* 0 00 0 <sup>−</sup> *G J <sup>E</sup>* 0 0 0 0 <sup>6</sup>*<sup>I</sup> <sup>L</sup>* 0 4*I* 0 <sup>0</sup> <sup>−</sup>6*<sup>I</sup> <sup>L</sup>* 0 0 04*I* ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ (6b) **K**2,2 *<sup>i</sup>* <sup>=</sup> *<sup>E</sup> L* ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ *A* 0 0 00 0 0 <sup>12</sup>*<sup>I</sup> <sup>L</sup>*<sup>2</sup> 0 00 <sup>−</sup>6*<sup>I</sup> L* 0 0 <sup>12</sup>*<sup>I</sup> <sup>L</sup>*<sup>2</sup> <sup>0</sup> <sup>6</sup>*<sup>I</sup> <sup>L</sup>* 0 00 0 *G J <sup>E</sup>* 0 0 0 0 <sup>6</sup>*<sup>I</sup> <sup>L</sup>* 0 4*I* 0 <sup>0</sup> <sup>−</sup>6*<sup>I</sup> <sup>L</sup>* 0 004*I* ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ (6c)

where the two diagonal block matrices, respectively, refer to the nodal displacements of the two end-nodes while the extra-diagonal blocks refer to the coupling among nodal displacements of different nodes: in fact, each entry of the generic 6 × 6 block-matrix **K***l*,*<sup>m</sup> <sup>i</sup>* can be thought as a force— torque—at the *l*th-node of the *i*th-beam when a unit displacement—rotation—is applied to the *m*th-node.

#### **3.3 Rigid body displacement**

Let us consider a rigid body with center of mass at the point *G* and a generic point *P* inside its volume, besides let **d***<sup>P</sup>* be the vector pointing from *G* to the point *P*. If the body can accomplish only small displacements/rotations, let **p** be the displacement of point *G* and **r**, be the axial vector of the small rotation matrix **R**, (9). Upon these assumptions, the displacement array **u***<sup>G</sup>* of *G* is defined as **u***<sup>G</sup>* = � **p***<sup>T</sup>* **r***<sup>T</sup>* �*<sup>T</sup>* .

the displacement array **u***<sup>P</sup>* can be expressed in terms of **u***<sup>G</sup>* by means of the following equation:

$$\mathbf{u}\_P = \mathbf{G}\_P \mathbf{u}\_G \quad \mathbf{G}\_P = \begin{bmatrix} \mathbf{1} & -\mathbf{D}\_P \\ \mathbf{O} & \mathbf{1} \end{bmatrix} \tag{7}$$

where **1** and **O**, respectively, are the 3 × 3 identity- and zero-matrices and **D***<sup>P</sup>* is the Cross-Product Matrix (C.P.M.) of **d***P*, (23).

Fig. 3. Case a): Rigid body-flexible beam.

#### **4. Stiffness matrix determination**

6 Will-be-set-by-IN-TECH

cross section area *A*, length of the beam *L*, torsional constant *J*, mass moment of inertia *I*, the Young module *E* and shear module *G*. For our purposes, it is convenient to divide **K***<sup>i</sup>* into

> **K**2,1 *<sup>i</sup>* **<sup>K</sup>**2,2 *i*

*A* 0 0 0 00

*<sup>L</sup>*<sup>2</sup> <sup>000</sup> <sup>6</sup>*<sup>I</sup>*

*<sup>L</sup>*<sup>2</sup> <sup>0</sup> <sup>−</sup>6*<sup>I</sup>*

*<sup>L</sup>* 0 0 04*I*

−*A* 0 0 0 00

�

*L*

⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦

*L*

⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦

*<sup>L</sup>* 0

*<sup>E</sup>* 0 0

⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦

*<sup>L</sup>* 0

*<sup>E</sup>* 0 0

*<sup>L</sup>*<sup>2</sup> 0 00 <sup>6</sup>*<sup>I</sup>*

*<sup>L</sup>* 0 0 04*I*

*<sup>L</sup>* 0

*<sup>E</sup>* 0 0

*<sup>L</sup>* 0 4*I* 0

*<sup>L</sup>* 0 004*I*

where the two diagonal block matrices, respectively, refer to the nodal displacements of the two end-nodes while the extra-diagonal blocks refer to the coupling among nodal displacements of different nodes: in fact, each entry of the generic 6 × 6 block-matrix

*<sup>i</sup>* can be thought as a force— torque—at the *l*th-node of the *i*th-beam when a unit

Let us consider a rigid body with center of mass at the point *G* and a generic point *P* inside its volume, besides let **d***<sup>P</sup>* be the vector pointing from *G* to the point *P*. If the body can accomplish only small displacements/rotations, let **p** be the displacement of point *G* and **r**, be the axial vector of the small rotation matrix **R**, (9). Upon these assumptions, the displacement array **u***<sup>G</sup>*

the displacement array **u***<sup>P</sup>* can be expressed in terms of **u***<sup>G</sup>* by means of the following equation:

where **1** and **O**, respectively, are the 3 × 3 identity- and zero-matrices and **D***<sup>P</sup>* is the

� **<sup>1</sup>** <sup>−</sup>**D***<sup>P</sup>* **O 1**

�

**u***<sup>P</sup>* = **G***P***u***G*, **G***<sup>P</sup>* =

*<sup>L</sup>*<sup>2</sup> <sup>0</sup> <sup>−</sup>6*<sup>I</sup>*

*<sup>L</sup>* 0 4*I* 0

*L*

*<sup>L</sup>* 0 4*I* 0

(5)

(6a)

(6b)

(6c)

(7)

� **K**1,1 *<sup>i</sup>* **<sup>K</sup>**1,2 *i*

**K***<sup>i</sup>* =

*<sup>i</sup>* )*<sup>T</sup>* and the other blocks are defined as

0 <sup>12</sup>*<sup>I</sup>*

0 <sup>6</sup>*<sup>I</sup>*

0 0 <sup>12</sup>*<sup>I</sup>*

0 0 <sup>−</sup>6*<sup>I</sup>*

<sup>0</sup> <sup>−</sup>12*<sup>I</sup>*

<sup>0</sup> <sup>−</sup>6*<sup>I</sup>*

0 0 <sup>12</sup>*<sup>I</sup>*

0 0 <sup>6</sup>*<sup>I</sup>*

<sup>0</sup> <sup>−</sup>6*<sup>I</sup>*

00 0 *G J*

0 <sup>12</sup>*<sup>I</sup>*

0 0 <sup>−</sup>12*<sup>I</sup>*

0 0 <sup>6</sup>*<sup>I</sup>*

00 0 <sup>−</sup> *G J*

*A* 0 0 00 0

*<sup>L</sup>*<sup>2</sup> 0 00 <sup>−</sup>6*<sup>I</sup>*

*<sup>L</sup>*<sup>2</sup> <sup>0</sup> <sup>6</sup>*<sup>I</sup>*

00 0 *G J*

⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣

⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣

⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣

blocks, i.e.

in which **K**2,1

**K***l*,*<sup>m</sup>*

*<sup>i</sup>* = (**K**1,2

**K**1,1 *<sup>i</sup>* <sup>=</sup> *<sup>E</sup> L*

**K**1,2 *<sup>i</sup>* <sup>=</sup> *<sup>E</sup> L*

**K**2,2 *<sup>i</sup>* <sup>=</sup> *<sup>E</sup> L*

displacement—rotation—is applied to the *m*th-node.

**p***<sup>T</sup>* **r***<sup>T</sup>* �*<sup>T</sup>* .

**3.3 Rigid body displacement**

of *G* is defined as **u***<sup>G</sup>* = �

Cross-Product Matrix (C.P.M.) of **d***P*, (23).

Let us consider a linkage of flexible beams and rigid bodies connected by means of joints. The strain energy *V*(**u**<sup>1</sup> <sup>1</sup>, **<sup>u</sup>**<sup>2</sup> <sup>1</sup>, **<sup>u</sup>**<sup>1</sup> <sup>2</sup>,..., **<sup>u</sup>***G*1, **<sup>u</sup>***G*2,..., *<sup>θ</sup>*1, *<sup>θ</sup>*2,...) of this system depends on nodal displacement arrays of both rigid and flexible parts and on joint displacement arrays. In the following, we will show how to express *V* in terms of a reduced set of independent nodal coordinates. In this process, we will start from three elementary blocks composed of rigid and flexible bodies that will be combined to build a generic linkage, for instance, the limb of a PKM. The first and the third case pertain the coupling between a rigid body and a flexible beam by means of a joint; the second case describes the coupling of two beams belonging to two different links. The choice to use two cases to describe the rigid body-flexible body connection is only due to convenience to follow the order of bodies from the base to the moving platform of a PKM. The reader might recur to a unique case to simplify the treatment. The last part of the section is devoted to some insight on joints' stiffness and feasible application to flexures.

#### **4.1 Case a) Rigid body-flexible beam**

Figure 3 describes a rigid body R coupled to a flexible beam F by means of a joint. The strain energy *Va* of the beam is function of the nodal displacement arrays **u**<sup>1</sup> <sup>2</sup> and **<sup>u</sup>**<sup>2</sup> <sup>2</sup>, i.e.

$$V\_a = \frac{1}{2} \begin{bmatrix} \mathbf{u}\_2^1 \\ \mathbf{u}\_2^2 \end{bmatrix}^T \begin{bmatrix} \mathbf{K}\_2^{1,1} \ \mathbf{K}\_2^{1,2} \\ \mathbf{K}\_2^{2,1} \ \mathbf{K}\_2^{2,2} \end{bmatrix} \begin{bmatrix} \mathbf{u}\_2^1 \\ \mathbf{u}\_2^2 \end{bmatrix} \tag{8}$$

The array **u**<sup>1</sup> can be expressed in terms of **<sup>u</sup>**<sup>2</sup> and *<sup>θ</sup>* recalling eq.(1), while **<sup>u</sup>**<sup>2</sup> , in turn, is tied to **u**1 from eq.(7): upon combining both the equations, we obtain

$$\mathbf{u}\_2^1 = \mathbf{G}\mathbf{u}\_1^1 + \mathbf{H}\boldsymbol{\theta} \tag{9}$$

where **G** depends on **d**. By substituting the previous equation into eq.(8) we find that *Va* = *V*(**u**<sup>1</sup> <sup>1</sup>, **<sup>u</sup>**<sup>2</sup> <sup>2</sup>, *θ*). If the joint is passive, its displacement *θ* depends on the elastic properties of the system and, therefore, on the two displacements **u**<sup>1</sup> <sup>1</sup> and **<sup>u</sup>**<sup>2</sup> <sup>2</sup>: thus, it implies that *Va* is only function of the two mentioned array, i.e. *Va* = *V*(**u**<sup>1</sup> <sup>1</sup>, **<sup>u</sup>**<sup>2</sup> <sup>2</sup>). In order to obtain the law for *θ* we minimize the strain energy *Va* w.r.t. *θ*:

$$dV\_a / d\mathbf{\theta} = \mathbf{0}\_6^T \tag{10}$$

where **0**<sup>6</sup> is the six-dimensional zero array. After rearrangements and simplifications, we find

$$\boldsymbol{\theta} = \mathbf{Y}^1 \mathbf{u}\_1^1 + \mathbf{Y}^2 \mathbf{u}\_2^2 \tag{11}$$

Fig. 4. Case b): Flexible beam-flexible beam.

where

$$\mathbf{Y}^{1} = -(\mathbf{H}^{T}\mathbf{K}\_{2}^{1,1}\mathbf{H})^{-1}\mathbf{H}^{T}\mathbf{K}\_{2}^{1,1}\mathbf{G} \tag{12a}$$

Now, by imposing that the derivative of *Vb* w.r.t. *θ* vanishes, we obtain

Following the previous condition, even the dependent nodal-array **u**<sup>1</sup>

*∂***u**<sup>2</sup> 1 *∂***u**<sup>1</sup> 2 + *∂Vb ∂***u**<sup>1</sup> 2

where the 6 <sup>×</sup> 6 matrices **<sup>G</sup>**<sup>1</sup> and **<sup>G</sup>**<sup>2</sup> have the following expressions:

<sup>1</sup> <sup>+</sup> **<sup>K</sup>**1,1

<sup>1</sup> <sup>+</sup> **<sup>K</sup>**2,2

<sup>1</sup> <sup>+</sup> **<sup>K</sup>**1,1

**u**2

<sup>1</sup> <sup>=</sup> **<sup>X</sup>**1**u**<sup>1</sup>

*dVb d***u**<sup>1</sup> 2 <sup>≡</sup> *<sup>∂</sup>Vb ∂***u**<sup>2</sup> 1

<sup>2</sup> is then written in terms of **<sup>u</sup>**<sup>1</sup>

(**K**2,1

**<sup>G</sup>**<sup>1</sup> <sup>=</sup> <sup>−</sup> (**K**2,2

**<sup>G</sup>**<sup>2</sup> <sup>=</sup> <sup>−</sup> (**K**2,2

The joint-arrays *θ<sup>j</sup>* becomes,

The dependent nodal-array **u**<sup>2</sup>

the final expression of *Vb* in terms ˜**u***<sup>b</sup>* is

 **1 O X**<sup>1</sup> **X**<sup>2</sup>

where **K***<sup>b</sup>* is the 12 × 12 stiffness matrix.

*<sup>T</sup>* **K**1,1 <sup>1</sup> **<sup>K</sup>**1,2 1

**K**2,1 <sup>1</sup> **<sup>K</sup>**2,2 1

Let ˜**u***<sup>b</sup>* =

 **u**1 1 *T* **u**2 2 *<sup>T</sup> <sup>T</sup>*

*Vb* <sup>=</sup> <sup>1</sup> 2 **u**˜ *T*

**K***<sup>b</sup>* =

where

*V*(**u**<sup>1</sup> <sup>1</sup>, **<sup>u</sup>**<sup>2</sup> 1(**u**<sup>1</sup> <sup>2</sup>, *<sup>θ</sup>*(**u**<sup>1</sup> <sup>1</sup>, **<sup>u</sup>**<sup>1</sup> <sup>2</sup>)), **<sup>u</sup>**<sup>1</sup> <sup>2</sup>, **<sup>u</sup>**<sup>2</sup>

The array **u**<sup>1</sup>

*θ* = **F**1**u**<sup>1</sup>

**<sup>F</sup>**<sup>1</sup> <sup>=</sup> <sup>−</sup>(**H***T***K**2,2

**<sup>F</sup>**<sup>2</sup> <sup>=</sup> <sup>−</sup>(**H***T***K**2,2

<sup>≡</sup> *<sup>∂</sup>Vb ∂***u**<sup>2</sup> 1

**H***T***K**2,2

**H***T***K**2,2

<sup>1</sup> <sup>+</sup> **<sup>Y</sup>**2**u**<sup>2</sup>

<sup>1</sup> <sup>+</sup> **<sup>X</sup>**2**u**<sup>2</sup>

 **1 O X**<sup>1</sup> **X**<sup>2</sup>  + **G**<sup>1</sup> **G**<sup>2</sup> **O 1** *<sup>T</sup>*

<sup>1</sup>, **<sup>u</sup>**<sup>2</sup> <sup>2</sup>, i.e.

<sup>2</sup> <sup>=</sup> **<sup>G</sup>**1**u**<sup>1</sup>

**u**1

<sup>2</sup> <sup>+</sup> **<sup>F</sup>**2*<sup>T</sup>*

<sup>2</sup> <sup>+</sup> **<sup>F</sup>**2*<sup>T</sup>*

*θ* = **Y**1**u**<sup>1</sup>

<sup>1</sup> **HF**<sup>1</sup> <sup>+</sup> **<sup>F</sup>**2*<sup>T</sup>*

<sup>1</sup> <sup>+</sup> **<sup>F</sup>**2**u**<sup>1</sup>

On the Stiffness Analysis and Elastodynamics of Parallel Kinematic Machines 93

<sup>1</sup> **<sup>H</sup>**)−1**H***T***K**2,1

<sup>1</sup> **<sup>H</sup>**)−1**H***T***K**2,2

<sup>2</sup>); therefore, by applying the chain-rule, we obtain

<sup>1</sup> <sup>+</sup> **<sup>G</sup>**2**u**<sup>2</sup>

<sup>1</sup> <sup>+</sup> **<sup>F</sup>**2*<sup>T</sup>*

<sup>1</sup> <sup>+</sup> **<sup>F</sup>**2*<sup>T</sup>*

**H***T***K**2,1

(**1**<sup>6</sup> <sup>+</sup> **HF**2) + *<sup>∂</sup>Vb*

**H***T***K**2,2

**H***T***K**2,2

**Y**<sup>1</sup> = **F**<sup>1</sup> + **F**2**G**1, **Y**<sup>2</sup> = **F**2**G**<sup>2</sup> (24)

**X**<sup>1</sup> = **G**<sup>1</sup> + **HY**1, **X**<sup>2</sup> = **G**<sup>2</sup> + **HY**<sup>2</sup> (26)

**K**1,1 <sup>2</sup> **<sup>K</sup>**1,2 2

**G**<sup>1</sup> **G**<sup>2</sup> **O 1**

(28)

**K**2,1 <sup>2</sup> **<sup>K</sup>**2,2 2

<sup>1</sup> is obtained by substituting eq.(21) and eq.(23) into eq.(17):

be the 12-dimensional array of independent nodal displacements, then

*<sup>b</sup>* **K***b***u**˜ *<sup>b</sup>* (27)

<sup>1</sup> <sup>+</sup> **<sup>F</sup>**2*<sup>T</sup>*

*∂***u**<sup>1</sup> 2

<sup>1</sup> **HF**2)−<sup>1</sup>

<sup>1</sup> **HF**2)−1**K**1,2

<sup>2</sup> (23)

<sup>2</sup> (25)

**H***T***K**2,2

= **0***<sup>T</sup>*

<sup>2</sup> (21)

<sup>1</sup> **HF**1) (22a)

<sup>1</sup> (22b)

<sup>2</sup> (18)

<sup>1</sup> (19a)

<sup>1</sup> (19b)

<sup>2</sup> must minimize *Vb* =

<sup>6</sup> (20)

$$\mathbf{Y}^2 = -(\mathbf{H}^T \mathbf{K}\_2^{1,1} \mathbf{H})^{-1} \mathbf{H}^T \mathbf{K}\_2^{1,2} \tag{12b}$$

Then, by substituting eq.(11) into eq.(9), we obtain

$$\mathbf{u}\_2^1 = \mathbf{X}^1 \mathbf{u}\_1^1 + \mathbf{X}^2 \mathbf{u}\_2^2 \tag{13}$$

$$\mathbf{X}^1 = \mathbf{G} + \mathbf{H} \mathbf{Y}^1, \quad \mathbf{X}^2 = \mathbf{H} \mathbf{Y}^2 \tag{14}$$

Let us define the 12-dimensional array ˜**u***<sup>a</sup>* = **u**1 1 *T* **u**2 2 *<sup>T</sup> <sup>T</sup>* and let us substitute the above expression into *Va*, thereby obtaining

$$V\_{a} = \frac{1}{2} \tilde{\mathbf{u}}\_{a}^{T} \mathbf{K}\_{d} \tilde{\mathbf{u}}\_{d\prime} \mathbf{K}\_{d} = \begin{bmatrix} \mathbf{X}^{1} \mathbf{X}^{2} \\ \mathbf{O} \mathbf{ 1} \end{bmatrix}^{T} \begin{bmatrix} \mathbf{K}\_{2}^{1,1} \mathbf{K}\_{2}^{1,2} \\ \mathbf{K}\_{2}^{2,1} \mathbf{K}\_{2}^{2,2} \end{bmatrix} \begin{bmatrix} \mathbf{X}^{1} \mathbf{X}^{2} \\ \mathbf{O} \mathbf{ 1} \end{bmatrix} \tag{15}$$

where **K***<sup>a</sup>* is the 12 × 12 stiffness matrix sought.

#### **4.2 Case b) Flexible body-flexible body**

For the case b) of Fig.4 two beams are coupled by a joint. The strain energy *Vb* is function of the nodal displacement arrays of the two flexible bodies, i.e.

$$V\_b = \frac{1}{2} \begin{bmatrix} \mathbf{u}\_1^1 \\ \mathbf{u}\_1^2 \end{bmatrix}^T \begin{bmatrix} \mathbf{K}\_1^{1,1} \mathbf{K}\_1^{1,2} \\ \mathbf{K}\_1^{2,1} \mathbf{K}\_1^{2,2} \end{bmatrix} \begin{bmatrix} \mathbf{u}\_1^1 \\ \mathbf{u}\_1^2 \end{bmatrix} + \frac{1}{2} \begin{bmatrix} \mathbf{u}\_2^1 \\ \mathbf{u}\_2^2 \end{bmatrix}^T \begin{bmatrix} \mathbf{K}\_2^{1,1} \mathbf{K}\_2^{1,2} \\ \mathbf{K}\_2^{2,1} \mathbf{K}\_2^{2,2} \end{bmatrix} \begin{bmatrix} \mathbf{u}\_2^1 \\ \mathbf{u}\_2^2 \end{bmatrix} \tag{16}$$

The four arrays are not all independent as the following equation stands:

$$\mathbf{u}\_1^2 = \mathbf{u}\_2^1 + \mathbf{H}\theta \tag{17}$$

The strain energy is, thus, dependent on **u**<sup>1</sup> <sup>1</sup>, **<sup>u</sup>**<sup>1</sup> <sup>2</sup>, **<sup>u</sup>**<sup>2</sup> <sup>2</sup> and *<sup>θ</sup>*: *Vb* <sup>=</sup> *<sup>V</sup>*(**u**<sup>1</sup> <sup>1</sup>, **<sup>u</sup>**<sup>2</sup> <sup>1</sup>, **<sup>u</sup>**<sup>1</sup> <sup>2</sup>, **<sup>u</sup>**<sup>2</sup> <sup>2</sup>) ≡ *V*(**u**<sup>1</sup> <sup>1</sup>, **<sup>u</sup>**<sup>1</sup> <sup>2</sup>, **<sup>u</sup>**<sup>2</sup> <sup>2</sup>, *<sup>θ</sup>*). Here, we choose to minimize *Vb* w.r.t. *<sup>θ</sup>*, as for the case a), and **<sup>u</sup>**<sup>1</sup> <sup>2</sup>. The reader should notice that our choice is not unique, it is only a particular reduction process necessary for our purposes, it means that the reader might develop a treatment in which **u**<sup>1</sup> <sup>2</sup> is not dependent anymore.

Now, by imposing that the derivative of *Vb* w.r.t. *θ* vanishes, we obtain

$$\boldsymbol{\theta} = \mathbf{F}^1 \mathbf{u}\_1^1 + \mathbf{F}^2 \mathbf{u}\_2^1 \tag{18}$$

where

8 Will-be-set-by-IN-TECH

**u**1 1

1

Fig. 4. Case b): Flexible beam-flexible beam.

Then, by substituting eq.(11) into eq.(9), we obtain

Let us define the 12-dimensional array ˜**u***<sup>a</sup>* =

*Va* <sup>=</sup> <sup>1</sup> 2 **u**˜ *T*

where **K***<sup>a</sup>* is the 12 × 12 stiffness matrix sought.

The strain energy is, thus, dependent on **u**<sup>1</sup>

expression into *Va*, thereby obtaining

**4.2 Case b) Flexible body-flexible body**

*Vb* <sup>=</sup> <sup>1</sup> 2 **u**1 1 **u**2 1

not dependent anymore.

where

*V*(**u**<sup>1</sup> <sup>1</sup>, **<sup>u</sup>**<sup>1</sup> <sup>2</sup>, **<sup>u</sup>**<sup>2</sup> **u**1 2

F

**<sup>Y</sup>**<sup>1</sup> <sup>=</sup> <sup>−</sup>(**H***T***K**1,1

**<sup>Y</sup>**<sup>2</sup> <sup>=</sup> <sup>−</sup>(**H***T***K**1,1

<sup>2</sup> <sup>=</sup> **<sup>X</sup>**1**u**<sup>1</sup>

 **X**<sup>1</sup> **X**<sup>2</sup> **O 1**

 **u**1 1 **u**2 1 + 1 2 **u**1 2 **u**2 2

**u**2 <sup>1</sup> <sup>=</sup> **<sup>u</sup>**<sup>1</sup>

**u**1

*<sup>a</sup>* **K***a***u**˜ *<sup>a</sup>*, **K***<sup>a</sup>* =

the nodal displacement arrays of the two flexible bodies, i.e.

**K**2,1 <sup>1</sup> **<sup>K</sup>**2,2 1

The four arrays are not all independent as the following equation stands:

*<sup>T</sup>* **K**1,1 <sup>1</sup> **<sup>K</sup>**1,2 1

**u**2 1

<sup>1</sup> <sup>+</sup> **<sup>X</sup>**2**u**<sup>2</sup>

 **u**1 1 *T* **u**2 2 *<sup>T</sup> <sup>T</sup>*

For the case b) of Fig.4 two beams are coupled by a joint. The strain energy *Vb* is function of

<sup>1</sup>, **<sup>u</sup>**<sup>1</sup> <sup>2</sup>, **<sup>u</sup>**<sup>2</sup>

<sup>2</sup>, *<sup>θ</sup>*). Here, we choose to minimize *Vb* w.r.t. *<sup>θ</sup>*, as for the case a), and **<sup>u</sup>**<sup>1</sup>

reader should notice that our choice is not unique, it is only a particular reduction process necessary for our purposes, it means that the reader might develop a treatment in which **u**<sup>1</sup>

*<sup>T</sup>* **K**1,1 <sup>2</sup> **<sup>K</sup>**1,2 2

**K**2,1 <sup>2</sup> **<sup>K</sup>**2,2 2

> *<sup>T</sup>* **K**1,1 <sup>2</sup> **<sup>K</sup>**1,2 2

**K**2,1 <sup>2</sup> **<sup>K</sup>**2,2 2

F

<sup>2</sup> **<sup>H</sup>**)−1**H***T***K**1,1

<sup>2</sup> **<sup>H</sup>**)−1**H***T***K**1,2

**u**2 2

<sup>2</sup> **G** (12a)

<sup>2</sup> (12b)

and let us substitute the above

(15)

(16)

<sup>2</sup>. The

<sup>2</sup> is

<sup>2</sup> (13)

**X**<sup>1</sup> = **G** + **HY**1, **X**<sup>2</sup> = **HY**<sup>2</sup> (14)

**X**<sup>1</sup> **X**<sup>2</sup> **O 1**

 **u**1 2 **u**2 2 

> <sup>1</sup>, **<sup>u</sup>**<sup>2</sup> <sup>1</sup>, **<sup>u</sup>**<sup>1</sup> <sup>2</sup>, **<sup>u</sup>**<sup>2</sup> <sup>2</sup>) ≡

<sup>2</sup> + **H***θ* (17)

<sup>2</sup> and *<sup>θ</sup>*: *Vb* <sup>=</sup> *<sup>V</sup>*(**u**<sup>1</sup>

2

**e**

*θ*

$$\mathbf{F}^1 = -(\mathbf{H}^T \mathbf{K}\_1^{2,2} \mathbf{H})^{-1} \mathbf{H}^T \mathbf{K}\_1^{2,1} \tag{19a}$$

$$\mathbf{F}^2 = -(\mathbf{H}^T \mathbf{K}\_1^{2,2} \mathbf{H})^{-1} \mathbf{H}^T \mathbf{K}\_1^{2,2} \tag{19b}$$

Following the previous condition, even the dependent nodal-array **u**<sup>1</sup> <sup>2</sup> must minimize *Vb* = *V*(**u**<sup>1</sup> <sup>1</sup>, **<sup>u</sup>**<sup>2</sup> 1(**u**<sup>1</sup> <sup>2</sup>, *<sup>θ</sup>*(**u**<sup>1</sup> <sup>1</sup>, **<sup>u</sup>**<sup>1</sup> <sup>2</sup>)), **<sup>u</sup>**<sup>1</sup> <sup>2</sup>, **<sup>u</sup>**<sup>2</sup> <sup>2</sup>); therefore, by applying the chain-rule, we obtain

$$\frac{dV\_b}{d\mathbf{u}\_2^1} \equiv \frac{\partial V\_b}{\partial \mathbf{u}\_1^2} \frac{\partial \mathbf{u}\_1^2}{\partial \mathbf{u}\_2^1} + \frac{\partial V\_b}{\partial \mathbf{u}\_2^1} \equiv \frac{\partial V\_b}{\partial \mathbf{u}\_1^2} (\mathbf{1}\_6 + \mathbf{H} \mathbf{F}^2) + \frac{\partial V\_b}{\partial \mathbf{u}\_2^1} = \mathbf{0}\_6^T \tag{20}$$

The array **u**<sup>1</sup> <sup>2</sup> is then written in terms of **<sup>u</sup>**<sup>1</sup> <sup>1</sup>, **<sup>u</sup>**<sup>2</sup> <sup>2</sup>, i.e.

$$\mathbf{u}\_2^1 = \mathbf{G}^1 \mathbf{u}\_1^1 + \mathbf{G}^2 \mathbf{u}\_2^2 \tag{21}$$

where the 6 <sup>×</sup> 6 matrices **<sup>G</sup>**<sup>1</sup> and **<sup>G</sup>**<sup>2</sup> have the following expressions:

$$\mathbf{G}^{1} = -\begin{pmatrix} \mathbf{K}\_{1}^{2,2} + \mathbf{K}\_{2}^{1,1} + \mathbf{F}^{2}\mathbf{H}^{T}\mathbf{K}\_{1}^{2,2} + \mathbf{F}^{2}\mathbf{H}^{T}\mathbf{K}\_{1}^{2,2}\mathbf{H}\mathbf{F}^{2} \end{pmatrix}^{-1}$$

$$(\mathbf{K}\_{1}^{2,1} + \mathbf{K}\_{1}^{2,2}\mathbf{H}\mathbf{F}^{1} + \mathbf{F}^{2}\mathbf{H}^{T}\mathbf{K}\_{1}^{2,1} + \mathbf{F}^{2}\mathbf{H}^{T}\mathbf{K}\_{1}^{2,2}\mathbf{H}\mathbf{F}^{1})\tag{22a}$$

$$\mathbf{G}^2 = -\left(\mathbf{K}\_1^{2,2} + \mathbf{K}\_2^{1,1} + \mathbf{F}^{2^T}\mathbf{H}^T\mathbf{K}\_1^{2,2} + \mathbf{F}^{2^T}\mathbf{H}^T\mathbf{K}\_1^{2,2}\mathbf{H}\mathbf{F}^2\right)^{-1}\mathbf{K}\_1^{1,2}\tag{22b}$$

The joint-arrays *θ<sup>j</sup>* becomes,

$$\boldsymbol{\theta} = \mathbf{Y}^1 \mathbf{u}\_1^1 + \mathbf{Y}^2 \mathbf{u}\_2^2 \tag{23}$$

$$\mathbf{Y}^1 = \mathbf{F}^1 + \mathbf{F}^2 \mathbf{G}^1, \mathbf{Y}^2 = \mathbf{F}^2 \mathbf{G}^2 \tag{24}$$

The dependent nodal-array **u**<sup>2</sup> <sup>1</sup> is obtained by substituting eq.(21) and eq.(23) into eq.(17):

$$\mathbf{u}\_1^2 = \mathbf{X}^1 \mathbf{u}\_1^1 + \mathbf{X}^2 \mathbf{u}\_2^2 \tag{25}$$

$$\mathbf{X}^1 = \mathbf{G}^1 + \mathbf{H}\mathbf{Y}^1,\\\mathbf{X}^2 = \mathbf{G}^2 + \mathbf{H}\mathbf{Y}^2 \tag{26}$$

Let ˜**u***<sup>b</sup>* = **u**1 1 *T* **u**2 2 *<sup>T</sup> <sup>T</sup>* be the 12-dimensional array of independent nodal displacements, then the final expression of *Vb* in terms ˜**u***<sup>b</sup>* is

$$V\_b = \frac{1}{2} \tilde{\mathbf{u}}\_b^T \mathbf{K}\_b \tilde{\mathbf{u}}\_b \tag{27}$$

$$\mathbf{K}\_{b} = \begin{bmatrix} \mathbf{1} & \mathbf{O} \\ \mathbf{X}^{1} \mathbf{X}^{2} \end{bmatrix}^{T} \begin{bmatrix} \mathbf{K}\_{1}^{1,1} \mathbf{K}\_{1}^{1,2} \\ \mathbf{K}\_{1}^{2,1} \mathbf{K}\_{1}^{2,2} \end{bmatrix} \begin{bmatrix} \mathbf{1} & \mathbf{O} \\ \mathbf{X}^{1} \mathbf{X}^{2} \end{bmatrix} + \begin{bmatrix} \mathbf{G}^{1} \mathbf{G}^{2} \\ \mathbf{O} \mathbf{1} \end{bmatrix}^{T} \begin{bmatrix} \mathbf{K}\_{2}^{1,1} \mathbf{K}\_{2}^{1,2} \\ \mathbf{K}\_{2}^{2,1} \mathbf{K}\_{2}^{2,2} \end{bmatrix} \begin{bmatrix} \mathbf{G}^{1} \mathbf{G}^{2} \\ \mathbf{O} \mathbf{1} \end{bmatrix} \tag{28}$$

where **K***<sup>b</sup>* is the 12 × 12 stiffness matrix.

Fig. 5. Case c): Flexible beam-rigid body.

#### **4.3 Case c) Flexible body-rigid body**

The case c) describes a beam coupled to a rigid body by means of a joint, as shown in Figure 5. The expressions involving the case c) can be easily obtained by extension of those of the case a), hence, we write only the final results for brevity:

$$\mathbf{Y}^1 = -(\mathbf{H}^T \mathbf{K}\_1^{2,2} \mathbf{H})^{-1} \mathbf{H}^T \mathbf{K}\_1^{2,1} \tag{29a}$$

$$\mathbf{Y}^2 = -(\mathbf{H}^T \mathbf{K}\_1^{2,2} \mathbf{H})^{-1} \mathbf{H}^T \mathbf{K}\_1^{2,2} \mathbf{G} \tag{29b}$$

$$\mathbf{X}^1 = \mathbf{H} \mathbf{Y}^1 \tag{30a}$$

for the real case **H** ≡ **1**6. The joint displacement array is now denoted with *δθ*.

configurations, we can write

**5. Mass matrix determination**

joints with mass and inertia.

**u**2 1

Fig. 6. Flexure revolute joint.

**5.1 Lumped approach**

simply reduces to

Referring to Figure 7, let us consider three different configurations of the flexure joints: an undeformed and unloaded configuration in which the joint-array is *θo*; a preloaded initial configuration with the joint-array *θ<sup>i</sup>* and a final configuration in which *θ<sup>f</sup>* = *θ<sup>i</sup>* + *δθ*, being *δθ* an array of small displacements around the initial joint-array *θi*. For the following explanation, we refer only to the ideal case of Figure 6. Let **K***<sup>θ</sup>* be the *m*(*j*) × *m*(*j*) joint stiffness matrix and let **w** be the generic wrench acting on the flexure joint; then, for the three

On the Stiffness Analysis and Elastodynamics of Parallel Kinematic Machines 95

where Δ*θ<sup>f</sup>* = Δ*θ<sup>i</sup>* + *δθ*. The strain energy *V<sup>θ</sup> <sup>f</sup>* of the flexure joint in its final configuration

The above expressions can be simplified considering *θ<sup>o</sup>* = **0**<sup>6</sup> for the unloaded configuration. We do not further discuss on flexure joints, leaving the reader to derive three new cases,

In this section two ways to include masses/inertias are discussed: the *lumped approach* and the *distributed approach*. The former concentrates masses on nodes of both rigid and flexible parts; the latter considers the real distribution of masses inside beams. As will be explained in the text, the two methods produce good results, particularly when an accurate degree of partitioning is chosen for flexible bodies. Finally, we focus our attention on a way to consider

Reducing mass and inertia of a rigid body to a particular point, the center of mass, without changing dynamic properties of the system, is a common procedure. On the contrary, for flexible bodies every reduction process is an approximation generating mistakes. Let us refer to Fig. 8 in which a link is divided into four beams. In the lumped approach the mass

> **u**2 1

**u**1 2

*<sup>V</sup><sup>θ</sup> <sup>f</sup>* <sup>=</sup> <sup>1</sup> 2 Δ*θ<sup>T</sup>*

similar to those discussed above, taking into account joint stiffness.

**h**(**e**, R)

**e**

*θ*

(a) Ideal flexure joint

**w***<sup>o</sup>* = **0**<sup>6</sup> (33a) **w***<sup>i</sup>* = **K***<sup>θ</sup>* (*θ<sup>i</sup>* − *θo*) ≡ **K***θ*Δ*θ<sup>i</sup>* (33b) **w***<sup>f</sup>* = **K***<sup>θ</sup>* (*θ<sup>f</sup>* − *θo*) ≡ **K***θ*Δ*θ<sup>f</sup>* (33c)

*<sup>f</sup>* **K***θ*Δ*θ<sup>f</sup>* (34)

*<sup>θ</sup><sup>x</sup> <sup>θ</sup><sup>y</sup>*

(b) Real flexure joint

*dz*

*θz dy dx* **u**1 2

$$\mathbf{X}^2 = \mathbf{G} + \mathbf{H}\mathbf{Y}^2\tag{30b}$$

The strain energy *Vc* associated to the beam becomes

$$V\_c = \frac{1}{2} \tilde{\mathbf{u}}\_c^T \mathbf{K}\_c \tilde{\mathbf{u}}\_{c\prime} \quad \mathbf{K}\_c = \begin{bmatrix} \mathbf{1} & \mathbf{O} \\ \mathbf{X}^1 \ \mathbf{X}^2 \end{bmatrix}^T \begin{bmatrix} \mathbf{K}\_1^{1,1} \ \mathbf{K}\_1^{1,2} \\ \mathbf{K}\_1^{2,1} \ \mathbf{K}\_1^{2,2} \end{bmatrix} \begin{bmatrix} \mathbf{1} & \mathbf{O} \\ \mathbf{X}^1 \ \mathbf{X}^2 \end{bmatrix} \tag{31}$$

where ˜**u***<sup>c</sup>* = **u**1 1 *T* **u**2 2 *<sup>T</sup> <sup>T</sup>* is the 12-dimensional array of independent nodal displacements and **K***<sup>c</sup>* is the 12 × 12 generalized stiffness matrix of the case c).

#### **4.4 Joint's stiffness and flexure joints**

The final part of the present section is devoted to show some feasible extension of the formulation to flexure mechanisms. Similar concepts can be applied even to ordinary joints to take into account joint stiffness. In order to reproduce the counterparts of mechanical joints, in a continuous structure it is a common strategy to recur to flexure joints, in fact, zones where the geometry and shape are designed to increase the compliance along specified degrees of freedom (*dofs*). An ideal flexure joint should allow for only motions along its *dofs*, while withstanding to remaining motions along its degrees of constraint (*docs*).

Figure 6 shows the cases of ideal and real flexure revolute joints. While for the ideal case the nodal displacements **u**<sup>2</sup> <sup>1</sup> and **<sup>u</sup>**<sup>1</sup> <sup>2</sup> of the two flexible beams (1) and (2), coupled by the joint, are tied by the usual constraint equation

$$\mathbf{u}\_1^2 = \mathbf{u}\_2^1 + \mathbf{H}\delta\theta \tag{32}$$

for the real case **H** ≡ **1**6. The joint displacement array is now denoted with *δθ*. Referring to Figure 7, let us consider three different configurations of the flexure joints: an undeformed and unloaded configuration in which the joint-array is *θo*; a preloaded initial configuration with the joint-array *θ<sup>i</sup>* and a final configuration in which *θ<sup>f</sup>* = *θ<sup>i</sup>* + *δθ*, being *δθ* an array of small displacements around the initial joint-array *θi*. For the following explanation, we refer only to the ideal case of Figure 6. Let **K***<sup>θ</sup>* be the *m*(*j*) × *m*(*j*) joint stiffness matrix and let **w** be the generic wrench acting on the flexure joint; then, for the three configurations, we can write

$$\mathbf{w}\_{\varnothing} = \mathbf{0}\_{\varnothing} \tag{33a}$$

$$\mathbf{w}\_{\dot{l}} = \mathbf{K}\_{\theta} (\theta\_{\dot{l}} - \theta\_{o}) \equiv \mathbf{K}\_{\theta} \Delta \theta\_{\dot{l}} \tag{33b}$$

$$\mathbf{w}\_f = \mathbf{K}\_\theta(\theta\_f - \theta\_o) \equiv \mathbf{K}\_\theta \Delta \theta\_f \tag{33c}$$

where Δ*θ<sup>f</sup>* = Δ*θ<sup>i</sup>* + *δθ*. The strain energy *V<sup>θ</sup> <sup>f</sup>* of the flexure joint in its final configuration simply reduces to

$$V\_{\theta\_f} = \frac{1}{2} \Delta \theta\_f^T \mathbf{K}\_{\theta} \Delta \theta\_f \tag{34}$$

The above expressions can be simplified considering *θ<sup>o</sup>* = **0**<sup>6</sup> for the unloaded configuration. We do not further discuss on flexure joints, leaving the reader to derive three new cases, similar to those discussed above, taking into account joint stiffness.

#### **5. Mass matrix determination**

In this section two ways to include masses/inertias are discussed: the *lumped approach* and the *distributed approach*. The former concentrates masses on nodes of both rigid and flexible parts; the latter considers the real distribution of masses inside beams. As will be explained in the text, the two methods produce good results, particularly when an accurate degree of partitioning is chosen for flexible bodies. Finally, we focus our attention on a way to consider joints with mass and inertia.

#### **5.1 Lumped approach**

10 Will-be-set-by-IN-TECH

**u**1 2 *θ*

2

The case c) describes a beam coupled to a rigid body by means of a joint, as shown in Figure 5. The expressions involving the case c) can be easily obtained by extension of those of the case

> **1 O X**<sup>1</sup> **X**<sup>2</sup>

The final part of the present section is devoted to show some feasible extension of the formulation to flexure mechanisms. Similar concepts can be applied even to ordinary joints to take into account joint stiffness. In order to reproduce the counterparts of mechanical joints, in a continuous structure it is a common strategy to recur to flexure joints, in fact, zones where the geometry and shape are designed to increase the compliance along specified degrees of freedom (*dofs*). An ideal flexure joint should allow for only motions along its *dofs*, while

Figure 6 shows the cases of ideal and real flexure revolute joints. While for the ideal case the

<sup>1</sup> **<sup>H</sup>**)−1**H***T***K**2,1

<sup>1</sup> **<sup>H</sup>**)−1**H***T***K**2,2

*<sup>T</sup>* **K**1,1 <sup>1</sup> **<sup>K</sup>**1,2 1

**K**2,1 <sup>1</sup> **<sup>K</sup>**2,2 1

is the 12-dimensional array of independent nodal displacements and

<sup>2</sup> of the two flexible beams (1) and (2), coupled by the joint, are

<sup>2</sup> + **H***δθ* (32)

<sup>1</sup> (29a)

<sup>1</sup> **G** (29b)

(31)

**X**<sup>1</sup> = **HY**<sup>1</sup> (30a) **X**<sup>2</sup> = **G** + **HY**<sup>2</sup> (30b)

> **1 O X**<sup>1</sup> **X**<sup>2</sup>

**<sup>Y</sup>**<sup>1</sup> <sup>=</sup> <sup>−</sup>(**H***T***K**2,2

**<sup>Y</sup>**<sup>2</sup> <sup>=</sup> <sup>−</sup>(**H***T***K**2,2

**u**2 1

F

**u**1 1

1

Fig. 5. Case c): Flexible beam-rigid body.

a), hence, we write only the final results for brevity:

The strain energy *Vc* associated to the beam becomes

**K***<sup>c</sup>* is the 12 × 12 generalized stiffness matrix of the case c).

<sup>1</sup> and **<sup>u</sup>**<sup>1</sup>

*<sup>c</sup>* **K***c***u**˜ *<sup>c</sup>*, **K***<sup>c</sup>* =

withstanding to remaining motions along its degrees of constraint (*docs*).

**u**2 <sup>1</sup> <sup>=</sup> **<sup>u</sup>**<sup>1</sup>

*Vc* <sup>=</sup> <sup>1</sup> 2 **u**˜ *T*

**4.4 Joint's stiffness and flexure joints**

tied by the usual constraint equation

where ˜**u***<sup>c</sup>* =

 **u**1 1 *T* **u**2 2 *<sup>T</sup> <sup>T</sup>*

nodal displacements **u**<sup>2</sup>

**4.3 Case c) Flexible body-rigid body**

**u**2 2

**d**

R

**e**

Reducing mass and inertia of a rigid body to a particular point, the center of mass, without changing dynamic properties of the system, is a common procedure. On the contrary, for flexible bodies every reduction process is an approximation generating mistakes. Let us refer to Fig. 8 in which a link is divided into four beams. In the lumped approach the mass

Fig. 6. Flexure revolute joint.

*m*/2

4

<sup>210</sup> <sup>−</sup> *Iy*

<sup>3</sup>*ρAi* 0 0

<sup>210</sup> <sup>+</sup> *Iz* 10*ρAiLi* ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦

> ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦

> > ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦

(39)

(40)

(37)

(38)

<sup>10</sup>*ρAiLi* 0

<sup>15</sup>*ρAi* 0

<sup>10</sup>*ρAiLi* 0

<sup>30</sup>*ρAi* 0

<sup>10</sup>*ρAiLi* 0

<sup>15</sup>*ρAi* 0

*i* <sup>105</sup> <sup>+</sup> <sup>2</sup>*Iz* 15*ρAi*

*i* <sup>140</sup> <sup>−</sup> *Iz* 30*ρAi*

<sup>210</sup> <sup>−</sup> *Iz* 10*ρAiLi*

*i* <sup>105</sup> <sup>+</sup> <sup>2</sup>*Iz* 15*ρAi*

<sup>420</sup> <sup>+</sup> *Iz* 10*ρAiLi*

*m*/2

On the Stiffness Analysis and Elastodynamics of Parallel Kinematic Machines 97

3

3

<sup>3</sup> 0 000 0

<sup>10</sup>*ρAiLi* <sup>0</sup> *<sup>L</sup>*<sup>2</sup>

<sup>10</sup>*ρAiLi* <sup>000</sup> *<sup>L</sup>*<sup>2</sup>

<sup>10</sup>*ρAiLi* <sup>0</sup> <sup>−</sup>*L*<sup>2</sup>

<sup>10</sup>*ρAiLi* 0 00 <sup>−</sup>*L*<sup>2</sup>

<sup>10</sup>*ρAiLi* <sup>0</sup> *<sup>L</sup>*<sup>2</sup>

<sup>10</sup>*ρAiLi* <sup>000</sup> *<sup>L</sup>*<sup>2</sup>

<sup>3</sup> 0 000 0

<sup>6</sup> 0 0 00 0

<sup>35</sup> <sup>+</sup> <sup>6</sup>*Iz* 5*ρAiL*<sup>2</sup> *i*

<sup>210</sup> <sup>−</sup> *Iy*

<sup>70</sup> <sup>−</sup> <sup>6</sup>*Iy* 5*ρAiL*<sup>2</sup> *i*

<sup>420</sup> <sup>+</sup> *Iz*

<sup>35</sup> <sup>+</sup> <sup>6</sup>*Iy* 5*ρAiL*<sup>2</sup> *i*

<sup>210</sup> <sup>+</sup> *Iy*

*<sup>i</sup>* of the nodal displacement arrays **<sup>u</sup>**<sup>1</sup>

�*<sup>T</sup>* �

where *ρ*, *Li*, *Ai*, *Ix*, *Iy* and *Iz*, respectively, are the density, the length, the cross section area and the mass moments of inertia for unit of length of the *i*th-beam. The matrix **M***<sup>i</sup>* is associated to the kinetic energy *Ti* of a beam, the latter being defined as a quadratic forms into the

> **M**1,1 *<sup>i</sup>* **<sup>M</sup>**1,2 *i*

> **M**2,1 *<sup>i</sup>* **<sup>M</sup>**2,2 *i*

00 0 *Ix*

00 0 *Ix*

00 0 *Ix*

4

000 <sup>11</sup>*Li*

*i* <sup>105</sup> <sup>+</sup> <sup>2</sup>*Iy*

0 00 <sup>−</sup>13*Li*

<sup>420</sup> <sup>−</sup> *Iy*

*i* <sup>140</sup> <sup>−</sup> *Iy*

<sup>000</sup> <sup>−</sup>11*Li*

*i* <sup>105</sup> <sup>+</sup> <sup>2</sup>*Iy*

*<sup>i</sup>* and **<sup>u</sup>**<sup>2</sup> *i* :

�

� � **u**˙ <sup>1</sup> *i* **u**˙ 2 *i*

<sup>210</sup> <sup>+</sup> *Iy*

<sup>3</sup>*ρAi* 0 0

0 <sup>11</sup>*Li*

<sup>3</sup>*ρAi* 0 0

0 <sup>13</sup>*<sup>L</sup>*

0 <sup>−</sup>11*Li*

1 *m*

2

2

1

Fig. 8. Lumped distribution of mass.

⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ 1

0 <sup>13</sup>

0 <sup>11</sup>*Li*

0 <sup>9</sup>

0 <sup>13</sup>*Li*

0 <sup>13</sup>

<sup>0</sup> <sup>−</sup>11*Li*

*<sup>i</sup>* , ˙**u**<sup>2</sup>

⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ 1

⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ 1

<sup>35</sup> <sup>+</sup> <sup>6</sup>*Iz* 5*ρAiL*<sup>2</sup> *i*

<sup>210</sup> <sup>+</sup> *Iz*

<sup>70</sup> <sup>−</sup> <sup>6</sup>*Iz* 5*ρAiL*<sup>2</sup> *i*

<sup>420</sup> <sup>−</sup> *Iz*

<sup>35</sup> <sup>+</sup> <sup>6</sup>*Iz* 5*ρAiL*<sup>2</sup> *i*

0 0 <sup>13</sup>

0 0 <sup>11</sup>*Li*

*Ti* <sup>=</sup> <sup>1</sup> 2 � **u**˙ 1 *i* **u**˙ 2 *i*

<sup>210</sup> <sup>−</sup> *Iz*

0 0 <sup>9</sup>

0 0 <sup>−</sup>13*Li*

0 0 <sup>13</sup>

0 0 <sup>−</sup>11*Li*

**M**1,1

**M**1,2

**M**2,2

*<sup>i</sup>* = *ρAiLi*

*<sup>i</sup>* = *ρAiLi*

*<sup>i</sup>* = *ρAiLi*

time-derivatives ˙**u**<sup>1</sup>

Fig. 7. Flexure joint: a) Undeformed configuration; b) starting preloaded configuration; c) small rotated configuration about the starting preloaded configuration.

and inertia of each beam is concentrated on its end-nodes. Here, we choose a symmetric distribution in which the mass *m* is divided by two, but the reader can use another distribution according to the case to be examined, as instance beams with varying cross section. The first node has a mass of *m*/2 while any other node, but the fourth, bears a mass *m* because it receives half a mass from each of the two contiguous beams coupled at its section. The fourth beam of the link is attached to a joint. According to what explained in the previous section, even the mass is concentrated only on the independent joints: it means that the mass of the last beam has to be concentrated on the fourth node, thereby the latter carrying a mass equal to 1/2*m* + *m* = 3/2*m*. Analogous arguments, not reported here for conciseness, may be used to describe inertias.

The lumped approach concentrates masses and inertias on all the independent nodes, belonging to both rigid and flexible bodies. Mathematically, a 6× 6 mass dyad **M***<sup>i</sup>* is associated at the *i*th-independent node, defined as

$$
\widetilde{\mathbf{M}}\_i = \begin{bmatrix} m\_i \mathbf{1} \bullet \mathbf{O} \\ \mathbf{O} \quad \mathbf{J}\_i \end{bmatrix} \tag{35}
$$

where *mi* and **J***<sup>i</sup>* are the mass and the 3 × 3 inertia matrix, respectively, of either a rigid body, if the independent node is at the center of mass of the said body, or of a beam's end-node. The generalized inertia matrix **M***PKM* of a PKM is readily derived upon assembling in diagonal blocks the previous inertia dyads, following the order chosen to enumerate all independent nodes inside the global array of independent nodal coordinates, hence

$$\mathbf{M}\_{PKM} = \operatorname{diag}(\mathbf{M}\_1, \mathbf{M}\_2, \dots, \mathbf{M}\_n) \tag{36}$$

#### **5.2 Distributed approach**

The distributed approach considers the true distribution of mass inside a flexible beam. Particularly, it is possible to find a 12 × 12 matrix **M***<sup>i</sup>* referred to the twelve nodal coordinates of a beam's end-nodes. This matrix can be divided into blocks, whose entries are reported below, as already done for the stiffness matrix **K***i*, i.e.

Fig. 8. Lumped distribution of mass.

12 Will-be-set-by-IN-TECH

**w***f*

**w***f*

(35)

**w***i*

Fig. 7. Flexure joint: a) Undeformed configuration; b) starting preloaded configuration; c)

and inertia of each beam is concentrated on its end-nodes. Here, we choose a symmetric distribution in which the mass *m* is divided by two, but the reader can use another distribution according to the case to be examined, as instance beams with varying cross section. The first node has a mass of *m*/2 while any other node, but the fourth, bears a mass *m* because it receives half a mass from each of the two contiguous beams coupled at its section. The fourth beam of the link is attached to a joint. According to what explained in the previous section, even the mass is concentrated only on the independent joints: it means that the mass of the last beam has to be concentrated on the fourth node, thereby the latter carrying a mass equal to 1/2*m* + *m* = 3/2*m*. Analogous arguments, not reported here for conciseness, may be used

The lumped approach concentrates masses and inertias on all the independent nodes, belonging to both rigid and flexible bodies. Mathematically, a 6× 6 mass dyad **M***<sup>i</sup>* is associated

where *mi* and **J***<sup>i</sup>* are the mass and the 3 × 3 inertia matrix, respectively, of either a rigid body, if the independent node is at the center of mass of the said body, or of a beam's end-node. The generalized inertia matrix **M***PKM* of a PKM is readily derived upon assembling in diagonal blocks the previous inertia dyads, following the order chosen to enumerate all independent

The distributed approach considers the true distribution of mass inside a flexible beam. Particularly, it is possible to find a 12 × 12 matrix **M***<sup>i</sup>* referred to the twelve nodal coordinates of a beam's end-nodes. This matrix can be divided into blocks, whose entries are reported

**M***PKM* = *diag*(**M**1, **M**2,..., **M***n*) (36)

**M***<sup>i</sup>* = *mi***1 O O J***<sup>i</sup>*

nodes inside the global array of independent nodal coordinates, hence

below, as already done for the stiffness matrix **K***i*, i.e.

small rotated configuration about the starting preloaded configuration.

*<sup>θ</sup><sup>o</sup> <sup>θ</sup><sup>i</sup> <sup>θ</sup><sup>f</sup>*

**w***<sup>o</sup>* **w***<sup>i</sup>*

**w***o*

to describe inertias.

**5.2 Distributed approach**

at the *i*th-independent node, defined as

**M**1,1 *<sup>i</sup>* = *ρAiLi* ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ 1 <sup>3</sup> 0 000 0 0 <sup>13</sup> <sup>35</sup> <sup>+</sup> <sup>6</sup>*Iz* 5*ρAiL*<sup>2</sup> *i* 000 <sup>11</sup>*Li* <sup>210</sup> <sup>+</sup> *Iz* 10*ρAiLi* 0 0 <sup>13</sup> <sup>35</sup> <sup>+</sup> <sup>6</sup>*Iz* 5*ρAiL*<sup>2</sup> *i* 0 <sup>−</sup>11*Li* <sup>210</sup> <sup>−</sup> *Iy* <sup>10</sup>*ρAiLi* 0 00 0 *Ix* <sup>3</sup>*ρAi* 0 0 0 0 <sup>−</sup>11*Li* <sup>210</sup> <sup>−</sup> *Iy* <sup>10</sup>*ρAiLi* <sup>0</sup> *<sup>L</sup>*<sup>2</sup> *i* <sup>105</sup> <sup>+</sup> <sup>2</sup>*Iy* <sup>15</sup>*ρAi* 0 0 <sup>11</sup>*Li* <sup>210</sup> <sup>+</sup> *Iz* <sup>10</sup>*ρAiLi* <sup>000</sup> *<sup>L</sup>*<sup>2</sup> *i* <sup>105</sup> <sup>+</sup> <sup>2</sup>*Iz* 15*ρAi* ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ (37) **M**1,2 *<sup>i</sup>* = *ρAiLi* ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ 1 <sup>6</sup> 0 0 00 0 0 <sup>9</sup> <sup>70</sup> <sup>−</sup> <sup>6</sup>*Iz* 5*ρAiL*<sup>2</sup> *i* 0 00 <sup>−</sup>13*Li* <sup>420</sup> <sup>+</sup> *Iz* 10*ρAiLi* 0 0 <sup>9</sup> <sup>70</sup> <sup>−</sup> <sup>6</sup>*Iy* 5*ρAiL*<sup>2</sup> *i* 0 <sup>13</sup>*<sup>L</sup>* <sup>420</sup> <sup>−</sup> *Iy* <sup>10</sup>*ρAiLi* 0 00 0 *Ix* <sup>3</sup>*ρAi* 0 0 0 0 <sup>−</sup>13*Li* <sup>420</sup> <sup>+</sup> *Iz* <sup>10</sup>*ρAiLi* <sup>0</sup> <sup>−</sup>*L*<sup>2</sup> *i* <sup>140</sup> <sup>−</sup> *Iy* <sup>30</sup>*ρAi* 0 0 <sup>13</sup>*Li* <sup>420</sup> <sup>−</sup> *Iz* <sup>10</sup>*ρAiLi* 0 00 <sup>−</sup>*L*<sup>2</sup> *i* <sup>140</sup> <sup>−</sup> *Iz* 30*ρAi* ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ (38) **M**2,2 *<sup>i</sup>* = *ρAiLi* ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ 1 <sup>3</sup> 0 000 0 0 <sup>13</sup> <sup>35</sup> <sup>+</sup> <sup>6</sup>*Iz* 5*ρAiL*<sup>2</sup> *i* <sup>000</sup> <sup>−</sup>11*Li* <sup>210</sup> <sup>−</sup> *Iz* 10*ρAiLi* 0 0 <sup>13</sup> <sup>35</sup> <sup>+</sup> <sup>6</sup>*Iy* 5*ρAiL*<sup>2</sup> *i* 0 <sup>11</sup>*Li* <sup>210</sup> <sup>+</sup> *Iy* <sup>10</sup>*ρAiLi* 0 00 0 *Ix* <sup>3</sup>*ρAi* 0 0 0 0 <sup>11</sup>*Li* <sup>210</sup> <sup>+</sup> *Iy* <sup>10</sup>*ρAiLi* <sup>0</sup> *<sup>L</sup>*<sup>2</sup> *i* <sup>105</sup> <sup>+</sup> <sup>2</sup>*Iy* <sup>15</sup>*ρAi* 0 <sup>0</sup> <sup>−</sup>11*Li* <sup>210</sup> <sup>−</sup> *Iz* <sup>10</sup>*ρAiLi* <sup>000</sup> *<sup>L</sup>*<sup>2</sup> *i* <sup>105</sup> <sup>+</sup> <sup>2</sup>*Iz* 15*ρAi* ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ (39)

where *ρ*, *Li*, *Ai*, *Ix*, *Iy* and *Iz*, respectively, are the density, the length, the cross section area and the mass moments of inertia for unit of length of the *i*th-beam. The matrix **M***<sup>i</sup>* is associated to the kinetic energy *Ti* of a beam, the latter being defined as a quadratic forms into the time-derivatives ˙**u**<sup>1</sup> *<sup>i</sup>* , ˙**u**<sup>2</sup> *<sup>i</sup>* of the nodal displacement arrays **<sup>u</sup>**<sup>1</sup> *<sup>i</sup>* and **<sup>u</sup>**<sup>2</sup> *i* :

$$T\_i = \frac{1}{2} \begin{bmatrix} \dot{\mathbf{u}}\_i^1 \\ \dot{\mathbf{u}}\_i^2 \end{bmatrix}^T \begin{bmatrix} \mathbf{M}\_i^{1,1} \ \mathbf{M}\_i^{1,2} \\ \mathbf{M}\_i^{2,1} \ \mathbf{M}\_i^{2,2} \end{bmatrix} \begin{bmatrix} \dot{\mathbf{u}}\_i^1 \\ \dot{\mathbf{u}}\_i^2 \end{bmatrix} \tag{40}$$

In the previous section we have found how dependent nodal displacement arrays are expressed in terms of independent ones. Let us consider the time-derivative of both sides of eq.(13), we write

$$
\dot{\mathbf{u}}\_2^1 = \mathbf{X}^1 \dot{\mathbf{u}}\_1^1 + \mathbf{X}^2 \dot{\mathbf{u}}\_2^2 \tag{41}
$$

Then, upon recalling eqs.(21) and (26), *TJ* can be expressed in terms of ˜

independent nodal displacements, respectively, defined as

**M***<sup>a</sup>* =

**M***<sup>b</sup>* =

**M***<sup>c</sup>* =

for the three cases in exam.

**u**<sup>1</sup> **u**<sup>2</sup> ... **u***<sup>n</sup>*

*B*1(*i*) =

**q** =

**6. Linearized elastodynamics equations**

**u***<sup>i</sup>* = **B**1(*i*)**q**, ˜**u** ≡

final dimension of the problem, defined as

**q**. The 6 × 6*n* matrix **B**<sup>1</sup> and the 12 × 6*n* matrix **B**<sup>2</sup> are defined as:

 **u***i* **u***j* 

**O**<sup>6</sup> **O**<sup>6</sup> ... **1**6(*i*) ... **O**<sup>6</sup>

*TJ* <sup>=</sup> <sup>1</sup> 2 (˜ **u**˙ )*T***M***<sup>J</sup>* ˜

**<sup>X</sup>**1*<sup>T</sup>***<sup>M</sup>***R***X**<sup>1</sup> <sup>+</sup> **<sup>G</sup>***<sup>T</sup>***<sup>M</sup>***L***G X**1*<sup>T</sup>***<sup>M</sup>***R***X**<sup>2</sup> **<sup>X</sup>**2*<sup>T</sup>***<sup>M</sup>***R***X**<sup>1</sup> **<sup>X</sup>**2*<sup>T</sup>***<sup>M</sup>***R***X**<sup>2</sup>

**<sup>X</sup>**1*<sup>T</sup>***<sup>M</sup>***L***X**<sup>1</sup> **<sup>X</sup>**1*<sup>T</sup>***<sup>M</sup>***L***X**<sup>2</sup> **<sup>X</sup>**2*<sup>T</sup>***<sup>M</sup>***L***X**<sup>1</sup> **<sup>X</sup>**2*<sup>T</sup>***<sup>M</sup>***L***X**<sup>2</sup> <sup>+</sup> **<sup>G</sup>***<sup>T</sup>***<sup>M</sup>***R***<sup>G</sup>**

where **M***<sup>J</sup>* is the 12 × 12 generalized inertia matrix of the joint expressed in terms of

On the Stiffness Analysis and Elastodynamics of Parallel Kinematic Machines 99

**<sup>X</sup>**1*<sup>T</sup>***<sup>M</sup>***L***X**<sup>1</sup> <sup>+</sup> **<sup>G</sup>**1*<sup>T</sup>***<sup>M</sup>***R***G**<sup>1</sup> **<sup>X</sup>**1*<sup>T</sup>***<sup>M</sup>***L***X**<sup>2</sup> <sup>+</sup> **<sup>G</sup>**1*<sup>T</sup>***<sup>M</sup>***R***G**<sup>2</sup> **<sup>X</sup>**2*<sup>T</sup>***<sup>M</sup>***L***X**<sup>1</sup> <sup>+</sup> **<sup>G</sup>**2*<sup>T</sup>***<sup>M</sup>***R***G**<sup>1</sup> **<sup>X</sup>**2*<sup>T</sup>***<sup>M</sup>***L***X**<sup>2</sup> <sup>+</sup> **<sup>G</sup>**2*<sup>T</sup>***<sup>M</sup>***R***G**<sup>2</sup>

In this section we will derive the generalized inertia and stiffness matrices necessary to write the linearized elastodynamics equations. As described in the two previous sections stiffness and inertia matrices of rigid bodies and flexible beams are referred to the independent nodal displacements of the system. Now, one might ask how to combine different matrices to build the global stiffness and inertia matrices. In order to solve this issue let us introduce two Boolean matrices **B**<sup>1</sup> and **B**<sup>2</sup> able to map a 6-dimensional displacement array **u***<sup>i</sup>* or a 12-dimensional displacement array ˜**u** in terms of a 6*n*-dimensional global array

important that the reader would define the order in which every single array **u***<sup>i</sup>* appears inside

The above expressions allow us to convert each nodal displacement array in term of **q**. As instance, let us take into exam the strain energy *Vb* of eq.(27), it simply turns into *Vb* <sup>=</sup> <sup>1</sup>

From *Vb* we can extract a new 6*n* × 6*n* matrix **K***b*, that is a stiffness matrix expanded to the

Likewise, the generic expanded 6*n* × 6*n* mass matrix **M***<sup>i</sup>* of **M***i*, see eq.(35), can be written as **<sup>M</sup>***<sup>i</sup>* <sup>=</sup> **<sup>B</sup>**1(*i*)*<sup>T</sup>***<sup>M</sup>***i***B**1(*i*), where *<sup>i</sup>* is the position index of **<sup>u</sup>***i*, i.e. the displacement array of the *i*th-rigid body's center of mass, inside **q**. By the same strategy, it is possible to expand every stiffness and inertia matrix. This operation is essential as, *only when referred to the same set*

, *<sup>B</sup>*2(*i*, *<sup>j</sup>*) =

<sup>2</sup>**q***T***B**2(*i*, *<sup>j</sup>*)*T***K***b***B**2(*i*, *<sup>j</sup>*)**q**, where *<sup>i</sup>* and *<sup>j</sup>* indicate the position indices of **<sup>u</sup>**<sup>1</sup>

containing all the independent displacement arrays of the robot. It is

= **B**2(*i*, *j*)**q** (47)

**K***<sup>b</sup>* = **B**2(*i*, *j*)*T***K***b***B**2(*i*, *j*) (49)

**O**<sup>6</sup> **O**<sup>6</sup> ... **1**6(*i*) ... **O**<sup>6</sup> **O**<sup>6</sup> **1**6(*j*) ... **O**<sup>6</sup> ... **O**<sup>6</sup> <sup>1</sup> and **<sup>u</sup>**<sup>2</sup>

(48)

<sup>2</sup> inside **q**.

**u**˙ = **u**1 1 *T* **u**2 2 *<sup>T</sup> <sup>T</sup>* , i.e.

**u**˙ (45)

(46a)

(46b)

(46c)

in which the matrices **X**<sup>1</sup> and **X**<sup>2</sup> are not dependent on time. Similar expressions can be obtained for eqs.(21) and (26) for the case b) and for the case c). It means that the same expressions standing for dependent and independent displacements can be extended at velocity and acceleration level. Therefore, the following matrices **M***a*, **M***<sup>b</sup>* and **M***<sup>c</sup>* can be written with perfect analogy to their counterparts **K***a*, **K***<sup>b</sup>* and **K***c*:

$$\mathbf{M}\_{d} = \begin{bmatrix} \mathbf{X}^{1} \mathbf{X}^{2} \\ \mathbf{O} \mathbf{ 1} \end{bmatrix}^{T} \begin{bmatrix} \mathbf{M}\_{2}^{1,1} \mathbf{M}\_{2}^{1,2} \\ \mathbf{M}\_{2}^{2,1} \mathbf{M}\_{2}^{2,2} \end{bmatrix} \begin{bmatrix} \mathbf{X}^{1} \mathbf{X}^{2} \\ \mathbf{O} \cdot \mathbf{1} \end{bmatrix} \tag{42a}$$

$$\mathbf{M}\_{b} = \begin{bmatrix} \mathbf{1} & \mathbf{O} \\ \mathbf{X}^{1} \mathbf{X}^{2} \end{bmatrix}^{T} \begin{bmatrix} \mathbf{M}\_{1}^{1,1} \mathbf{M}\_{1}^{1,2} \\ \mathbf{M}\_{1}^{2,1} \mathbf{M}\_{1}^{2,2} \end{bmatrix} \begin{bmatrix} \mathbf{1} & \mathbf{O} \\ \mathbf{X}^{1} \mathbf{X}^{2} \end{bmatrix} + \begin{bmatrix} \mathbf{G}^{1} \mathbf{G}^{2} \\ \mathbf{O} \mathbf{1} \end{bmatrix}^{T} \begin{bmatrix} \mathbf{M}\_{2}^{1,1} \mathbf{M}\_{2}^{1,2} \\ \mathbf{M}\_{2}^{2,1} \mathbf{M}\_{2}^{2,2} \end{bmatrix} \begin{bmatrix} \mathbf{G}^{1} \mathbf{G}^{2} \\ \mathbf{O} \mathbf{1} \end{bmatrix} \tag{42b}$$

$$\mathbf{M}\_{\mathbb{C}} = \begin{bmatrix} \mathbf{1} & \mathbf{O} \\ \mathbf{X}^{1} \mathbf{X}^{2} \end{bmatrix}^{T} \begin{bmatrix} \mathbf{M}\_{1}^{1,1} \mathbf{M}\_{1}^{1,2} \\ \mathbf{M}\_{1}^{2,1} \mathbf{M}\_{1}^{2,2} \end{bmatrix} \begin{bmatrix} \mathbf{1} & \mathbf{O} \\ \mathbf{X}^{1} \mathbf{X}^{2} \end{bmatrix} \tag{42c}$$

with obvious meaning of all terms. The three mass matrices, above defined, are referred to the time-derivatives ˜ **u**˙ *<sup>a</sup>*, ˜ **u**˙ *<sup>b</sup>* and ˜ **u**˙ *<sup>c</sup>* of the independent nodal displacement arrays ˜**u***a*, ˜**u***<sup>b</sup>* and ˜**u***c*. To clarify some doubt that might arise let us refer to Fig. 3. In this case the mass matrix of the rigid body, see eq.(35), is referred to its center of mass with displacement array **u**<sup>1</sup> <sup>1</sup>. The distributed approach first allows us to find the 12 × 12 mass matrix **M**<sup>2</sup> of the beam, then the latter is expressed in terms of only independent displacements by means of **M***a*. Obviously, **M**<sup>2</sup> and **M***<sup>a</sup>* refer to the same object, the beam; but **M***<sup>a</sup> distributes* **M**<sup>2</sup> into the independent nodes with displacement arrays **u**<sup>1</sup> <sup>1</sup> and **<sup>u</sup>**<sup>2</sup> <sup>2</sup>. This implies that the center of mass of the rigid body carries its own mass/inertia and part of the mass/inertia of the beam. Similar conclusions may be sought for the cases b) and c).

#### **5.3 Joint's mass and inertia**

In this final subsection we describe a method to consider mass and inertia of joints. For convenience, let us refer to Fig. 4. The joint is split into two parts, one belonging to the first beam, the remaining to the second beam. Let **M***<sup>L</sup>* and **M***R*, where capital letters stand for left and right, be the mass dyads of the two half-parts of the joint, defined as

$$
\widetilde{\mathbf{M}}\_L = \begin{bmatrix} m\_L \mathbf{1} \bullet \mathbf{O} \\ \mathbf{O} \quad \mathbf{J}\_L \end{bmatrix} \quad \widetilde{\mathbf{M}}\_R = \begin{bmatrix} m\_R \mathbf{1} \bullet \mathbf{O} \\ \mathbf{O} \quad \mathbf{J}\_R \end{bmatrix} \tag{43}
$$

The kinetic energy *TJ* of the joint can be written in terms of the above matrices **M***<sup>L</sup>* and **M***<sup>R</sup>* and of the time-derivatives of the nodal arrays ˙**u**<sup>2</sup> 1, ˙**u**<sup>1</sup> <sup>2</sup>, thus

$$T\_I = \frac{1}{2} (\dot{\mathbf{u}}\_1^2)^T \widetilde{\mathbf{M}}\_L \dot{\mathbf{u}}\_1^2 + \frac{1}{2} (\dot{\mathbf{u}}\_2^1)^T \widetilde{\mathbf{M}}\_R \dot{\mathbf{u}}\_2^1 \tag{44}$$

Then, upon recalling eqs.(21) and (26), *TJ* can be expressed in terms of ˜ **u**˙ = **u**1 1 *T* **u**2 2 *<sup>T</sup> <sup>T</sup>* , i.e.

$$T\_I = \frac{1}{2} (\mathbf{\check{u}})^T \mathbf{M}\_I \mathbf{\check{u}} \tag{45}$$

where **M***<sup>J</sup>* is the 12 × 12 generalized inertia matrix of the joint expressed in terms of independent nodal displacements, respectively, defined as

$$\mathbf{M}\_{\mathfrak{d}} = \begin{bmatrix} \mathbf{X}^{1} \widetilde{\mathbf{M}}\_{\mathbb{R}} \mathbf{X}^{1} + \mathbf{G}^{T} \widetilde{\mathbf{M}}\_{L} \mathbf{G} & \mathbf{X}^{1} \widetilde{\mathbf{M}}\_{\mathbb{R}} \mathbf{X}^{2} \\\ \mathbf{X}^{2} \widetilde{\mathbf{M}}\_{\mathbb{R}} \mathbf{X}^{1} & \mathbf{X}^{2} \widetilde{\mathbf{M}}\_{\mathbb{R}} \mathbf{X}^{2} \end{bmatrix} \tag{46a}$$

$$\mathbf{M}\_{b} = \begin{bmatrix} \mathbf{X}^{1}\widetilde{\mathbf{M}}\_{L}\mathbf{X}^{1} + \mathbf{G}^{1}\widetilde{\mathbf{M}}\_{R}\mathbf{G}^{1} & \mathbf{X}^{1}\widetilde{\mathbf{M}}\_{L}\mathbf{X}^{2} + \mathbf{G}^{1}\widetilde{\mathbf{M}}\_{R}\mathbf{G}^{2} \\ \mathbf{X}^{2}\widetilde{\mathbf{M}}\_{L}\mathbf{X}^{1} + \mathbf{G}^{2}\widetilde{\mathbf{M}}\_{R}\mathbf{G}^{1} & \mathbf{X}^{2}\widetilde{\mathbf{M}}\_{L}\mathbf{X}^{2} + \mathbf{G}^{2}\widetilde{\mathbf{M}}\_{R}\mathbf{G}^{2} \end{bmatrix} \tag{46b}$$

$$\mathbf{M}\_{\mathbb{C}} = \begin{bmatrix} \mathbf{X}^{1} \widetilde{\mathbf{M}}\_{L} \mathbf{X}^{1} & \mathbf{X}^{1} \widetilde{\mathbf{M}}\_{L} \mathbf{X}^{2} \\ \mathbf{X}^{2} \widetilde{\mathbf{M}}\_{L} \mathbf{X}^{1} & \mathbf{X}^{2} \widetilde{\mathbf{M}}\_{L} \mathbf{X}^{2} + \mathbf{G}^{T} \widetilde{\mathbf{M}}\_{R} \mathbf{G} \end{bmatrix} \tag{46c}$$

for the three cases in exam.

14 Will-be-set-by-IN-TECH

In the previous section we have found how dependent nodal displacement arrays are expressed in terms of independent ones. Let us consider the time-derivative of both sides

in which the matrices **X**<sup>1</sup> and **X**<sup>2</sup> are not dependent on time. Similar expressions can be obtained for eqs.(21) and (26) for the case b) and for the case c). It means that the same expressions standing for dependent and independent displacements can be extended at velocity and acceleration level. Therefore, the following matrices **M***a*, **M***<sup>b</sup>* and **M***<sup>c</sup>* can be

 + **G**<sup>1</sup> **G**<sup>2</sup> **O 1**

with obvious meaning of all terms. The three mass matrices, above defined, are referred to the

clarify some doubt that might arise let us refer to Fig. 3. In this case the mass matrix of the rigid

approach first allows us to find the 12 × 12 mass matrix **M**<sup>2</sup> of the beam, then the latter is expressed in terms of only independent displacements by means of **M***a*. Obviously, **M**<sup>2</sup> and **M***<sup>a</sup>* refer to the same object, the beam; but **M***<sup>a</sup> distributes* **M**<sup>2</sup> into the independent nodes with

its own mass/inertia and part of the mass/inertia of the beam. Similar conclusions may be

In this final subsection we describe a method to consider mass and inertia of joints. For convenience, let us refer to Fig. 4. The joint is split into two parts, one belonging to the first beam, the remaining to the second beam. Let **M***<sup>L</sup>* and **M***R*, where capital letters stand for left

<sup>1</sup>)*T***M***L***u**˙ <sup>2</sup>

The kinetic energy *TJ* of the joint can be written in terms of the above matrices **M***<sup>L</sup>* and **M***<sup>R</sup>*

**M***<sup>R</sup>* =

1, ˙**u**<sup>1</sup> <sup>2</sup>, thus

<sup>1</sup> + 1 2 (**u**˙ <sup>1</sup>

 *mR***1 O O J***<sup>R</sup>*

<sup>2</sup>)*T***M***R***u**˙ <sup>1</sup>

*<sup>T</sup>*

**u**˙ *<sup>c</sup>* of the independent nodal displacement arrays ˜**u***a*, ˜**u***<sup>b</sup>* and ˜**u***c*. To

<sup>2</sup>. This implies that the center of mass of the rigid body carries

**M**1,1 <sup>2</sup> **<sup>M</sup>**1,2 2

**M**2,1 <sup>2</sup> **<sup>M</sup>**2,2 2

<sup>1</sup> <sup>+</sup> **<sup>X</sup>**2**u**˙ <sup>2</sup>

<sup>2</sup> (41)

**G**<sup>1</sup> **G**<sup>2</sup> **O 1**

<sup>1</sup>. The distributed

(42a)

(42b)

(42c)

(43)

<sup>2</sup> (44)

**u**˙ 1

**X**<sup>1</sup> **X**<sup>2</sup> **O 1**

 **1 O X**<sup>1</sup> **X**<sup>2</sup>

 **1 O X**<sup>1</sup> **X**<sup>2</sup>

body, see eq.(35), is referred to its center of mass with displacement array **u**<sup>1</sup>

and right, be the mass dyads of the two half-parts of the joint, defined as

 *mL***1 O O J***<sup>L</sup>*

**M***<sup>L</sup>* =

*TJ* <sup>=</sup> <sup>1</sup> 2 (**u**˙ <sup>2</sup>

and of the time-derivatives of the nodal arrays ˙**u**<sup>2</sup>

written with perfect analogy to their counterparts **K***a*, **K***<sup>b</sup>* and **K***c*:

*<sup>T</sup>*

*<sup>T</sup>*

*<sup>T</sup>*

**u**˙ *<sup>b</sup>* and ˜

<sup>1</sup> and **<sup>u</sup>**<sup>2</sup>

**M**1,1 <sup>2</sup> **<sup>M</sup>**1,2 2

**M**2,1 <sup>2</sup> **<sup>M</sup>**2,2 2

**M**1,1 <sup>1</sup> **<sup>M</sup>**1,2 1

**M**2,1 <sup>1</sup> **<sup>M</sup>**2,2 1

**M**1,1 <sup>1</sup> **<sup>M</sup>**1,2 1

**M**2,1 <sup>1</sup> **<sup>M</sup>**2,2 1

<sup>2</sup> <sup>=</sup> **<sup>X</sup>**1**u**˙ <sup>1</sup>

of eq.(13), we write

**M***<sup>a</sup>* = **X**<sup>1</sup> **X**<sup>2</sup> **O 1**

**M***<sup>b</sup>* =

**M***<sup>c</sup>* =

time-derivatives ˜

displacement arrays **u**<sup>1</sup>

sought for the cases b) and c).

**5.3 Joint's mass and inertia**

 **1 O X**<sup>1</sup> **X**<sup>2</sup>

 **1 O X**<sup>1</sup> **X**<sup>2</sup>

**u**˙ *<sup>a</sup>*, ˜

#### **6. Linearized elastodynamics equations**

In this section we will derive the generalized inertia and stiffness matrices necessary to write the linearized elastodynamics equations. As described in the two previous sections stiffness and inertia matrices of rigid bodies and flexible beams are referred to the independent nodal displacements of the system. Now, one might ask how to combine different matrices to build the global stiffness and inertia matrices. In order to solve this issue let us introduce two Boolean matrices **B**<sup>1</sup> and **B**<sup>2</sup> able to map a 6-dimensional displacement array **u***<sup>i</sup>* or a 12-dimensional displacement array ˜**u** in terms of a 6*n*-dimensional global array **q** = **u**<sup>1</sup> **u**<sup>2</sup> ... **u***<sup>n</sup>* containing all the independent displacement arrays of the robot. It is important that the reader would define the order in which every single array **u***<sup>i</sup>* appears inside **q**. The 6 × 6*n* matrix **B**<sup>1</sup> and the 12 × 6*n* matrix **B**<sup>2</sup> are defined as:

$$\mathbf{u}\_{i} = \mathbf{B}\_{1}(i)\mathbf{q}\_{\prime} \quad \tilde{\mathbf{u}} \equiv \begin{bmatrix} \mathbf{u}\_{i} \\ \mathbf{u}\_{j} \end{bmatrix} = \mathbf{B}\_{2}(i,j)\mathbf{q} \tag{47}$$

$$\begin{aligned} \mathbf{B}\_1(i) &= \begin{bmatrix} \mathbf{O}\_6 \ \mathbf{O}\_6 \ \dots \ \mathbf{1}\_6(i) \ \dots \ \mathbf{O}\_6 \end{bmatrix}, \quad \mathbf{B}\_2(i,j) = \begin{bmatrix} \mathbf{O}\_6 \ \mathbf{O}\_6 \ \dots \ \mathbf{1}\_6(i) \ \dots \ \mathbf{O}\_6 \\ \mathbf{O}\_6 \ \mathbf{1}\_6(j) \ \dots \ \mathbf{O}\_6 \ \dots \ \mathbf{O}\_6 \end{bmatrix} \tag{48}$$

The above expressions allow us to convert each nodal displacement array in term of **q**. As instance, let us take into exam the strain energy *Vb* of eq.(27), it simply turns into *Vb* <sup>=</sup> <sup>1</sup> <sup>2</sup>**q***T***B**2(*i*, *<sup>j</sup>*)*T***K***b***B**2(*i*, *<sup>j</sup>*)**q**, where *<sup>i</sup>* and *<sup>j</sup>* indicate the position indices of **<sup>u</sup>**<sup>1</sup> <sup>1</sup> and **<sup>u</sup>**<sup>2</sup> <sup>2</sup> inside **q**. From *Vb* we can extract a new 6*n* × 6*n* matrix **K***b*, that is a stiffness matrix expanded to the final dimension of the problem, defined as

$$\overline{\mathbf{K}}\_b = \mathbf{B}\_2(i, j)^T \mathbf{K}\_b \mathbf{B}\_2(i, j) \tag{49}$$

Likewise, the generic expanded 6*n* × 6*n* mass matrix **M***<sup>i</sup>* of **M***i*, see eq.(35), can be written as **<sup>M</sup>***<sup>i</sup>* <sup>=</sup> **<sup>B</sup>**1(*i*)*<sup>T</sup>***<sup>M</sup>***i***B**1(*i*), where *<sup>i</sup>* is the position index of **<sup>u</sup>***i*, i.e. the displacement array of the *i*th-rigid body's center of mass, inside **q**. By the same strategy, it is possible to expand every stiffness and inertia matrix. This operation is essential as, *only when referred to the same set*

and the dynamics of a mechanical system, can only amplify or reduce, without creating or

On the Stiffness Analysis and Elastodynamics of Parallel Kinematic Machines 101

Different considerations can be made for the boundary surface of the workspace. For the latter case, a robot undergoes poses in which one or more legs are completely extended, thus the robot loses dofs associated to those directions normal to the tangent plane at the boundary surface of the workspace. If the displacement of the MP for the corresponding eigenmode, when evaluated at a pose lying on the boundary surface, is normal to the said tangent plane the robot is virtually locked, it means that its stiffness along that direction is high, thus reflecting into an increment of frequency. On the contrary, if the displacement of the MP lies onto the tangent plane its in-plane stiffness, and its natural frequency as well, in theory goes to zero. All remaining cases are included between zero and a maximum value depending on the projection of displacements on the tangent plane. This important conclusion can be summarized by the following statement: *the value of a given natural frequency of a PKM on the boundary surface of its workspace can be interpreted through the projection of the moving platform's displacements, for the corresponding eigenmode, on the tangent plane to the said surface*. The previous statement is strict only when a three dimensional workspace, e.g. the constant orientation workspace, can be considered. In other cases, one should recur

to concepts of projection based on twist theory, beyond the scopes of this chapter.

In this section we apply the method to a complex PKM with 4-dofs. The robot, shown in Fig. 9, is composed of four legs connecting the base frame to a moving platform. Due to the particular kind of constraints generated by the limbs, the moving platform can accomplish three translations and a rotation around the vertical axis. This motion, even referred to as Schönflies motion, is typical of SCARA robots and it is largely used in industry for pick and place applications. Each leg is composed of a distal link connected to the base frame by means of an actuated revolute joint and to a parallelogram linkage, i.e. a Π-joint, by means of a revolute joint. The parallelogram, in turn, is linked to a proximal link by another revolute. Finally, the proximal link is coupled to the MP by a vertical axis revolute joint. All revolute joints, apart from the one fixed to the inertial frame, are passive. We do not further go inside the analysis of the kinematics, citing (24) for more explanations and for any detail pertaining

In order to apply the proposed method each leg is decomposed into three modules. The first module includes the BP and the proximal link, the second the parallelogram linkage, while the last includes the distal link and the MP. We recur to some simplifications that do not change the meaning of the treatment. In general, the two bases of a PKM are one order stiffer than the links composing the structure and can be modeled as rigid without loss of accuracy. Even the short input and output links of a Π-joint can be considered rigid when compared to the slender coupler links. As verified by FEM, the said assumptions are good approximations that simplify the final model introducing rigid-flexible and flexible-rigid connections inside

Let us analyze the first module shown in Fig. 10. The distal link is split into three flexible bodies, the choice of discretization being absolutely arbitrary. In this simple case the end-bodies *F*<sup>1</sup> and *F*<sup>3</sup> of the link are coupled to two rigid bodies, i.e., the base *BP* and the input link *I* of the Π-joint, by means of two revolute joints of axes parallel to the unit vector

nullifying, the effects of the kinematics.

**7. Case study**

the robot's geometry.

and between each module.

**7.1 Application of the method to a generic leg**

#### Fig. 9. CAD model of the robot.

*of nodal displacement coordinates*, these matrices can be summed and combined to obtain the global generalized matrices **K***PKM* and **M***PKM* of the robot. The way to use **B**<sup>1</sup> and **B**<sup>2</sup> will be shown in detail in the case study of the next section.

Let us introduce the 6*n*-dimensional array **f** of generalized nodal wrenches, i.e. nodal forces and torques, then, the system of linearized elastodynamics equations becomes

$$\mathbf{M}\_{PKM}\ddot{\mathbf{q}} + \mathbf{K}\_{PKM}\mathbf{q} = \mathbf{f} \tag{50}$$

The previous system may be used to solve statics around a starting posture of the robot. The homogeneous part of eq.(50), i.e. the left side, may be used to find eigenfrequencies and eigenmodes of the robot, or the zero-input response. As already pointed out, the natural frequencies of the robot change with regard to the pose that the robot is attaining. To determine how stiffness or natural frequencies vary inside the workspace local and global performance indices may be introduced to investigate elastodynamic behavior of PKMs. As instance, let us consider the first natural frequency *f*<sup>1</sup> mapped inside a robot's workspace. The latter can be analyzed to differentiate areas with high range of frequency from areas in which the elastodynamic performances worsen. The multidimensional integral Ω<sup>1</sup> of *f*1, extended to the whole workspace W or part of it, can be a good global index to show how global performances increase or decrease changing some geometrical, structural or inertial parameter, i.e.

$$
\Omega\_1 = \frac{\int\_W f\_1 d\mathcal{W}}{\mathcal{W}} \approx \frac{\sum\_{i=1}^{n\_w} \mathcal{W}\_i f\_{1i}}{\sum\_{i=1}^{n\_w} \mathcal{W}\_i} \tag{51}
$$

where, in the approximated discrete formula of the right side, *W* is divided into *nw* hypercubes *Wi*, while the frequency *f*1*<sup>i</sup>* is calculated at the center point of *Wi*.

We conclude this section considering an useful application of the method to find singularities. It is well known that as a robot approaches to a singularity configuration it gains mobility. As a consequence, its generalized stiffness matrix becomes singular and its first eigenfrequency goes to zero. Analyzing the variation of *f*<sup>1</sup> inside the workspace may be useful to avoid zones close to singularity or with low values of frequency. We stress that singularities directly come from kinematics: stiffness and inertial parameters, respectively influencing the elasticity

and the dynamics of a mechanical system, can only amplify or reduce, without creating or nullifying, the effects of the kinematics.

Different considerations can be made for the boundary surface of the workspace. For the latter case, a robot undergoes poses in which one or more legs are completely extended, thus the robot loses dofs associated to those directions normal to the tangent plane at the boundary surface of the workspace. If the displacement of the MP for the corresponding eigenmode, when evaluated at a pose lying on the boundary surface, is normal to the said tangent plane the robot is virtually locked, it means that its stiffness along that direction is high, thus reflecting into an increment of frequency. On the contrary, if the displacement of the MP lies onto the tangent plane its in-plane stiffness, and its natural frequency as well, in theory goes to zero. All remaining cases are included between zero and a maximum value depending on the projection of displacements on the tangent plane. This important conclusion can be summarized by the following statement: *the value of a given natural frequency of a PKM on the boundary surface of its workspace can be interpreted through the projection of the moving platform's displacements, for the corresponding eigenmode, on the tangent plane to the said surface*. The previous statement is strict only when a three dimensional workspace, e.g. the constant orientation workspace, can be considered. In other cases, one should recur to concepts of projection based on twist theory, beyond the scopes of this chapter.

## **7. Case study**

16 Will-be-set-by-IN-TECH

distal link

shown in detail in the case study of the next section.

Fig. 9. CAD model of the robot.

extended to the whole workspace

parameter, i.e.

MP

proximal link

**M***PKM***q**¨ + **K***PKM***q** = **f** (50)

or part of it, can be a good global index to show how

(51)

Parallelogram

BP

*of nodal displacement coordinates*, these matrices can be summed and combined to obtain the global generalized matrices **K***PKM* and **M***PKM* of the robot. The way to use **B**<sup>1</sup> and **B**<sup>2</sup> will be

Let us introduce the 6*n*-dimensional array **f** of generalized nodal wrenches, i.e. nodal forces

The previous system may be used to solve statics around a starting posture of the robot. The homogeneous part of eq.(50), i.e. the left side, may be used to find eigenfrequencies and eigenmodes of the robot, or the zero-input response. As already pointed out, the natural frequencies of the robot change with regard to the pose that the robot is attaining. To determine how stiffness or natural frequencies vary inside the workspace local and global performance indices may be introduced to investigate elastodynamic behavior of PKMs. As instance, let us consider the first natural frequency *f*<sup>1</sup> mapped inside a robot's workspace. The latter can be analyzed to differentiate areas with high range of frequency from areas in which the elastodynamic performances worsen. The multidimensional integral Ω<sup>1</sup> of *f*1,

global performances increase or decrease changing some geometrical, structural or inertial

*<sup>W</sup>* <sup>≈</sup> <sup>∑</sup>*nw*

where, in the approximated discrete formula of the right side, *W* is divided into *nw* hypercubes

We conclude this section considering an useful application of the method to find singularities. It is well known that as a robot approaches to a singularity configuration it gains mobility. As a consequence, its generalized stiffness matrix becomes singular and its first eigenfrequency goes to zero. Analyzing the variation of *f*<sup>1</sup> inside the workspace may be useful to avoid zones close to singularity or with low values of frequency. We stress that singularities directly come from kinematics: stiffness and inertial parameters, respectively influencing the elasticity

*<sup>i</sup>*=<sup>1</sup> *Wi f*1*<sup>i</sup>* ∑*nw <sup>i</sup>*=<sup>1</sup> *Wi*

and torques, then, the system of linearized elastodynamics equations becomes

W

 *<sup>W</sup> f*1*dW*

Ω<sup>1</sup> =

*Wi*, while the frequency *f*1*<sup>i</sup>* is calculated at the center point of *Wi*.

In this section we apply the method to a complex PKM with 4-dofs. The robot, shown in Fig. 9, is composed of four legs connecting the base frame to a moving platform. Due to the particular kind of constraints generated by the limbs, the moving platform can accomplish three translations and a rotation around the vertical axis. This motion, even referred to as Schönflies motion, is typical of SCARA robots and it is largely used in industry for pick and place applications. Each leg is composed of a distal link connected to the base frame by means of an actuated revolute joint and to a parallelogram linkage, i.e. a Π-joint, by means of a revolute joint. The parallelogram, in turn, is linked to a proximal link by another revolute. Finally, the proximal link is coupled to the MP by a vertical axis revolute joint. All revolute joints, apart from the one fixed to the inertial frame, are passive. We do not further go inside the analysis of the kinematics, citing (24) for more explanations and for any detail pertaining the robot's geometry.

#### **7.1 Application of the method to a generic leg**

In order to apply the proposed method each leg is decomposed into three modules. The first module includes the BP and the proximal link, the second the parallelogram linkage, while the last includes the distal link and the MP. We recur to some simplifications that do not change the meaning of the treatment. In general, the two bases of a PKM are one order stiffer than the links composing the structure and can be modeled as rigid without loss of accuracy. Even the short input and output links of a Π-joint can be considered rigid when compared to the slender coupler links. As verified by FEM, the said assumptions are good approximations that simplify the final model introducing rigid-flexible and flexible-rigid connections inside and between each module.

Let us analyze the first module shown in Fig. 10. The distal link is split into three flexible bodies, the choice of discretization being absolutely arbitrary. In this simple case the end-bodies *F*<sup>1</sup> and *F*<sup>3</sup> of the link are coupled to two rigid bodies, i.e., the base *BP* and the input link *I* of the Π-joint, by means of two revolute joints of axes parallel to the unit vector

*I*

*I*

parallelogram linkage

**y***i* **z***i*

*O*

link-moving platform

*<sup>T</sup>***K***M*1**q***M*<sup>1</sup>

independent nodal-arrays in terms of the global array **q**.

**e**1

can be calculated as

*Vm*<sup>1</sup> <sup>=</sup> <sup>1</sup> 2 **q***M*<sup>1</sup> *Fi*

**<sup>x</sup>***<sup>i</sup>* **<sup>y</sup>***<sup>i</sup>*

**z***i*

**d***i*

*Fi* **x***i*

(a) Drawing of a the second module:

**e**2

(b) Drawing of a the third module: proximal

from the center of mass of *I* to the center of the revolute joint, in this case, is the zero vector. Then, by introducing the 12 × *nqm*<sup>1</sup> binary-entry matrix **B**2(•, •) mapping the independent nodal-arrays in terms of the array **q***M*1, the expression of the generalized stiffness matrix **K***M*<sup>1</sup>

Here, the matrix **K***M*<sup>1</sup> is expressed in terms of **q***M*1. The reader must notice that to obtain **K***PKM* each matrix has to be expressed in terms of a final *nq*-dimensional array **q** including *all* the independent nodal displacements. In order to define **q**, as the modules of a leg are in series, the last node of a module coincides with the first one of the next module and it must be denoted by a unique nodal displacement array. Extension to the final dimension of **q** is readily obtained by substituting into eq.(54) a new 12 × *nq* binary-entry matrix **B**2(•, •) mapping the

The second module is shown in Fig. 11(a). In this case two rigid-flexible and two flexible-rigid connections must be used to describe the couplings between flexible couplers and input-output rigid bodies *I* and *O*. The joint array **h***R*(**e**) remains the same for all the four revolute joints as a matter of fact of the parallelogram architecture. The third module, shown in Fig. 11(b), includes one rigid-flexible and one flexible-rigid connection, respectively, between the distal link and the output link *O* of the Π-joint and the distal link and the MP. As displayed in the same figure, the two revolute joints to be used to find the joint arrays **h***<sup>R</sup>* have axes of unit vectors **e**<sup>1</sup> and **e**2. We do not describe the calculations to find the generalized

**K***M*<sup>1</sup> = **B**2(1, 2)*T***K**1**B**2(1, 2) + **B**2(2, 3)*T***K**2**B**2(2, 3) + **B**2(3, 4)*T***K***c*3**B**2(3, 4) (54)

**d***i*

**z***<sup>M</sup>*

*MP P*

**x***<sup>M</sup>*

**y***<sup>M</sup>*

**e**

On the Stiffness Analysis and Elastodynamics of Parallel Kinematic Machines 103

*O*

Fig. 10. Drawing of a the first module: Base platform-distal link

**e**. The remaining flexible body *F*<sup>2</sup> is internal to the proximal link and coupled to the others bodies by means of fixed connections. Following what explained in the previous sections, the strain energy *Vm*<sup>1</sup> of the first module is the sum of three components, i.e. the strain energies *VFi* of the three beams in which is decomposed, namely

$$V\_{m1} = \sum\_{i=1}^{3} V\_{F\_i} = \sum\_{i=1}^{3} \frac{1}{2} \tilde{\mathbf{u}}\_i^T \mathbf{K}\_i \tilde{\mathbf{u}}\_i \tag{52}$$

with ˜**u**<sup>1</sup> = **u**1 1 *T* **u**2 1 *<sup>T</sup> <sup>T</sup>* , ˜**u**<sup>2</sup> = **u**2 2 *T* **u**3 2 *<sup>T</sup> <sup>T</sup>* and ˜**u**<sup>3</sup> = **u**3 3 *T* **u**4 3 *<sup>T</sup> <sup>T</sup>* . The first actuated revolute joint is considered locked to perform the elastodynamics analysis, thereby, **u**<sup>1</sup> <sup>1</sup> <sup>≡</sup> **<sup>u</sup>**<sup>1</sup> *BP* <sup>≡</sup> **<sup>u</sup>**1. In this way the generalized stiffness matrix of the robot is not singular as rigid body motions of the robot are deleted and the analysis is performed at a frozen configuration. Even for the other three fixed connections at points 2 and 3, similar equations stand: **u**<sup>2</sup> <sup>1</sup> <sup>≡</sup> **<sup>u</sup>**<sup>2</sup> <sup>2</sup> <sup>≡</sup> **<sup>u</sup>**<sup>2</sup> and **u**3 <sup>2</sup> <sup>≡</sup> **<sup>u</sup>**<sup>3</sup> <sup>3</sup> <sup>≡</sup> **<sup>u</sup>**3. Finally, for the displacement array of the input link *<sup>R</sup>*, we can write: **<sup>u</sup>***<sup>R</sup>* <sup>≡</sup> **<sup>u</sup>**4. By substituting into the previous expressions, we derive ˜**u**<sup>1</sup> = **<sup>u</sup>**1*<sup>T</sup>* **<sup>u</sup>**2*<sup>T</sup> <sup>T</sup>* , ˜**u**<sup>2</sup> = **<sup>u</sup>**2*<sup>T</sup>* **<sup>u</sup>**3*<sup>T</sup> <sup>T</sup>* and ˜**u**<sup>3</sup> = **<sup>u</sup>**3*<sup>T</sup>* **u**4 3 *<sup>T</sup> <sup>T</sup>* . The local stiffness matrix **K***<sup>i</sup>* of the *i*th-body *Fi* is expressed into the global reference frame. In Figure 10 the first local frame {**x**1, **y**1, **z**1}, with origin at point 1, and the base global reference frame {**x**, **y**, **z**}, with origin at point *O*, are displayed. The set of independent displacements are defined inside the independent displacement array **q***M*<sup>1</sup> of the first module: **<sup>q</sup>***M*<sup>1</sup> = **u**<sup>1</sup> **u**<sup>2</sup> **u**<sup>3</sup> **u**<sup>4</sup> , where **u**<sup>4</sup> is the nodal displacement array of the center of mass of the rigid body *I*. While the first and the second term, i.e. *VF*<sup>1</sup> and *VF*<sup>2</sup> , of the strain energy *Vm*<sup>1</sup> contain only independent displacements, the last term *VF*<sup>3</sup> is function of the dependent nodal array **u**<sup>4</sup> <sup>3</sup>. The flexible body *F*3, indeed, is coupled by a revolute joint to the rigid body *I*, thus, the case c), flexible-rigid, can be applied to express its strain energy only in terms of independent displacements. Following the notation above introduced, we write:

$$V\_{\mathbf{F}\_3} = \frac{1}{2} \tilde{\mathbf{u}}\_3^T \mathbf{K}\_3 \tilde{\mathbf{u}}\_3 \equiv \frac{1}{2} \tilde{\mathbf{u}}\_{c3}^T \mathbf{K}\_{c3} \tilde{\mathbf{u}}\_{c3} \tag{53}$$

where ˜**u***c*<sup>3</sup> ≡ **<sup>u</sup>**3*<sup>T</sup>* **<sup>u</sup>**4*<sup>T</sup> <sup>T</sup>* and **K***c*<sup>3</sup> has been defined in eq.(31). Notice that in order to find **K***c*<sup>3</sup> we have to use the 6-dimensional joint array **h***R*(**e**); besides, the position vector **d** going 18 Will-be-set-by-IN-TECH

**x**

**e**. The remaining flexible body *F*<sup>2</sup> is internal to the proximal link and coupled to the others bodies by means of fixed connections. Following what explained in the previous sections, the strain energy *Vm*<sup>1</sup> of the first module is the sum of three components, i.e. the strain energies

*VFi* =

and ˜**u**<sup>3</sup> =

In this way the generalized stiffness matrix of the robot is not singular as rigid body motions of the robot are deleted and the analysis is performed at a frozen configuration. Even for the

global reference frame. In Figure 10 the first local frame {**x**1, **y**1, **z**1}, with origin at point 1, and the base global reference frame {**x**, **y**, **z**}, with origin at point *O*, are displayed. The set of independent displacements are defined inside the independent displacement array **q***M*<sup>1</sup>

center of mass of the rigid body *I*. While the first and the second term, i.e. *VF*<sup>1</sup> and *VF*<sup>2</sup> , of the strain energy *Vm*<sup>1</sup> contain only independent displacements, the last term *VF*<sup>3</sup> is function of the

rigid body *I*, thus, the case c), flexible-rigid, can be applied to express its strain energy only in terms of independent displacements. Following the notation above introduced, we write:

<sup>3</sup> **<sup>K</sup>**3**u**˜ <sup>3</sup> <sup>≡</sup> <sup>1</sup>

**K***c*<sup>3</sup> we have to use the 6-dimensional joint array **h***R*(**e**); besides, the position vector **d** going

2 **u**˜ *T*

<sup>3</sup> <sup>≡</sup> **<sup>u</sup>**3. Finally, for the displacement array of the input link *<sup>R</sup>*, we can write: **<sup>u</sup>***<sup>R</sup>* <sup>≡</sup> **<sup>u</sup>**4.

3 ∑ *i*=1

1 2 **u**˜ *T*

> **u**3 3 *T* **u**4 3 *<sup>T</sup> <sup>T</sup>*

> > **<sup>u</sup>**1*<sup>T</sup>*

. The local stiffness matrix **K***<sup>i</sup>* of the *i*th-body *Fi* is expressed into the

<sup>3</sup>. The flexible body *F*3, indeed, is coupled by a revolute joint to the

and **K***c*<sup>3</sup> has been defined in eq.(31). Notice that in order to find

**<sup>u</sup>**2*<sup>T</sup> <sup>T</sup>*

, where **u**<sup>4</sup> is the nodal displacement array of the

*<sup>c</sup>*3**K***c*3**u**˜ *<sup>c</sup>*<sup>3</sup> (53)

1

**<sup>x</sup>**<sup>1</sup> **<sup>y</sup>**<sup>1</sup> **z**1

*F*1

2

*F*2

3

*I*

*<sup>i</sup>* **K***i***u**˜*<sup>i</sup>* (52)

. The first actuated revolute

<sup>1</sup> <sup>≡</sup> **<sup>u</sup>**<sup>2</sup>

, ˜**u**<sup>2</sup> = **<sup>u</sup>**2*<sup>T</sup>*

<sup>1</sup> <sup>≡</sup> **<sup>u</sup>**<sup>1</sup>

*BP* <sup>≡</sup> **<sup>u</sup>**1.

<sup>2</sup> <sup>≡</sup> **<sup>u</sup>**<sup>2</sup> and

**<sup>u</sup>**3*<sup>T</sup> <sup>T</sup>*

*F*3

4

*BP*

3 ∑ *i*=1

joint is considered locked to perform the elastodynamics analysis, thereby, **u**<sup>1</sup>

other three fixed connections at points 2 and 3, similar equations stand: **u**<sup>2</sup>

**u**<sup>1</sup> **u**<sup>2</sup> **u**<sup>3</sup> **u**<sup>4</sup>

*VF*<sup>3</sup> <sup>=</sup> <sup>1</sup> 2 **u**˜ *T*

*O*

**y z**

Fig. 10. Drawing of a the first module: Base platform-distal link

*Vm*<sup>1</sup> =

*VFi* of the three beams in which is decomposed, namely

, ˜**u**<sup>2</sup> = **u**2 2 *T* **u**3 2 *<sup>T</sup> <sup>T</sup>*

By substituting into the previous expressions, we derive ˜**u**<sup>1</sup> =

with ˜**u**<sup>1</sup> =

and ˜**u**<sup>3</sup> =

where ˜**u***c*<sup>3</sup> ≡

**u**3 <sup>2</sup> <sup>≡</sup> **<sup>u</sup>**<sup>3</sup>  **u**1 1 *T* **u**2 1 *<sup>T</sup> <sup>T</sup>*

 **<sup>u</sup>**3*<sup>T</sup>* **u**4 3 *<sup>T</sup> <sup>T</sup>*

of the first module: **<sup>q</sup>***M*<sup>1</sup> =

dependent nodal array **u**<sup>4</sup>

 **<sup>u</sup>**3*<sup>T</sup>* **<sup>u</sup>**4*<sup>T</sup> <sup>T</sup>*

link-moving platform

from the center of mass of *I* to the center of the revolute joint, in this case, is the zero vector. Then, by introducing the 12 × *nqm*<sup>1</sup> binary-entry matrix **B**2(•, •) mapping the independent nodal-arrays in terms of the array **q***M*1, the expression of the generalized stiffness matrix **K***M*<sup>1</sup> can be calculated as

$$\begin{aligned} V\_{m1} &= \frac{1}{2} \mathbf{q}\_{M1}{}^T \mathbf{K}\_{M1} \mathbf{q}\_{M1} \\ \mathbf{K}\_{M1} &= \mathbf{B}\_2(1,2) {}^T \mathbf{K}\_1 \mathbf{B}\_2(1,2) + \mathbf{B}\_2(2,3) {}^T \mathbf{K}\_2 \mathbf{B}\_2(2,3) + \mathbf{B}\_2(3,4) {}^T \mathbf{K}\_{c3} \mathbf{B}\_2(3,4) \end{aligned} \tag{54}$$

Here, the matrix **K***M*<sup>1</sup> is expressed in terms of **q***M*1. The reader must notice that to obtain **K***PKM* each matrix has to be expressed in terms of a final *nq*-dimensional array **q** including *all* the independent nodal displacements. In order to define **q**, as the modules of a leg are in series, the last node of a module coincides with the first one of the next module and it must be denoted by a unique nodal displacement array. Extension to the final dimension of **q** is readily obtained by substituting into eq.(54) a new 12 × *nq* binary-entry matrix **B**2(•, •) mapping the independent nodal-arrays in terms of the global array **q**.

The second module is shown in Fig. 11(a). In this case two rigid-flexible and two flexible-rigid connections must be used to describe the couplings between flexible couplers and input-output rigid bodies *I* and *O*. The joint array **h***R*(**e**) remains the same for all the four revolute joints as a matter of fact of the parallelogram architecture. The third module, shown in Fig. 11(b), includes one rigid-flexible and one flexible-rigid connection, respectively, between the distal link and the output link *O* of the Π-joint and the distal link and the MP. As displayed in the same figure, the two revolute joints to be used to find the joint arrays **h***<sup>R</sup>* have axes of unit vectors **e**<sup>1</sup> and **e**2. We do not describe the calculations to find the generalized

(a) First mode

On the Stiffness Analysis and Elastodynamics of Parallel Kinematic Machines 105

(b) Second mode

(c) Third mode

(d) Fourth mode

As in the case of statics, the elastodynamics model of the 3T1R robot is used to compare natural modes and frequencies to FEA software. We report only the comparison to Ansys©, as more indicative of a feasible application of the method to a real case. Figure 12 shows the first four natural modes obtained when the robot is attaining its home pose. In Table 2 the first ten natural frequency of the robot, at the same pose, are reported. Results show good agreement

Fig. 12. Modes comparison: Algorithm *vs*. Ansys©.

**7.2.2 Natural modes and frequencies**

Fig. 11. Static deformation of the robot.


Table 1. Translational displacements of the end-effector's center node.

stiffness matrices of the two modules as can be readily obtained following procedures similar to the case of the first module.

#### **7.2 Comparison to FEA software and validation**

In this subsection the proposed elastodynamics model is compared to FEA results for validation. Two FEA models, with increasing complexity, are implemented in order to establish the nearness of the examined case to the real case. The first model, developed by the commercial software Marc©, considers only 3D-beams and rigid bodies, while joints are modeled by means of relative degrees of freedom between common nodes belonging to two coupled bodies. This model perfectly fits the simplifications of the method in exam. The second model, developed by the commercial software Ansys©, describes a robot with a complex structure closer to reality. Links are solids while joints have finite burdens and provide their function by means of the coupling of surfaces and screws.

#### **7.2.1 Statics**

We first compared the static deformation of the robot when an external force of 1000 N, along the x-direction, is applied on the end-effector, when the latter is positioned at its home pose with angle of rotation of MP equal to zero. Figure 11 shows the displacements of our model when compared to Ansys© model, while in Tab. 1 the translational displacements of the end-effector center node are reported.

It can be observed how the first MARC© model perfectly fit to the results of our method. The relative error grows, but still remaining limited, when displacements are compared to Ansys© model. The reason is well understood as it comes from the use of solid bodies and joints to simulate the robot's structure.

20 Will-be-set-by-IN-TECH

(c) Algorithm model (d) Ansys© model

*Disp. Algorithm Marc© Ansys© err% Marc© err% Ansys© ux* 0.068066 0.068065 0.067436 −0.0008815 −0.933626 *uy* −0.030545 −0.030545 −0.029890 0.003274 −2.191368 *uz* 0 0 0 0 0

stiffness matrices of the two modules as can be readily obtained following procedures similar

In this subsection the proposed elastodynamics model is compared to FEA results for validation. Two FEA models, with increasing complexity, are implemented in order to establish the nearness of the examined case to the real case. The first model, developed by the commercial software Marc©, considers only 3D-beams and rigid bodies, while joints are modeled by means of relative degrees of freedom between common nodes belonging to two coupled bodies. This model perfectly fits the simplifications of the method in exam. The second model, developed by the commercial software Ansys©, describes a robot with a complex structure closer to reality. Links are solids while joints have finite burdens and

We first compared the static deformation of the robot when an external force of 1000 N, along the x-direction, is applied on the end-effector, when the latter is positioned at its home pose with angle of rotation of MP equal to zero. Figure 11 shows the displacements of our model when compared to Ansys© model, while in Tab. 1 the translational displacements of the

It can be observed how the first MARC© model perfectly fit to the results of our method. The relative error grows, but still remaining limited, when displacements are compared to Ansys© model. The reason is well understood as it comes from the use of solid bodies and joints to

Table 1. Translational displacements of the end-effector's center node.

provide their function by means of the coupling of surfaces and screws.

u

Fig. 11. Static deformation of the robot.

to the case of the first module.

end-effector center node are reported.

simulate the robot's structure.

**7.2.1 Statics**

**7.2 Comparison to FEA software and validation**

(d) Fourth mode

Fig. 12. Modes comparison: Algorithm *vs*. Ansys©.

#### **7.2.2 Natural modes and frequencies**

As in the case of statics, the elastodynamics model of the 3T1R robot is used to compare natural modes and frequencies to FEA software. We report only the comparison to Ansys©, as more indicative of a feasible application of the method to a real case. Figure 12 shows the first four natural modes obtained when the robot is attaining its home pose. In Table 2 the first ten natural frequency of the robot, at the same pose, are reported. Results show good agreement

**8. Conclusions**

**9. References**

This chapter has discussed a method, based on the Matrix Structural Analysis, to study the linearized elastodynamics of PKMs. Base and moving platforms are considered rigid, while links can be modeled as rigid or flexible parts, the latter being decomposed into two or more flexible bodies. Here, we used 3D Euler beams but the method can be extended to superelements with two end-nodes. Joints are directly included, without recurring to Lagrange multipliers, by means of kinematic constraints between nodal displacement arrays. Three cases have been taken into account to model the rigid-flexible, flexible-flexible and flexible-rigid coupling of bodies by means joints. Each case yields equations, linking dependent, independent nodal displacement arrays and joint displacements as well, to be used to find generalized stiffness and inertia matrices. The latter are then combined as elementary blocks to find the global matrices of the whole system. Some useful extension to compliant mechanism has been introduced, while two strategy, the lumped and the distributed one, have been explained to include mass/inertia into the model. Feasible applications of the method pertain: the study of natural frequencies inside the robot's workspace by means of local and global indices, the singularity finding, the optimization of

On the Stiffness Analysis and Elastodynamics of Parallel Kinematic Machines 107

elastodynamic performances varying geometric, structural or inertial parameters.

and time-consuming FEM meshing routines.

*Automation*, May 11-15, Washington, US.

Machines, McGill University, Montreal, Canada.

Heinemann, London.

Book Company.

New York.

Finally, the method has been applied to an articulated four-dofs PKMs with Schönflies motions. A modular approach is used to split each of the four legs into three modules. Results, compared to commercial software, revealed good accuracy in determining natural frequency range, while drastically reducing the time of computation avoiding the annoying

[1] Tyapin, I., Hovland, G. & Brogårdh, T. (2008). Kinematic and Elastodynamic Design Optimisation of the 3-DOF Gantry-Tau Parallel Kinematic Manipulator, In: *Proceedings of the Second International Workshop on Fundamental Issues and Future Research Directions*

[3] Martin, H. C. (1966). *Introduction to matrix methods of structural analysis*, McGraw-Hill

[4] Wang, C.K. (1966). *Matrix methods of structural analysis*, International Textbook Company. [5] Przemieniecki, J. S. (1985). *Theory of Matrix Structural Analysis*, Dover Publications, Inc,

[6] Huang, T., Zhao, X. & Withehouse, D.J. (2002). Stiffness estimation of a tripod-based parallel kinematic machine, *IEEE Trans. on Robotics and Automation*, Vol. 18(1). [7] Li, Y.W., Wang, J.S. & Wang, L.P. (2002). Stiffness analysis of a Stewart platform-based parallel kinematic machine, In: *Proceedings of IEEE ICRA: Int. Conf. On Robotics and*

[8] Clinton, C.M., Zhang, G. & Wavering, A.J. (1997). Stiffness modeling of a Stewart-platform-based milling machine, In: *Trans. of the North America Manufacturing*

*Research Institution of SME*, May 20-23, Vol. XXV, pp 335-340, Lincoln, NB, US. [9] Al Bassit, L., Angeles, J., Al-Wydyan, K. & Morozov., A. (2002). The elastodynamics of a Schönflies -Motion generator, *Technical Report TR-CIM-02-06*, Centre for Intelligent

*for Parallel Mechanisms and Manipulators*, September 21-22, Montpellier, France. [2] Zienkiewicz, O. C. & Taylor, R. L. (2000). *Solid mechanics - Volume 2*, Butterworth


Table 2. The first ten natural frequencies of the robot at the home pose.

with Ansys© model. Some discrepancies still occur due to the same reasons explained for the static case and to the use of flexible bodies to simulate all links and platforms in Ansys©. This choice has been taken in order to avoid asymmetric contacts between rigid and flexible solids leading to convergence mistakes. In turn, we used a stiffer material to approximate rigid behavior of the MP and rigid links.

#### **7.2.3 Distribution of frequencies inside the workspace**

In order to provide a further useful extension of the proposed method, the first natural frequency is calculated and plotted inside the constant orientation workspace of the robot (25). We have chosen elementary cubes of side 5 cm to discretize the workspace of the robot while the angle of rotation of the MP is *θ<sup>z</sup>* = 0. It can be observed that two privileged diagonals divide the workspace into four symmetric areas. These directions coincide with the two diagonals of the squared BP of Fig. 9, meaning that, at the home pose, geometric symmetry reflects itself into the elastic behavior of the robot. The boundary of the constant orientation workspace shows areas with high range of frequency along with other areas in which the first natural frequency reaches values close to zero (Hz). As already outlined in the text, this behavior is due to the MP's displacement (deformation) in correspondence of the first eigenmode: if at a certain pose, in which the analysis is performed, the displacement is along a *doc* of the MP, the ensuing frequency will be high; conversely, if it is along a *dof* of the MP, the frequency will come near zero.

Fig. 13. Distribution of the first natural frequency inside the robot's constant orientation workspace.

#### **8. Conclusions**

22 Will-be-set-by-IN-TECH

**Freq. Algorithm [Hz] Ansys© [Hz] err%** 2.920 3.019 3.28% 7.062 7.204 1.97% 11.275 11.487 1.85% 20.834 21.554 3.34% 25.480 25.722 0.94% 25.545 25.804 1.00% 26.028 26.181 0.58% 28.306 28.830 1.82% 30.368 31.295 2.96% 31.060 32.286 3.80%

with Ansys© model. Some discrepancies still occur due to the same reasons explained for the static case and to the use of flexible bodies to simulate all links and platforms in Ansys©. This choice has been taken in order to avoid asymmetric contacts between rigid and flexible solids leading to convergence mistakes. In turn, we used a stiffer material to approximate

In order to provide a further useful extension of the proposed method, the first natural frequency is calculated and plotted inside the constant orientation workspace of the robot (25). We have chosen elementary cubes of side 5 cm to discretize the workspace of the robot while the angle of rotation of the MP is *θ<sup>z</sup>* = 0. It can be observed that two privileged diagonals divide the workspace into four symmetric areas. These directions coincide with the two diagonals of the squared BP of Fig. 9, meaning that, at the home pose, geometric symmetry reflects itself into the elastic behavior of the robot. The boundary of the constant orientation workspace shows areas with high range of frequency along with other areas in which the first natural frequency reaches values close to zero (Hz). As already outlined in the text, this behavior is due to the MP's displacement (deformation) in correspondence of the first eigenmode: if at a certain pose, in which the analysis is performed, the displacement is along a *doc* of the MP, the ensuing frequency will be high; conversely, if it is along a *dof* of the

[*Hz*]

Fig. 13. Distribution of the first natural frequency inside the robot's constant orientation

Table 2. The first ten natural frequencies of the robot at the home pose.

rigid behavior of the MP and rigid links.

MP, the frequency will come near zero.

workspace.

**7.2.3 Distribution of frequencies inside the workspace**

This chapter has discussed a method, based on the Matrix Structural Analysis, to study the linearized elastodynamics of PKMs. Base and moving platforms are considered rigid, while links can be modeled as rigid or flexible parts, the latter being decomposed into two or more flexible bodies. Here, we used 3D Euler beams but the method can be extended to superelements with two end-nodes. Joints are directly included, without recurring to Lagrange multipliers, by means of kinematic constraints between nodal displacement arrays. Three cases have been taken into account to model the rigid-flexible, flexible-flexible and flexible-rigid coupling of bodies by means joints. Each case yields equations, linking dependent, independent nodal displacement arrays and joint displacements as well, to be used to find generalized stiffness and inertia matrices. The latter are then combined as elementary blocks to find the global matrices of the whole system. Some useful extension to compliant mechanism has been introduced, while two strategy, the lumped and the distributed one, have been explained to include mass/inertia into the model. Feasible applications of the method pertain: the study of natural frequencies inside the robot's workspace by means of local and global indices, the singularity finding, the optimization of elastodynamic performances varying geometric, structural or inertial parameters.

Finally, the method has been applied to an articulated four-dofs PKMs with Schönflies motions. A modular approach is used to split each of the four legs into three modules. Results, compared to commercial software, revealed good accuracy in determining natural frequency range, while drastically reducing the time of computation avoiding the annoying and time-consuming FEM meshing routines.

#### **9. References**


**6** 

**Parallel, Serial and** 

*1UAE University* 

*3Tawazun Holding 1,3United Arab Emirates* 

*2Japan* 

*2Kitami Institute of Technology* 

**Hybrid Machine Tools and** 

Khalifa H. Harib1, Kamal A.F. Moustafa1, A.M.M. Sharif Ullah2 and Salah Zenieh3

**Robotics Structures: Comparative** 

**Study on Optimum Kinematic Designs** 

After their inception in the past two decades as possible alternatives to conventional Computer Numerical Controlled (CNC) machine tools structures that dominantly adapt serial structures, Parallel Kinematic Machines (PKM) were anticipated to form a basis for a new generation of future machining centers. However this hope quickly faded out as most problems associated with this type of structures still persist and could not be completely solved satisfactorily. This especially becomes more apparent in machining applications where accuracy, rigidity, dexterity and large workspace are important requirements. Although the PKMs possess superior mechanical characteristics to serial structures, particularly in terms of high rigidity, accuracy and dynamic response, however the PKMs have their own drawbacks including singularity problems, inconsistent dexterity, irregular

To alleviate the PKMs' limitations, considerable research efforts were directed to solve these problems. Optimum design methods are among the various methods that are attempted to improve the dexterity as well as to maximize the workspace (Stoughton and Arai, 1993, Huang et al., 2000). Various methods to evaluate the workspace were suggested (Gosselin, 1990; Luh et al., 1996; Conti et al., 1998; Tsai et al., 2006). Workspace optimization is also addressed (Wang and Hsieh, 1998). A new shift in tackling the aforementioned problems came when researchers start to look at hybrid structures, consisting of parallel and serial linkages as a compromise to exploit the advantageous characteristics of the serial and parallel structures. This shift creates new research and development needs and founded

Among the early hybrid kinematic designs, the Tricept was considered as the first commercially successful hybrid machine tools. This hybrid machine which was developed by Neos Robotics, has a three-degrees-of-freedom parallel kinematic structure and a

workspace, and limited range of motion, particularly rotational motion.

**1. Introduction** 

new ideas.

