# Kineto-Elastodynamic Characteristics of the Six-Degree-of-Freedom Parallel Structure Seismic Simulator

Authors: Yongjie Zhao
Journal: Journal of Robotics
 (2011)

---

## Abstract

Based on the kineto-elastodynamic assumptions, the dynamic model of the six-degree-of-freedom parallel structure seismic simulator is developed by virtue of the finite element method and the substructure synthesis technique. The kineto-elastodynamic characteristics represented by the natural frequency, the sensitivity analysis, the energy ratios, and the displacement response of the moving platform are investigated. It is shown that the second-order natural frequency is much higher than the first-order natural frequency, and the first-order natural frequency is sensitive to the radius of the strut and the radius of the lead screw. In order to improve the dynamic characteristic of the manipulator, the mass of the moving platform should be reduced or the stiffness of the strut should be increased especially for the sixth strut. For the investigated trajectory, the displacement response of the moving platform along thex direction is smaller than these displacement responses along the y direction and along the z direction. The angular displacement response of the moving platform rotating about z-axis is slightly larger than those angular displacement responses rotating about the x-axis and about the y-axis.


---

## Body

## 1. Introduction

A seismic simulator is one of the most important equipments in the earthquake resistance testing.Due to the requirement of the large and variable load capability, these kinds of equipments are usually developed with the parallel structure manipulators [1–4]. The parallel manipulator is a closed-loop kinematic chain mechanism whose end effector is linked to the base by several independent kinematic chains [5–7]. For this type of manipulators, there are some potential advantages such as high accuracy, rigidity, and speed. They have been successfully used in the motion simulators, robotic end effectors, and other circumstances like fast pick-and-place operation. Many investigations have been carried out on the parallel manipulators since the concept was introduced. However, there is not many works on the flexible dynamics of the parallel manipulator [8, 9] compared with the vast of papers on the kinematics and rigid dynamics due to the following facts: (i) computational cost; (ii) geometrical complexity; (iii) unidentified mechanics property. For the 6-PSS (prismatic-spherical-spherical joint) flexible parallel manipulator under consideration in this paper, which is developed for the six-degree-of-freedom seismic simulator, the dynamics considering the structure flexibility is fundamental for the modeling, design, and control.The demands of high speed, high load, high precision, or lightweight structure from industry make it necessary to consider the deformation, stiffness, and other dynamic characteristics for the parallel manipulator [10–23]. Mathematical modeling of a general flexible parallel manipulator is a challenging task since there is no availability of closed-form solutions to the inverse kinematic model for the flexible parallel manipulator. The nominal motion of the manipulator involves changing geometries resulting in varying system parameters. The equations of motion are usually configuration dependent and need to be computed at each configuration of the manipulator [10]. The equations of motion of a flexible five-bar manipulator were developed by means of the instantaneous structural approach, and it had been found that the mode shapes and natural frequencies of this particular manipulator are invariant throughout most of the workspace [11]. The design, dynamic modeling, and experiment validation of a three-degree-of-freedom flexible arm were presented in [12] on the assumption that all the arm mass is concentrated at the tip and at the base. So the dynamic of the arm becomes a lumped single mass model instead of the usual distributed mass model. The finite element method and the Euler-Lagrange formulation were used in [13] to model the flexible link of a three-degree-of-freedom parallel manipulator by assuming that the influence of flexible motion on rigid motion is negligible. With the piston being modeled as a mass-spring damper, a set of twelve Lagrange equations for flexible Stewart manipulator was derived by using tensor representation in [14]. The dynamic model of the 3-PRR planar parallel manipulator with flexible links was formulated by using the Lagrange equations of the first type on the assumption that the intermediate links being modeled with pinned-free boundary conditions [15]. The Lagrange finite element formulation was used to derive such a dynamic model for the flexible planar linkage with two translational and one rotational degrees of freedom, and then the dynamic model was applied to the flexible link planar parallel manipulator based on standard Kineto-elastodynamic assumptions [16]. Based on the model, strain rate feedback control using PZT transducers was used to simulate the active control of Kineto-elastodynamic responses. The dynamic finite element analysis of the flexible planar parallel manipulator was presented in [17] including the convergence analysis of the natural frequencies and the mapping of the first-order natural frequency with respect to the robot configuration. It had also been found that the geometric stiffness and the dynamic terms have a negligible effect on the response for this particular manipulator. A substructure modeling procedure was presented to develop the dynamic model for the flexible planar parallel manipulator in [18]. The Craig-Bampton method was used to reduce the model order and assemble the complete dynamic model. On the assumption that the deformations of the intermediate links are small relative to the length of the links, a procedure for the development of structural dynamic model for the 3-PRR flexible parallel manipulator was presented in [19] based on the assumed mode method. Without considering the effect of nominal motion, reference [20] provided the stationary vibration model of the sliding-leg parallel kinematic machine where the links were modeled as finite elements and the joint as virtual spring/dampers. Then, the nonstationary model of the same mechanism was developed with the elastodynamic method [21]. In the researches cited above, there is little investigation on the Kineto-elastodynamic characteristics of the six-degree-of-freedom parallel manipulator while considering the natural frequency, the sensitivity analysis, the energy ratios, and the displacement response.This paper presents the Kineto-elastodynamic modeling and the Kineto-elastodynamic characteristics analysis of the 6-PSS parallel structure seismic simulator. It is organized as follows: in Section2, the description of the seismic simulator and the rigid dynamic equations are presented. Section 3 gives the Kineto-elastodynamic model of the manipulator developed by virtue of the finite element method and the substructure synthesis technique. Section 4 investigates the Kineto-elastodynamic characteristics represented by the natural frequency, the sensitivity analysis, the energy ratios, and the displacement response of the moving platform through simulation. Section 5 gives the conclusions.


## 2. System Description and Rigid Dynamics

### 2.1. Description

The schematic diagram of the 6-dof parallel structure seismic simulator is shown in Figure1. As shown in Figure 1, the parallel manipulator is composed of a moving platform and six sliders. In each kinematic chain, the platform and the slider are connected via spherical ball-bearing joints by a strut of fixed length. Each slider is driven by DC motor via a linear ball screw. The lead screws of B1,B2, and B3 are vertical to the ground.Figure 1
Schematic diagram of the 6-dof parallel structure seismic simulator.For the purpose of analysis, the following coordinate systems are defined. As shown in Figure2, the coordinate system O-xyz  is attached to the fixed base; another moving coordinate frame O′-uvw  is located at the center of mass of the moving platform. The pose of the moving platform can be described by a position vector r and a rotation matrixRo′o. Let the rotation matrix be defined by the roll, pitch, and yaw angles, namely, a rotation of ϕx about the fixed x axis, followed by a rotation of ϕy about the fixed y axis, and a rotation of ϕz about the fixed z axis. Thus, the rotation matrix is(1)Ro′o=Rot(z,ϕz)Rot(y,ϕy)Rot(x,ϕx),
where sϕ denotes the sine of angle ϕ, and cϕ denotes the cosine of angle ϕ. In the hypothesis of small rotations, the angular velocity of the moving platform is given by [24, 25](2)ω=[ϕ̇xϕ̇yϕ̇z]T.Figure 2
Vector diagram of a PSS kinematic.The orientation of each kinematic strut with respect to the fixed base can be described by two Euler angles. As shown in Figure3, the local coordinate system of the ith strut can be thought of as a rotation of ϕi about the z axis resulting in a Ci-xi′yi′zi′ system followed by another rotation of φi about the rotated  yi′-axis. So the rotation matrix of the ith strut can be written as(3)Rio=Rot(z,ϕi)Rot(yi′,φi)=[cϕicφi-sϕicϕisφisϕicφicϕisϕisφi-sφi0cφi],i=1,2,…,6.
The unit vector along the strut in the coordinate system O-xyz is(4)wi=Riowii=Rio[001]=[cϕisφisϕisφicφi].
So the Euler angles ϕi and φi can be computed as follows:(5)cφi=wiz,sφi=wix2+wiy2,(0≤φi<π),sϕi=wiysφi,(φi≠0),cϕi=wixsφi,(φi≠0),ifφi=0,thenϕi=0.Figure 3
The local coordinate system of theith strut.


### 2.2. Rigid Dynamics

When the seismic simulator is not at a singular configuration, the rigid dynamic model can be formulated by means of the principle of virtual work and the concept link Jacobian matrices [25]. It can be expressed as(6)F=-J-T[fene]-J-T{[mpg0]+∑i=16JivωT[miRoig0]+JT[(mc1g)Te1(mc2g)Te2(mc3g)Te3(mc4g)Te4(mc5g)Te5(mc6g)Te6]T∑i=16JivωT}+J-T{[mpv̇Ipoω̇]+∑i=16JivωT[miv̇iiIiiω̇ii]+JT[mc1q̈1mc2q̈2mc3q̈3mc4q̈4mc5q̈5mc6q̈6]T}+J-T{[0ω×(Ipoω)]+∑i=16JivωT[0ωii×(Iiiωii)]},where 0=[000]T,Jivω is the link Jacobian matrix which maps the velocity of the moving platform into the velocity of the ith strut in the Ci-xiyizi coordinate system, mp, mci, and mi denote the mass of the moving platform, the mass of the slider, and the mass of the ith strut, respectively, Ipo is the inertia matrix of the moving platform taken about the center of mass expressed in the O-xyz coordinate system, Iii is the inertia matrix of the ith cylindrical strut about their respective centers of mass expressed in the Ci-xiyizi coordinate system, fe and ne are the external force and moment exerted at the center of mass of the moving platform, v̇ and ω̇ are the linear and angular acceleration of the moving platform,v̇ii,ωii, and ω̇ii are the linear acceleration, the angular velocity, and acceleration of the ith strut expressed in the Ci-xiyizi coordinate system, respectively. q̈i denotes the joint acceleration, g  is the gravity acceleration, J  is the Jacobian matrix which maps the velocity vector of the moving platform into the velocity vector of the actuating joint. F  is the input force vector exerted at the center of the slider. And(7)J=diag⁡(1w1Te11w2Te21w3Te31w4Te41w5Te51w6Te6)×[w1w2w3w4w5w6a1×w1a2×w2a3×w3a4×w4a5×w5a6×w6]T,Jivω=[[Roi-S(aii)Roi]+li2S(wii)Jiω1li{[S(wii)Roi-S(wii)S(aii)Roi]-(wii×eii)[wiTwiTei(ai×wi)TwiTei]}]=[JivJiω],S(wii)=[0-wiziwiyiwizi0-wixi-wiyiwixi0],S(aii)=[0-aiziaiyiaizi0-aixi-aiyiaixi0],iωi=Jiω[vω],q̈i=1wiTei(wiTv̇+(ai×wi)Tω̇+wiT(ω×(ω×ai))-wiT(ωi×(ωi×liwi))),iω̇i=Jiω[v̇ω̇]+1li(Δ1+Δ2),Δ1=-(wii×eii)wiTei((wiTω)(aiTω)-(wiTai)(ωTω)+li|ωi×wi|2),Δ2=(iωiTaii)(wii×iωi)-(iωTiω)(wii×aii),v̇ii=Jiv[v̇ω̇]+S(iωi)S(iωi)aii+12S(wii)(Δ1+Δ2)-li2S(iωi)S(iωi)wii.ei, ai, and wi are shown in Figure 2; they are the unit vector along the lead screw, the vector O′Ai, and the unit vector along strut CiAi, respectively.


## 2.1. Description

The schematic diagram of the 6-dof parallel structure seismic simulator is shown in Figure1. As shown in Figure 1, the parallel manipulator is composed of a moving platform and six sliders. In each kinematic chain, the platform and the slider are connected via spherical ball-bearing joints by a strut of fixed length. Each slider is driven by DC motor via a linear ball screw. The lead screws of B1,B2, and B3 are vertical to the ground.Figure 1
Schematic diagram of the 6-dof parallel structure seismic simulator.For the purpose of analysis, the following coordinate systems are defined. As shown in Figure2, the coordinate system O-xyz  is attached to the fixed base; another moving coordinate frame O′-uvw  is located at the center of mass of the moving platform. The pose of the moving platform can be described by a position vector r and a rotation matrixRo′o. Let the rotation matrix be defined by the roll, pitch, and yaw angles, namely, a rotation of ϕx about the fixed x axis, followed by a rotation of ϕy about the fixed y axis, and a rotation of ϕz about the fixed z axis. Thus, the rotation matrix is(1)Ro′o=Rot(z,ϕz)Rot(y,ϕy)Rot(x,ϕx),
where sϕ denotes the sine of angle ϕ, and cϕ denotes the cosine of angle ϕ. In the hypothesis of small rotations, the angular velocity of the moving platform is given by [24, 25](2)ω=[ϕ̇xϕ̇yϕ̇z]T.Figure 2
Vector diagram of a PSS kinematic.The orientation of each kinematic strut with respect to the fixed base can be described by two Euler angles. As shown in Figure3, the local coordinate system of the ith strut can be thought of as a rotation of ϕi about the z axis resulting in a Ci-xi′yi′zi′ system followed by another rotation of φi about the rotated  yi′-axis. So the rotation matrix of the ith strut can be written as(3)Rio=Rot(z,ϕi)Rot(yi′,φi)=[cϕicφi-sϕicϕisφisϕicφicϕisϕisφi-sφi0cφi],i=1,2,…,6.
The unit vector along the strut in the coordinate system O-xyz is(4)wi=Riowii=Rio[001]=[cϕisφisϕisφicφi].
So the Euler angles ϕi and φi can be computed as follows:(5)cφi=wiz,sφi=wix2+wiy2,(0≤φi<π),sϕi=wiysφi,(φi≠0),cϕi=wixsφi,(φi≠0),ifφi=0,thenϕi=0.Figure 3
The local coordinate system of theith strut.


## 2.2. Rigid Dynamics

When the seismic simulator is not at a singular configuration, the rigid dynamic model can be formulated by means of the principle of virtual work and the concept link Jacobian matrices [25]. It can be expressed as(6)F=-J-T[fene]-J-T{[mpg0]+∑i=16JivωT[miRoig0]+JT[(mc1g)Te1(mc2g)Te2(mc3g)Te3(mc4g)Te4(mc5g)Te5(mc6g)Te6]T∑i=16JivωT}+J-T{[mpv̇Ipoω̇]+∑i=16JivωT[miv̇iiIiiω̇ii]+JT[mc1q̈1mc2q̈2mc3q̈3mc4q̈4mc5q̈5mc6q̈6]T}+J-T{[0ω×(Ipoω)]+∑i=16JivωT[0ωii×(Iiiωii)]},where 0=[000]T,Jivω is the link Jacobian matrix which maps the velocity of the moving platform into the velocity of the ith strut in the Ci-xiyizi coordinate system, mp, mci, and mi denote the mass of the moving platform, the mass of the slider, and the mass of the ith strut, respectively, Ipo is the inertia matrix of the moving platform taken about the center of mass expressed in the O-xyz coordinate system, Iii is the inertia matrix of the ith cylindrical strut about their respective centers of mass expressed in the Ci-xiyizi coordinate system, fe and ne are the external force and moment exerted at the center of mass of the moving platform, v̇ and ω̇ are the linear and angular acceleration of the moving platform,v̇ii,ωii, and ω̇ii are the linear acceleration, the angular velocity, and acceleration of the ith strut expressed in the Ci-xiyizi coordinate system, respectively. q̈i denotes the joint acceleration, g  is the gravity acceleration, J  is the Jacobian matrix which maps the velocity vector of the moving platform into the velocity vector of the actuating joint. F  is the input force vector exerted at the center of the slider. And(7)J=diag⁡(1w1Te11w2Te21w3Te31w4Te41w5Te51w6Te6)×[w1w2w3w4w5w6a1×w1a2×w2a3×w3a4×w4a5×w5a6×w6]T,Jivω=[[Roi-S(aii)Roi]+li2S(wii)Jiω1li{[S(wii)Roi-S(wii)S(aii)Roi]-(wii×eii)[wiTwiTei(ai×wi)TwiTei]}]=[JivJiω],S(wii)=[0-wiziwiyiwizi0-wixi-wiyiwixi0],S(aii)=[0-aiziaiyiaizi0-aixi-aiyiaixi0],iωi=Jiω[vω],q̈i=1wiTei(wiTv̇+(ai×wi)Tω̇+wiT(ω×(ω×ai))-wiT(ωi×(ωi×liwi))),iω̇i=Jiω[v̇ω̇]+1li(Δ1+Δ2),Δ1=-(wii×eii)wiTei((wiTω)(aiTω)-(wiTai)(ωTω)+li|ωi×wi|2),Δ2=(iωiTaii)(wii×iωi)-(iωTiω)(wii×aii),v̇ii=Jiv[v̇ω̇]+S(iωi)S(iωi)aii+12S(wii)(Δ1+Δ2)-li2S(iωi)S(iωi)wii.ei, ai, and wi are shown in Figure 2; they are the unit vector along the lead screw, the vector O′Ai, and the unit vector along strut CiAi, respectively.


## 3. Kineto-Elastodynamic Model

The idea of substructure synthesis and the finite element method are adopted to develop the Kineto-elastodynamic model of the 6-PSS parallel structure seismic simulator. The finite element method used here is based on the basic assumptions [26] as follows. (1) The deflections of the links of the manipulator obey the small deflection theory. The small amplitude structural vibrations do not have a significant effect on its rigid-body motion and the coupling term between the elastic deformation and the rigid-body motion is neglected. The true motion is regarded as the sum of the rigid-body motion, and the elastic motion. (2) The instantaneous structural approach is adopted. At each instant, the manipulator is modeled as an instantaneous structure undergoing elastic deformations about its mean rigid body configuration. (3) The model is based on the Euler-Bernoulli beam theory. (4) The transverse deflections are modeled as a cubic polynomial of the nodal displacement, the longitudinal deflections and the torsional deflections are modeled as a first-order polynomial of the nodal displacement. The manipulator is divided into seven substructures, namely, one moving platform substructure and six kinematic chain substructures which are composed of the lead-screw assembly and the strut. Each strut is divided into three elements. The moving platform and the sliders are regarded as the rigid bodies since their deformations are small relative to the elastic deformations.


### 3.1. Strut Dynamic Equation

#### 3.1.1. Element Model

The nodal elastic displacement of the element is shown in Figure4. So the elastic displacement of the element can be expressed as (8)δ=[δiδj]T,
where(9)δi=[uiviwiθixθiyθiz]T,δj=[ujvjwjθjxθjyθjz]T.Figure 4
Nodal elastic displacement of the element.The elastic displacement vector of an arbitrary pointP within the element can be expressed by the nodal displacement of the element [27](10)p=[uvwθ]T=Nδ,
where N is relational matrix which maps the nodal elastic displacement vector of the element into that of the point P.The polynomial of the nodal displacement of the element is chosen to formulate the displacement of the pointP. The transverse displacement, the longitudinal displacement, and the torsional displacement of the point  P are modeled as a cubic polynomial and a linear function of the nodal displacement, respectively. So the longitudinal displacement and the torsional displacement of the point  P are expressed as(11)u=a0+a1x,θ=d0+d1x.
The transverse displacements of the point P can be expressed as (12)v=b0+b1x+b2x2+b3x3,w=c0+c1x+c2x2+c3x3.
Substituting (11) and (12) into (10) yields(13)p=[uvwθ]=[Hu(x)aHv(x)bHw(x)cHθ(x)d],
where(14)a=[a0a1]T,b=[b0b1b2b3]T,c=[c0c1c2c3]T,d=[d0d1]T,Hu(x)=Hθ(x)=[1x],Hv(x)=Hw(x)=[1xx2x3].Considering the nodei and the node j, where x=0 and x=l, respectively, yields(15)a=Auθ-1δu,b=Avw-1δv,c=Avw-1δw,d=Auθ-1δθ,
where(16)δu=[uiuj]T,δv=[viθizvjθjz]T,δw=[wiθiywjθjy]T,δθ=[θixθjx]T,Auθ-1=[10-1l1l],Avw-1=[10000100-3l2-2l3l2-1l2l31l2-2l31l2].
Substituting (15) into (13) yields
(17)p=[uvwθ]=[hu(x)hv(x)hw(x)hθ(x)]Aδ=Nδ,
where(18)hu(x)=[100000x00000],hv(x)=[01000x0x2000x3],hw(x)=[0010x000x20x30],hθ(x)=[000100000x00],A=[100000000000010000000000001000000000000100000000000010000000000001000000-1l000001l000000-3l2000-2l03l2000-1l00-3l20-2l0003l20-1l0000-1l000001l00002l301l2000-2l301l2002l30001l20-2l30001l2].Considering the knowledge of material mechanics,  the strain of the pointP is
(19)ε=[ε0εbyεbzεr]=[dudx-yd2vdx2-zd2wdx2GE16Jkπd3dθdx],
where ε0 is the axial strain, εby and εbz are the flexural strain in the plane Oxy and Oxz, respectively, εr is the torsional strain, y and z are the distance along the y and z direction from the axis of the element to the point P, G is the torsional modulus, E is the Young’s modulus, Jk is the polar moment of inertia of cross-section, d is the diameter of the strut.Substituting (17) into (19) yields(20)ε=[hu′(x)-yhv′′(x)-zhw′′(x)GE16Jkπd3hθ′(x)]Aδ=Bδ,
where(21)hu′(x)=[000000100000],hv′′(x)=[000000020006x],hw′′(x)=[00000000206x0],hθ′(x)=[000000000100].
So the stress is (22)σ=Eε.According to the knowledge of material mechanics, the strain energy of the element can be expressed as(23)U=12∭εTσdv.
Substituting (20) and (22) into (23) yields (24)U=12∭EδTBTBδdv=12δTE∭BTBdvδ=12δTkδ,
where(25)k=12E∭BTBdv
is the element stiffness matrix. Substituting (20) into (25) yields(26)k=[EAl00000-EAl00000012EIzl30006EIzl20-12EIzl30006EIzl20012EIyl30-6EIyl2000-12EIyl30-6EIyl20000GJkl00000-GJkl0000-6EIyl204EIyl0006EIyl202EIyl006EIzl20004EIzl0-6EIzl20002EIzl-EAl00000EAl000000-12EIzl3000-6EIzl2012EIzl3000-6EIzl200-12EIyl306EIyl200012EIyl306EIyl20000-GJkl00000GJkl0000-6EIyl202EIyl0006EIyl204EIyl006EIzl20002EIzl0-6EIzl20004EIzl],where  Iy  and  Iz  are the principal moments of inertia corresponding to the y axis and z axis, respectively, A  is the area of the cross-section of the uniform beam.Based on the presented assumption, there are kinematics relationships as follows:(27)ṗa=ṗr+ṗ,p̈a=p̈r+p̈,
where  ṗa  and  p̈a  are the absolute velocity and the absolute acceleration of a certain point within the element,  ṗr  and  p̈r  are the velocity and the acceleration of the rigid-body motion, and  ṗ  and  p̈  are the velocity and the acceleration of the elastic motion.It can be proved that (see [27])(28)ṗr=Nrδ̇r=Nδ̇r,p̈r=Nrδ̈r=Nδ̈r,
where  δ̇r  and  δ̈r  are the nodal velocity and the acceleration of the rigid-body motion.So the kinetic energy of the element can be expressed as(29)T=12∫0lm(x)ṗaTṗadx=12∫0lm(x)δ̇aTNTNδ̇adx=12δ̇aT∫0lm(x)NTNdxδ̇a=12δ̇aTmδ̇a,
where  δ̇a  is the absolute velocity of the node. And (30)m=12∫0lm(x)NTNdx
is the element mass matrix. For the uniform beam, the mass function is(31)m(x)=ρA,
where  ρ  is the density of the beam. So the element mass matrix is(32)m=ρAl420[140000007000000015600022l054000-13l001560-22l00054013l0000140JkA0000070JkA0000-22l04l2000-13l0-3l20022l0004l2013l000-3l270000001400000005400013l0156000-22l00540-13l000156022l000070JkA00000140JkA000013l0-3l200022l04l200-13l000-3l20-22l0004l2].Substituting (23) and (29) into Lagrange equation yields
(33)ddt(∂T∂δ̇)-∂T∂δ+∂U∂δ=f,
where f denotes the resultant of the applied and the internal forces exerted at the element. Thus, the kinematic differential equation of the strut element within ith kinematic chain can be achieved as follows:(34)mδ̈+kδ=f-mδ̈r.
Equation (34) can be expressed in the coordinate system O-xyz  as(35)MeiÜei+KeiUei=Fei-MeiÜeri,i=1,2,…,6
where(36)Uei=RiTδ,Üei=RiTδ̈,Üeri=RiTδ̈r,Fei=RiTf,Mei=RiTmRi,Kei=RiTkRi,Ri=diag⁡((RioRot(y,-π2))T(RioRot(y,-π2))T(RioRot(y,-π2))T(RioRot(y,-π2))T)


#### 3.1.2. Strut Dynamic Equation

For the strut withinith kinematic chain, the dynamic motion in the coordinate system O-xyz can be assembled as(37)Msi′Üsi′+Ksi′Usi′=Fsi′-Msi′Üsri′,
where
(38)Msi′=∑j=13AjTMeiAj,Ksi′=∑j=13AjTKeiAj,Fsi′=∑j=13AjTFei=[Fc2i*018×1Fc1i*],A1=[E12012×12],A2=[012×6E12012×6],A3=[012×12E12],
where Aj is the connectivity matrix which maps the total nodal coordinates Usi′ of the strut within the ith kinematics chain into the jth nodal coordinate within the strut. Üsi′ and Üsri′ denote the the nodal acceleration of the elastic motion and the rigid-body motion of the strut within the ith kinematics chain. E12 is the unit matrix of order twelve. Fc1i* is the internal force between the strut and the slider. Fc2i* is the internal force between the rigid moving platform and the strut. All the above coordinates are measured in the coordinate system O-xyz.


### 3.2. Slider Dynamic Equation

The oscillation of the slider along the axial direction of the lead-screw can be expressed as(39)mciÜci+kciUci=fi-mciq̈i-eiTFc1i*,
where Uci is the elastic displacement of the slider along the axial direction of the lead-screw.  kci  is the equivalent axial stiffness of the lead screw assembly that is composed of three serially connected component, that is, lead screw, ball nut which links slider with lead-screw and the two sets of support bearings at both ends. Let  kcsi,kcni, and  kcbi  denote their respective axial stiffness then the equivalent axial stiffness can be expressed as(40)1kci=1kcsi+1kcni+1kcbi,kcsi=AsEs(L1i+L2i)L1iL2i,
where As and Es stand for the cross-sectional area of the lead screw and its Young’s modular, and L1i and L2i are the distances between the nut and the two sets of support bearings located at each end of the lead screw.


### 3.3. Deformation Compatibility Condition

The compatibility of the deformations between the rigid moving platform and the flexible strut can be expressed as(41)Uc2i=[E3-S(ai)]Up,
where Up  denotes the generalized coordinates of vibration motion of the moving platform,  Uc2i=Usi′(1:3,1:1) denotes the elastic displacement of the corresponding node within the flexible strut, S(ai) is the screw matrix of ai, and E3 denotes the unit matrix of the order three.The compatibility of the deformations between the rigid slider and the flexible strut can be expressed as(42)eiTUc1i=Uci,
where Uc1i=Usi′(18:21,1:1) is the elastic displacement of the corresponding node within the flexible strut.


### 3.4. Substructure Motion Equation

#### 3.4.1. Moving Platform Substructure

The oscillation equation of the rigid moving platform substructure in the coordinate systemO-xyz is(43)[mpE303×303×3Ipo]Üp-[fe-mpv̇-∑i=16Fc2i*ne-Ipoω̇-ω×(Ipoω)-∑i=16ai×Fc2i*]=0.


#### 3.4.2. Kinematic Chain Substructure

Employing the deformation compatibility conditions between the flexible strut and the rigid slider and the boundary conditions of the slider, the motion equation of theith kinematic chain substructure can be assembled as(44)MkiÜki+KkiUki=Fkci+Fkdi,
where
(45)Uki=[Usi′(1:18,1:1)Uci],Mki=[Msi′(1:18,1:18)Msi′(1:18,21:21)Msi′(21:21,1:18)Msi′(21:21,21:21)+mci]i=1,2,3,Mki=[Msi′(1:18,1:18)Msi′(1:18,20:20)Msi′(20:20,1:18)Msi′(20:20,20:20)+mci]i=4,5,Mki=[Msi′(1:18,1:18)Msi′(1:18,19:19)Msi′(19:19,1:18)Msi′(19:19,19:19)+mci]i=6,Kki=[Ksi′(1:18,1:18)Ksi′(1:18,21:21)Ksi′(21:21,1:18)Ksi′(21:21,21:21)+kci]i=1,2,3,Kki=[Ksi′(1:18,1:18)Ksi′(1:18,20:20)Ksi′(20:20,1:18)Ksi′(20:20,20:20)+kci]i=4,5,Kki=[Ksi′(1:18,1:18)Ksi′(1:18,19:19)Ksi′(19:19,1:18)Ksi′(19:19,19:19)+kci]i=6,Fkci=[Fc2i*015×10],Fkdi=[018×1fi]-MkiÜkri,Ukri=[Usri′(1:18,1:1)q̈i]Fkci is the internal forces between the elements within the strut, Fkdi is the resultant force of the generalized inertial force and the outside force.


### 3.5. Kineto-elasticdynamic Model of the Manipulator

Gathering the dynamic equations of the substructures and employing the deformation compatibility conditions between the rigid moving platform and the flexible strut yields(46)DTM′DÜ+DTK′DU=DTFc′+DTFd′,
where(47)D=[E606×96E3-S(a1)03×96016×6E16016×80E3-S(a2)03×96016×22E16016×64E3-S(a3)03×96016×38E16016×48E3-S(a4)03×96016×54E16016×32E3-S(a5)03×96016×70E16016×16E3-S(a6)03×96016×86E16],U=[UpUk1(4:19,1)Uk2(4:19,1)Uk3(4:19,1)Uk4(4:19,1)Uk5(4:19,1)Uk6(4:19,1)],Fc′=[-[∑i=16Fc2i*]T-[∑i=16ai×Fc2i*]TFkc1TFkc2TFkc3TFkc4TFkc5TFkc6T]T,Fd′=[[fe-mpv̇]T[ne-Ipoω̇-ω×(Ipoω)]TFkd1TFkd2TFkd3TFkd4TFkd5TFkd6T[ne-Ipoω̇-ω×(Ipoω)]T]T Simplifying (46) yields the kineto-elasticdynamic model of the manipulator(48)MÜ+KU=Fd,
where(49)M=DTM′D,(50)K=DTK′D,(51)Fd=DTFd′.


## 3.1. Strut Dynamic Equation

### 3.1.1. Element Model

The nodal elastic displacement of the element is shown in Figure4. So the elastic displacement of the element can be expressed as (8)δ=[δiδj]T,
where(9)δi=[uiviwiθixθiyθiz]T,δj=[ujvjwjθjxθjyθjz]T.Figure 4
Nodal elastic displacement of the element.The elastic displacement vector of an arbitrary pointP within the element can be expressed by the nodal displacement of the element [27](10)p=[uvwθ]T=Nδ,
where N is relational matrix which maps the nodal elastic displacement vector of the element into that of the point P.The polynomial of the nodal displacement of the element is chosen to formulate the displacement of the pointP. The transverse displacement, the longitudinal displacement, and the torsional displacement of the point  P are modeled as a cubic polynomial and a linear function of the nodal displacement, respectively. So the longitudinal displacement and the torsional displacement of the point  P are expressed as(11)u=a0+a1x,θ=d0+d1x.
The transverse displacements of the point P can be expressed as (12)v=b0+b1x+b2x2+b3x3,w=c0+c1x+c2x2+c3x3.
Substituting (11) and (12) into (10) yields(13)p=[uvwθ]=[Hu(x)aHv(x)bHw(x)cHθ(x)d],
where(14)a=[a0a1]T,b=[b0b1b2b3]T,c=[c0c1c2c3]T,d=[d0d1]T,Hu(x)=Hθ(x)=[1x],Hv(x)=Hw(x)=[1xx2x3].Considering the nodei and the node j, where x=0 and x=l, respectively, yields(15)a=Auθ-1δu,b=Avw-1δv,c=Avw-1δw,d=Auθ-1δθ,
where(16)δu=[uiuj]T,δv=[viθizvjθjz]T,δw=[wiθiywjθjy]T,δθ=[θixθjx]T,Auθ-1=[10-1l1l],Avw-1=[10000100-3l2-2l3l2-1l2l31l2-2l31l2].
Substituting (15) into (13) yields
(17)p=[uvwθ]=[hu(x)hv(x)hw(x)hθ(x)]Aδ=Nδ,
where(18)hu(x)=[100000x00000],hv(x)=[01000x0x2000x3],hw(x)=[0010x000x20x30],hθ(x)=[000100000x00],A=[100000000000010000000000001000000000000100000000000010000000000001000000-1l000001l000000-3l2000-2l03l2000-1l00-3l20-2l0003l20-1l0000-1l000001l00002l301l2000-2l301l2002l30001l20-2l30001l2].Considering the knowledge of material mechanics,  the strain of the pointP is
(19)ε=[ε0εbyεbzεr]=[dudx-yd2vdx2-zd2wdx2GE16Jkπd3dθdx],
where ε0 is the axial strain, εby and εbz are the flexural strain in the plane Oxy and Oxz, respectively, εr is the torsional strain, y and z are the distance along the y and z direction from the axis of the element to the point P, G is the torsional modulus, E is the Young’s modulus, Jk is the polar moment of inertia of cross-section, d is the diameter of the strut.Substituting (17) into (19) yields(20)ε=[hu′(x)-yhv′′(x)-zhw′′(x)GE16Jkπd3hθ′(x)]Aδ=Bδ,
where(21)hu′(x)=[000000100000],hv′′(x)=[000000020006x],hw′′(x)=[00000000206x0],hθ′(x)=[000000000100].
So the stress is (22)σ=Eε.According to the knowledge of material mechanics, the strain energy of the element can be expressed as(23)U=12∭εTσdv.
Substituting (20) and (22) into (23) yields (24)U=12∭EδTBTBδdv=12δTE∭BTBdvδ=12δTkδ,
where(25)k=12E∭BTBdv
is the element stiffness matrix. Substituting (20) into (25) yields(26)k=[EAl00000-EAl00000012EIzl30006EIzl20-12EIzl30006EIzl20012EIyl30-6EIyl2000-12EIyl30-6EIyl20000GJkl00000-GJkl0000-6EIyl204EIyl0006EIyl202EIyl006EIzl20004EIzl0-6EIzl20002EIzl-EAl00000EAl000000-12EIzl3000-6EIzl2012EIzl3000-6EIzl200-12EIyl306EIyl200012EIyl306EIyl20000-GJkl00000GJkl0000-6EIyl202EIyl0006EIyl204EIyl006EIzl20002EIzl0-6EIzl20004EIzl],where  Iy  and  Iz  are the principal moments of inertia corresponding to the y axis and z axis, respectively, A  is the area of the cross-section of the uniform beam.Based on the presented assumption, there are kinematics relationships as follows:(27)ṗa=ṗr+ṗ,p̈a=p̈r+p̈,
where  ṗa  and  p̈a  are the absolute velocity and the absolute acceleration of a certain point within the element,  ṗr  and  p̈r  are the velocity and the acceleration of the rigid-body motion, and  ṗ  and  p̈  are the velocity and the acceleration of the elastic motion.It can be proved that (see [27])(28)ṗr=Nrδ̇r=Nδ̇r,p̈r=Nrδ̈r=Nδ̈r,
where  δ̇r  and  δ̈r  are the nodal velocity and the acceleration of the rigid-body motion.So the kinetic energy of the element can be expressed as(29)T=12∫0lm(x)ṗaTṗadx=12∫0lm(x)δ̇aTNTNδ̇adx=12δ̇aT∫0lm(x)NTNdxδ̇a=12δ̇aTmδ̇a,
where  δ̇a  is the absolute velocity of the node. And (30)m=12∫0lm(x)NTNdx
is the element mass matrix. For the uniform beam, the mass function is(31)m(x)=ρA,
where  ρ  is the density of the beam. So the element mass matrix is(32)m=ρAl420[140000007000000015600022l054000-13l001560-22l00054013l0000140JkA0000070JkA0000-22l04l2000-13l0-3l20022l0004l2013l000-3l270000001400000005400013l0156000-22l00540-13l000156022l000070JkA00000140JkA000013l0-3l200022l04l200-13l000-3l20-22l0004l2].Substituting (23) and (29) into Lagrange equation yields
(33)ddt(∂T∂δ̇)-∂T∂δ+∂U∂δ=f,
where f denotes the resultant of the applied and the internal forces exerted at the element. Thus, the kinematic differential equation of the strut element within ith kinematic chain can be achieved as follows:(34)mδ̈+kδ=f-mδ̈r.
Equation (34) can be expressed in the coordinate system O-xyz  as(35)MeiÜei+KeiUei=Fei-MeiÜeri,i=1,2,…,6
where(36)Uei=RiTδ,Üei=RiTδ̈,Üeri=RiTδ̈r,Fei=RiTf,Mei=RiTmRi,Kei=RiTkRi,Ri=diag⁡((RioRot(y,-π2))T(RioRot(y,-π2))T(RioRot(y,-π2))T(RioRot(y,-π2))T)


### 3.1.2. Strut Dynamic Equation

For the strut withinith kinematic chain, the dynamic motion in the coordinate system O-xyz can be assembled as(37)Msi′Üsi′+Ksi′Usi′=Fsi′-Msi′Üsri′,
where
(38)Msi′=∑j=13AjTMeiAj,Ksi′=∑j=13AjTKeiAj,Fsi′=∑j=13AjTFei=[Fc2i*018×1Fc1i*],A1=[E12012×12],A2=[012×6E12012×6],A3=[012×12E12],
where Aj is the connectivity matrix which maps the total nodal coordinates Usi′ of the strut within the ith kinematics chain into the jth nodal coordinate within the strut. Üsi′ and Üsri′ denote the the nodal acceleration of the elastic motion and the rigid-body motion of the strut within the ith kinematics chain. E12 is the unit matrix of order twelve. Fc1i* is the internal force between the strut and the slider. Fc2i* is the internal force between the rigid moving platform and the strut. All the above coordinates are measured in the coordinate system O-xyz.


## 3.1.1. Element Model

The nodal elastic displacement of the element is shown in Figure4. So the elastic displacement of the element can be expressed as (8)δ=[δiδj]T,
where(9)δi=[uiviwiθixθiyθiz]T,δj=[ujvjwjθjxθjyθjz]T.Figure 4
Nodal elastic displacement of the element.The elastic displacement vector of an arbitrary pointP within the element can be expressed by the nodal displacement of the element [27](10)p=[uvwθ]T=Nδ,
where N is relational matrix which maps the nodal elastic displacement vector of the element into that of the point P.The polynomial of the nodal displacement of the element is chosen to formulate the displacement of the pointP. The transverse displacement, the longitudinal displacement, and the torsional displacement of the point  P are modeled as a cubic polynomial and a linear function of the nodal displacement, respectively. So the longitudinal displacement and the torsional displacement of the point  P are expressed as(11)u=a0+a1x,θ=d0+d1x.
The transverse displacements of the point P can be expressed as (12)v=b0+b1x+b2x2+b3x3,w=c0+c1x+c2x2+c3x3.
Substituting (11) and (12) into (10) yields(13)p=[uvwθ]=[Hu(x)aHv(x)bHw(x)cHθ(x)d],
where(14)a=[a0a1]T,b=[b0b1b2b3]T,c=[c0c1c2c3]T,d=[d0d1]T,Hu(x)=Hθ(x)=[1x],Hv(x)=Hw(x)=[1xx2x3].Considering the nodei and the node j, where x=0 and x=l, respectively, yields(15)a=Auθ-1δu,b=Avw-1δv,c=Avw-1δw,d=Auθ-1δθ,
where(16)δu=[uiuj]T,δv=[viθizvjθjz]T,δw=[wiθiywjθjy]T,δθ=[θixθjx]T,Auθ-1=[10-1l1l],Avw-1=[10000100-3l2-2l3l2-1l2l31l2-2l31l2].
Substituting (15) into (13) yields
(17)p=[uvwθ]=[hu(x)hv(x)hw(x)hθ(x)]Aδ=Nδ,
where(18)hu(x)=[100000x00000],hv(x)=[01000x0x2000x3],hw(x)=[0010x000x20x30],hθ(x)=[000100000x00],A=[100000000000010000000000001000000000000100000000000010000000000001000000-1l000001l000000-3l2000-2l03l2000-1l00-3l20-2l0003l20-1l0000-1l000001l00002l301l2000-2l301l2002l30001l20-2l30001l2].Considering the knowledge of material mechanics,  the strain of the pointP is
(19)ε=[ε0εbyεbzεr]=[dudx-yd2vdx2-zd2wdx2GE16Jkπd3dθdx],
where ε0 is the axial strain, εby and εbz are the flexural strain in the plane Oxy and Oxz, respectively, εr is the torsional strain, y and z are the distance along the y and z direction from the axis of the element to the point P, G is the torsional modulus, E is the Young’s modulus, Jk is the polar moment of inertia of cross-section, d is the diameter of the strut.Substituting (17) into (19) yields(20)ε=[hu′(x)-yhv′′(x)-zhw′′(x)GE16Jkπd3hθ′(x)]Aδ=Bδ,
where(21)hu′(x)=[000000100000],hv′′(x)=[000000020006x],hw′′(x)=[00000000206x0],hθ′(x)=[000000000100].
So the stress is (22)σ=Eε.According to the knowledge of material mechanics, the strain energy of the element can be expressed as(23)U=12∭εTσdv.
Substituting (20) and (22) into (23) yields (24)U=12∭EδTBTBδdv=12δTE∭BTBdvδ=12δTkδ,
where(25)k=12E∭BTBdv
is the element stiffness matrix. Substituting (20) into (25) yields(26)k=[EAl00000-EAl00000012EIzl30006EIzl20-12EIzl30006EIzl20012EIyl30-6EIyl2000-12EIyl30-6EIyl20000GJkl00000-GJkl0000-6EIyl204EIyl0006EIyl202EIyl006EIzl20004EIzl0-6EIzl20002EIzl-EAl00000EAl000000-12EIzl3000-6EIzl2012EIzl3000-6EIzl200-12EIyl306EIyl200012EIyl306EIyl20000-GJkl00000GJkl0000-6EIyl202EIyl0006EIyl204EIyl006EIzl20002EIzl0-6EIzl20004EIzl],where  Iy  and  Iz  are the principal moments of inertia corresponding to the y axis and z axis, respectively, A  is the area of the cross-section of the uniform beam.Based on the presented assumption, there are kinematics relationships as follows:(27)ṗa=ṗr+ṗ,p̈a=p̈r+p̈,
where  ṗa  and  p̈a  are the absolute velocity and the absolute acceleration of a certain point within the element,  ṗr  and  p̈r  are the velocity and the acceleration of the rigid-body motion, and  ṗ  and  p̈  are the velocity and the acceleration of the elastic motion.It can be proved that (see [27])(28)ṗr=Nrδ̇r=Nδ̇r,p̈r=Nrδ̈r=Nδ̈r,
where  δ̇r  and  δ̈r  are the nodal velocity and the acceleration of the rigid-body motion.So the kinetic energy of the element can be expressed as(29)T=12∫0lm(x)ṗaTṗadx=12∫0lm(x)δ̇aTNTNδ̇adx=12δ̇aT∫0lm(x)NTNdxδ̇a=12δ̇aTmδ̇a,
where  δ̇a  is the absolute velocity of the node. And (30)m=12∫0lm(x)NTNdx
is the element mass matrix. For the uniform beam, the mass function is(31)m(x)=ρA,
where  ρ  is the density of the beam. So the element mass matrix is(32)m=ρAl420[140000007000000015600022l054000-13l001560-22l00054013l0000140JkA0000070JkA0000-22l04l2000-13l0-3l20022l0004l2013l000-3l270000001400000005400013l0156000-22l00540-13l000156022l000070JkA00000140JkA000013l0-3l200022l04l200-13l000-3l20-22l0004l2].Substituting (23) and (29) into Lagrange equation yields
(33)ddt(∂T∂δ̇)-∂T∂δ+∂U∂δ=f,
where f denotes the resultant of the applied and the internal forces exerted at the element. Thus, the kinematic differential equation of the strut element within ith kinematic chain can be achieved as follows:(34)mδ̈+kδ=f-mδ̈r.
Equation (34) can be expressed in the coordinate system O-xyz  as(35)MeiÜei+KeiUei=Fei-MeiÜeri,i=1,2,…,6
where(36)Uei=RiTδ,Üei=RiTδ̈,Üeri=RiTδ̈r,Fei=RiTf,Mei=RiTmRi,Kei=RiTkRi,Ri=diag⁡((RioRot(y,-π2))T(RioRot(y,-π2))T(RioRot(y,-π2))T(RioRot(y,-π2))T)


## 3.1.2. Strut Dynamic Equation

For the strut withinith kinematic chain, the dynamic motion in the coordinate system O-xyz can be assembled as(37)Msi′Üsi′+Ksi′Usi′=Fsi′-Msi′Üsri′,
where
(38)Msi′=∑j=13AjTMeiAj,Ksi′=∑j=13AjTKeiAj,Fsi′=∑j=13AjTFei=[Fc2i*018×1Fc1i*],A1=[E12012×12],A2=[012×6E12012×6],A3=[012×12E12],
where Aj is the connectivity matrix which maps the total nodal coordinates Usi′ of the strut within the ith kinematics chain into the jth nodal coordinate within the strut. Üsi′ and Üsri′ denote the the nodal acceleration of the elastic motion and the rigid-body motion of the strut within the ith kinematics chain. E12 is the unit matrix of order twelve. Fc1i* is the internal force between the strut and the slider. Fc2i* is the internal force between the rigid moving platform and the strut. All the above coordinates are measured in the coordinate system O-xyz.


## 3.2. Slider Dynamic Equation

The oscillation of the slider along the axial direction of the lead-screw can be expressed as(39)mciÜci+kciUci=fi-mciq̈i-eiTFc1i*,
where Uci is the elastic displacement of the slider along the axial direction of the lead-screw.  kci  is the equivalent axial stiffness of the lead screw assembly that is composed of three serially connected component, that is, lead screw, ball nut which links slider with lead-screw and the two sets of support bearings at both ends. Let  kcsi,kcni, and  kcbi  denote their respective axial stiffness then the equivalent axial stiffness can be expressed as(40)1kci=1kcsi+1kcni+1kcbi,kcsi=AsEs(L1i+L2i)L1iL2i,
where As and Es stand for the cross-sectional area of the lead screw and its Young’s modular, and L1i and L2i are the distances between the nut and the two sets of support bearings located at each end of the lead screw.


## 3.3. Deformation Compatibility Condition

The compatibility of the deformations between the rigid moving platform and the flexible strut can be expressed as(41)Uc2i=[E3-S(ai)]Up,
where Up  denotes the generalized coordinates of vibration motion of the moving platform,  Uc2i=Usi′(1:3,1:1) denotes the elastic displacement of the corresponding node within the flexible strut, S(ai) is the screw matrix of ai, and E3 denotes the unit matrix of the order three.The compatibility of the deformations between the rigid slider and the flexible strut can be expressed as(42)eiTUc1i=Uci,
where Uc1i=Usi′(18:21,1:1) is the elastic displacement of the corresponding node within the flexible strut.


## 3.4. Substructure Motion Equation

### 3.4.1. Moving Platform Substructure

The oscillation equation of the rigid moving platform substructure in the coordinate systemO-xyz is(43)[mpE303×303×3Ipo]Üp-[fe-mpv̇-∑i=16Fc2i*ne-Ipoω̇-ω×(Ipoω)-∑i=16ai×Fc2i*]=0.


### 3.4.2. Kinematic Chain Substructure

Employing the deformation compatibility conditions between the flexible strut and the rigid slider and the boundary conditions of the slider, the motion equation of theith kinematic chain substructure can be assembled as(44)MkiÜki+KkiUki=Fkci+Fkdi,
where
(45)Uki=[Usi′(1:18,1:1)Uci],Mki=[Msi′(1:18,1:18)Msi′(1:18,21:21)Msi′(21:21,1:18)Msi′(21:21,21:21)+mci]i=1,2,3,Mki=[Msi′(1:18,1:18)Msi′(1:18,20:20)Msi′(20:20,1:18)Msi′(20:20,20:20)+mci]i=4,5,Mki=[Msi′(1:18,1:18)Msi′(1:18,19:19)Msi′(19:19,1:18)Msi′(19:19,19:19)+mci]i=6,Kki=[Ksi′(1:18,1:18)Ksi′(1:18,21:21)Ksi′(21:21,1:18)Ksi′(21:21,21:21)+kci]i=1,2,3,Kki=[Ksi′(1:18,1:18)Ksi′(1:18,20:20)Ksi′(20:20,1:18)Ksi′(20:20,20:20)+kci]i=4,5,Kki=[Ksi′(1:18,1:18)Ksi′(1:18,19:19)Ksi′(19:19,1:18)Ksi′(19:19,19:19)+kci]i=6,Fkci=[Fc2i*015×10],Fkdi=[018×1fi]-MkiÜkri,Ukri=[Usri′(1:18,1:1)q̈i]Fkci is the internal forces between the elements within the strut, Fkdi is the resultant force of the generalized inertial force and the outside force.


## 3.4.1. Moving Platform Substructure

The oscillation equation of the rigid moving platform substructure in the coordinate systemO-xyz is(43)[mpE303×303×3Ipo]Üp-[fe-mpv̇-∑i=16Fc2i*ne-Ipoω̇-ω×(Ipoω)-∑i=16ai×Fc2i*]=0.


## 3.4.2. Kinematic Chain Substructure

Employing the deformation compatibility conditions between the flexible strut and the rigid slider and the boundary conditions of the slider, the motion equation of theith kinematic chain substructure can be assembled as(44)MkiÜki+KkiUki=Fkci+Fkdi,
where
(45)Uki=[Usi′(1:18,1:1)Uci],Mki=[Msi′(1:18,1:18)Msi′(1:18,21:21)Msi′(21:21,1:18)Msi′(21:21,21:21)+mci]i=1,2,3,Mki=[Msi′(1:18,1:18)Msi′(1:18,20:20)Msi′(20:20,1:18)Msi′(20:20,20:20)+mci]i=4,5,Mki=[Msi′(1:18,1:18)Msi′(1:18,19:19)Msi′(19:19,1:18)Msi′(19:19,19:19)+mci]i=6,Kki=[Ksi′(1:18,1:18)Ksi′(1:18,21:21)Ksi′(21:21,1:18)Ksi′(21:21,21:21)+kci]i=1,2,3,Kki=[Ksi′(1:18,1:18)Ksi′(1:18,20:20)Ksi′(20:20,1:18)Ksi′(20:20,20:20)+kci]i=4,5,Kki=[Ksi′(1:18,1:18)Ksi′(1:18,19:19)Ksi′(19:19,1:18)Ksi′(19:19,19:19)+kci]i=6,Fkci=[Fc2i*015×10],Fkdi=[018×1fi]-MkiÜkri,Ukri=[Usri′(1:18,1:1)q̈i]Fkci is the internal forces between the elements within the strut, Fkdi is the resultant force of the generalized inertial force and the outside force.


## 3.5. Kineto-elasticdynamic Model of the Manipulator

Gathering the dynamic equations of the substructures and employing the deformation compatibility conditions between the rigid moving platform and the flexible strut yields(46)DTM′DÜ+DTK′DU=DTFc′+DTFd′,
where(47)D=[E606×96E3-S(a1)03×96016×6E16016×80E3-S(a2)03×96016×22E16016×64E3-S(a3)03×96016×38E16016×48E3-S(a4)03×96016×54E16016×32E3-S(a5)03×96016×70E16016×16E3-S(a6)03×96016×86E16],U=[UpUk1(4:19,1)Uk2(4:19,1)Uk3(4:19,1)Uk4(4:19,1)Uk5(4:19,1)Uk6(4:19,1)],Fc′=[-[∑i=16Fc2i*]T-[∑i=16ai×Fc2i*]TFkc1TFkc2TFkc3TFkc4TFkc5TFkc6T]T,Fd′=[[fe-mpv̇]T[ne-Ipoω̇-ω×(Ipoω)]TFkd1TFkd2TFkd3TFkd4TFkd5TFkd6T[ne-Ipoω̇-ω×(Ipoω)]T]T Simplifying (46) yields the kineto-elasticdynamic model of the manipulator(48)MÜ+KU=Fd,
where(49)M=DTM′D,(50)K=DTK′D,(51)Fd=DTFd′.


## 4. Kineto-Elastodynamic Characteristics Analysis

In this section, the investigation on the Kineto-elastodynamic characteristics of the 6-PSS parallel structure seismic simulator is carried out through simulation. The program is developed by the MATLAB software. The parameters of the seismic simulator used for the simulation are given in Tables1, 2, 3, and 4.Table 1
The parameters of the base platform (m).
123456xBi0.4000000.000000−0.4000000.400000−0.400000−2.000000yBi−0.4000000.400000−0.400000−2.000000−2.0000000.000000zBi0.0000000.0000000.0000001.5000001.5000001.500000Table 2
The parameters of the moving platform which are measured in the coordinate frameO′-uvw (m).
123456xAi0.4000000.000000−0.4000000.400000−0.400000−0.681000yAi−0.4000000.400000−0.400000−0.681000−0.6810000.000000zAi−0.166000−0.166000−0.166000−0.037500−0.037500−0.037500Table 3
The length of the strutCiAi (m).
123456li1.0000001.0000001.0000001.0000001.0000001.000000Table 4
The mass parameters of the manipulator (kg).
123456mi202020202020mci50100505050100The mass of the moving platform ismp=200 kg. The inertia parameters used in the simulation are given as(52)Ipo′=[17.33300017.33300033.333]kg⋅m2,Iii=[1.2790001.2790000.005]kg⋅m2.
Other parameters used in the simulation are given asE=2.06×1011 Pa,  G=79.38×109 Pa, Es=2.06×1011 Pa, As=1.96×10-3 m2, L1i+L2i=1.1 m, ρ=7800 kg/m3, di=0.244 m, d=0.05 m, h1=2 m, h2=1.5 m, h=0.01 m, z0=1.744 m.


### 4.1. Natural Frequency

According to the vibration theory, the rigidity of the system may be represented by the natural frequency. The seismic simulator with the higher frequency would have the higher stiffness.From (48), we get(53)det⁡(-ω2M+K)=0,
where ω denotes the natural frequency. The distribution of the natural frequency is shown in Figure 5 when the pose of the moving platform is given as ϕx=ϕy=ϕz=0 and z=z0.Distributions of the natural frequencies in the workspace. (a) First-order natural frequency. (b) Second-order natural frequency.
(a)(b)It is shown in Figure5 that the second-order natural frequency is much higher than the first-order natural frequency.


### 4.2. Sensitivity Analysis

The sensitivity analysis is usually used to evaluate the effect of the structural design variables on the performance of the manipulator. From (48), we get(54)(-ωr2M+K)φr=0,
where φr and ωr are the mode shape value and the natural frequency of the vibration in the rth mode. Taking the derivative of (54) with respect to the structural design value pm such as the radius of the strut and the radius of the lead screw yields (55)(-2ωr∂ωr∂pmM-ωr2∂M∂pm+∂K∂pm)φr+(-ωr2M+K)∂φr∂pm=0.
Taking dot product of φr on both sides of the equation yields(56)φrT(-2ωr∂ωr∂pmM-ωr2∂M∂pm+∂K∂pm)φr+φrT(-ωr2M+K)∂φr∂pm=0.
Since (57)φrT(-ωr2M+K)=((-ωr2M+K)φr)T=0,φrTMφr=E.
give(58)-2ωr∂ωr∂pm-ωr2φrT∂M∂pmφr+φrT∂K∂pmφr=0.
so (59)∂ωr∂pm=-12ωr(ωr2φrT∂M∂pmφr-φrT∂K∂pmφr).
Figure 6 shows the sensitivity distribution of the manipulator when the pose of the moving platform is given as ϕx=ϕy=ϕz=0 and z=z0. It is shown that the first-order natural frequency is sensitive to the radius of the strut and the radius of the lead screw.Sensitivities of the first-order natural frequency to the structure parameters. (a) Radius of the struts. (b) Radius of the lead screws.
(a)(b)


### 4.3. Energy Ratio Distribution

The computation of the energy ratio is usually used to evaluate the allocation of the stiffness and the mass of the manipulator. Suppose thatTsr and Vsr are the maximum kinetic energy and elastic potential energy of the substructures vibrating in its rth mode. TAr and VAr denote the maximum kinetic energy and elastic potential energy of the system vibrating in the rth mode. Thus,(60)TAr=∑s=1NTsr,VAr=∑s=1NVsr,
where(61)Tsr=12ωr2AsrTmsAsr,Vsr=12AsrTksAsr.Asr is the oscillating amplitude array of the substructure vibrating in the rth mode. ms and ks denote the mass matrix and the stiffness matrix of the substructure, respectively.So the energy ratio of the substructure can be achieved as(62)TsrTAr=γsr,∑s=1Nγsr=1,VsrVAr=μsr,∑s=1Nμsr=1,
where γsr and μsr denote the kinetic energy ratio and the elastic potential energy ratio of the substructure, respectively. Figure 7 shows the distributions of the kinetic energy ratio and the elastic potential energy ratio, respectively, when the pose of the moving platform is given as ϕx=ϕy=ϕz=0 and z=z0. It is shown that the mass of the moving platform should be reduced or the stiffness of the strut should be increased in order to improve the dynamic characteristics of the manipulator and the stiffness of the sixth strut must be increased from the energy ratios computation.Distributions of the energy ratios in the workspace. (a) Kinetic energy ratios of the sliders.  (b) Kinetic energy ratio of the moving platform. (c) Kinetic energy ratios of the struts. (d) Elastic potential energy ratios of the sliders. (e) Elastic potential energy ratio of the moving platform. (f) Elastic potential energy ratios of the struts.
(a)(b)(c)(d)(e)(f)


### 4.4. Displacement Response Analysis

The displacement response analysis will be carried out by solving (48) subject to the initial conditions(63)U0=U(0),U̇0=U̇(0),Since the damping in the structure is a very complex subject [28], the modal damping ratios of ςr=0.1% are added to the Kineto-elastodynamic model of the manipulator.From (48), we get(64)(-ω2M+K)φ=0.
Neglecting higher-order terms, the displacement vector U of a multi-degree-of-freedom system can be expressed in terms of the four dominant modal contributions. Thus, the dynamic response of the system can be expressed as(65)U=φη,
where φ=[φ1φ2φ3φ4] is the modal matrix.Substituting (65) into (50) and adding the modal damping ratio yields(66)φTMφη̈+φTCφη̇+φTKφη=φTFd,
where(67)φTMφ=E4,φTCφ=C¯=diag⁡(2ς1ω12ς2ω22ς3ω32ς4ω4),φTKφ=Ω2=diag⁡(ω12ω22ω32ω42)E4 is the unit matrix of order four.Substituting (67) into (66) yields(68)η̈+C¯η̇+Ω2η=N,
where(69)N=φTFd.The stiffness matrix and the mass matrix of the Kineto-elastodynamic model of the parallel manipulator are time varying. The common strategy of solving this kind of problem is dividing the motion period into several small time internals and regarding the stiffness matrix and the mass matrix as constant in each small time interval [26].LetT denote the motion period which is divided into n intervals(70)Δt=Tn.
In the ith time interval (ti-1<t<ti), the motion equation of the manipulator is (71)η̈r+2ςrωr(i)η̇r+(ωr(i))2ηr=Nr,(r=1,2,…,N).
So the contribution of the rth mode to the displacement response is (72)ηr(t)=1ωdr(i)∫ti-1tNr(τ)e-ςrωr(i)(t-τ)sinωdr(i)(t-τ)dτ+ηr(ti-1)(1-ςr2)1/2e-ςrωr(i)(t-ti-1)cos⁡[ωdr(i)(t-ti-1)-ψr]+η̇r(ti-1)ωdr(i)e-ςrωr(i)(t-ti-1)sinωdr(i)(t-ti-1),(r=1,2,…,N),
where(73)ωdr(i)=(1-ςr2)1/2ωr(i),(74)ψr=arctanςr(1-ςr2)1/2.
Substituting t=ti into (73) yields(75)ηr(ti)=1ωdr(i)∫ti-1tiNr(τ)e-ςrωr(i)(ti-τ)sinωdr(i)(ti-τ)dτ+ηr(ti-1)(1-ςr2)1/2e-ςrωr(i)Δtcos⁡[ωdr(i)Δt-ψr]+η̇r(ti-1)ωdr(i)e-ςrωr(i)Δtsinωdr(i)Δt,(r=1,2,…,N).
Taking the derivative of (72) with respect to time and substituting t=ti into it yields
(76)η̇r(ti)=1ωdr(i)∫ti-1tiNr(τ)[-ςrωr(i)e-ςrωr(i)(ti-τ)sinωdr(i)(ti-τ)+ωdr(i)e-ςrωr(i)(ti-τ)cos⁡ωdr(i)(ti-τ)]dτ-ηr(ti-1)[ωdr(i)(1-ςr2)1/2e-ςrωr(i)Δtsin(ωdr(i)Δt-ψr)+ςrωr(i)(1-ςr2)1/2e-ςrωr(i)Δtcos⁡(ωdr(i)Δt-ψr)]-η̇r(ti-1)[ςrωr(i)ωdr(i)e-ςrωr(i)Δtsinωdr(i)Δt-e-ςrωr(i)Δtcos⁡ωdr(i)Δt](r=1,2,…,N).
It is shown from (75) and (76) that ηr(ti) and η̇r(ti) can be achieved when ηr(ti-1) and η̇r(ti-1) are given. As for ti=0,(77)ηr(0)=(φ(r))TMU(0),η̇r(0)=(φ(r))TMU̇(0).
So the total displacement response can be achieved by combining these modal contributions
(78)U(ti)=∑r=1Nηr(ti)φr(ti),
It is the sum of the steady-state response and the transient state response.Assuming that the investigated trajectory of the moving platform used in the simulation is expressed as(79)x=-0.1+amax⁡T22π(τ-12πsin(2πτ)),y=-0.1+amax⁡T22π(τ-12πsin(2πτ)),z=1.644+amax⁡T22π(τ-12πsin(2πτ)),ϕx=-0.1+amax⁡T22π(τ-12πsin(2πτ)),ϕy=-0.1+amax⁡T22π(τ-12πsin(2πτ)),ϕz=-0.1+amax⁡T22π(τ-12πsin(2πτ)),
where amax=9.8m/s2, τ=t/T, T=2πS/amax s in seconds, and S=0.2 m(rad). The motion period is divided into 512 intervals. The displacement response of the moving platform is shown in Figure 8. It is shown that the displacement response of the moving platform along the x direction is smaller than these displacement responses along the y direction and along the z direction. The angular displacement response of the moving platform rotating about z axis is slightly larger than those angular displacement responses rotating about the x axis and about the y axis.Displacement responses of the moving platform.(a)x direction,  (b) y direction,  (c) z direction, (d) ϕx direction, (e) ϕy direction, and  (f) ϕz direction.
(a)(b)(c)(d)(e)(f)


## 4.1. Natural Frequency

According to the vibration theory, the rigidity of the system may be represented by the natural frequency. The seismic simulator with the higher frequency would have the higher stiffness.From (48), we get(53)det⁡(-ω2M+K)=0,
where ω denotes the natural frequency. The distribution of the natural frequency is shown in Figure 5 when the pose of the moving platform is given as ϕx=ϕy=ϕz=0 and z=z0.Distributions of the natural frequencies in the workspace. (a) First-order natural frequency. (b) Second-order natural frequency.
(a)(b)It is shown in Figure5 that the second-order natural frequency is much higher than the first-order natural frequency.


## 4.2. Sensitivity Analysis

The sensitivity analysis is usually used to evaluate the effect of the structural design variables on the performance of the manipulator. From (48), we get(54)(-ωr2M+K)φr=0,
where φr and ωr are the mode shape value and the natural frequency of the vibration in the rth mode. Taking the derivative of (54) with respect to the structural design value pm such as the radius of the strut and the radius of the lead screw yields (55)(-2ωr∂ωr∂pmM-ωr2∂M∂pm+∂K∂pm)φr+(-ωr2M+K)∂φr∂pm=0.
Taking dot product of φr on both sides of the equation yields(56)φrT(-2ωr∂ωr∂pmM-ωr2∂M∂pm+∂K∂pm)φr+φrT(-ωr2M+K)∂φr∂pm=0.
Since (57)φrT(-ωr2M+K)=((-ωr2M+K)φr)T=0,φrTMφr=E.
give(58)-2ωr∂ωr∂pm-ωr2φrT∂M∂pmφr+φrT∂K∂pmφr=0.
so (59)∂ωr∂pm=-12ωr(ωr2φrT∂M∂pmφr-φrT∂K∂pmφr).
Figure 6 shows the sensitivity distribution of the manipulator when the pose of the moving platform is given as ϕx=ϕy=ϕz=0 and z=z0. It is shown that the first-order natural frequency is sensitive to the radius of the strut and the radius of the lead screw.Sensitivities of the first-order natural frequency to the structure parameters. (a) Radius of the struts. (b) Radius of the lead screws.
(a)(b)


## 4.3. Energy Ratio Distribution

The computation of the energy ratio is usually used to evaluate the allocation of the stiffness and the mass of the manipulator. Suppose thatTsr and Vsr are the maximum kinetic energy and elastic potential energy of the substructures vibrating in its rth mode. TAr and VAr denote the maximum kinetic energy and elastic potential energy of the system vibrating in the rth mode. Thus,(60)TAr=∑s=1NTsr,VAr=∑s=1NVsr,
where(61)Tsr=12ωr2AsrTmsAsr,Vsr=12AsrTksAsr.Asr is the oscillating amplitude array of the substructure vibrating in the rth mode. ms and ks denote the mass matrix and the stiffness matrix of the substructure, respectively.So the energy ratio of the substructure can be achieved as(62)TsrTAr=γsr,∑s=1Nγsr=1,VsrVAr=μsr,∑s=1Nμsr=1,
where γsr and μsr denote the kinetic energy ratio and the elastic potential energy ratio of the substructure, respectively. Figure 7 shows the distributions of the kinetic energy ratio and the elastic potential energy ratio, respectively, when the pose of the moving platform is given as ϕx=ϕy=ϕz=0 and z=z0. It is shown that the mass of the moving platform should be reduced or the stiffness of the strut should be increased in order to improve the dynamic characteristics of the manipulator and the stiffness of the sixth strut must be increased from the energy ratios computation.Distributions of the energy ratios in the workspace. (a) Kinetic energy ratios of the sliders.  (b) Kinetic energy ratio of the moving platform. (c) Kinetic energy ratios of the struts. (d) Elastic potential energy ratios of the sliders. (e) Elastic potential energy ratio of the moving platform. (f) Elastic potential energy ratios of the struts.
(a)(b)(c)(d)(e)(f)


## 4.4. Displacement Response Analysis

The displacement response analysis will be carried out by solving (48) subject to the initial conditions(63)U0=U(0),U̇0=U̇(0),Since the damping in the structure is a very complex subject [28], the modal damping ratios of ςr=0.1% are added to the Kineto-elastodynamic model of the manipulator.From (48), we get(64)(-ω2M+K)φ=0.
Neglecting higher-order terms, the displacement vector U of a multi-degree-of-freedom system can be expressed in terms of the four dominant modal contributions. Thus, the dynamic response of the system can be expressed as(65)U=φη,
where φ=[φ1φ2φ3φ4] is the modal matrix.Substituting (65) into (50) and adding the modal damping ratio yields(66)φTMφη̈+φTCφη̇+φTKφη=φTFd,
where(67)φTMφ=E4,φTCφ=C¯=diag⁡(2ς1ω12ς2ω22ς3ω32ς4ω4),φTKφ=Ω2=diag⁡(ω12ω22ω32ω42)E4 is the unit matrix of order four.Substituting (67) into (66) yields(68)η̈+C¯η̇+Ω2η=N,
where(69)N=φTFd.The stiffness matrix and the mass matrix of the Kineto-elastodynamic model of the parallel manipulator are time varying. The common strategy of solving this kind of problem is dividing the motion period into several small time internals and regarding the stiffness matrix and the mass matrix as constant in each small time interval [26].LetT denote the motion period which is divided into n intervals(70)Δt=Tn.
In the ith time interval (ti-1<t<ti), the motion equation of the manipulator is (71)η̈r+2ςrωr(i)η̇r+(ωr(i))2ηr=Nr,(r=1,2,…,N).
So the contribution of the rth mode to the displacement response is (72)ηr(t)=1ωdr(i)∫ti-1tNr(τ)e-ςrωr(i)(t-τ)sinωdr(i)(t-τ)dτ+ηr(ti-1)(1-ςr2)1/2e-ςrωr(i)(t-ti-1)cos⁡[ωdr(i)(t-ti-1)-ψr]+η̇r(ti-1)ωdr(i)e-ςrωr(i)(t-ti-1)sinωdr(i)(t-ti-1),(r=1,2,…,N),
where(73)ωdr(i)=(1-ςr2)1/2ωr(i),(74)ψr=arctanςr(1-ςr2)1/2.
Substituting t=ti into (73) yields(75)ηr(ti)=1ωdr(i)∫ti-1tiNr(τ)e-ςrωr(i)(ti-τ)sinωdr(i)(ti-τ)dτ+ηr(ti-1)(1-ςr2)1/2e-ςrωr(i)Δtcos⁡[ωdr(i)Δt-ψr]+η̇r(ti-1)ωdr(i)e-ςrωr(i)Δtsinωdr(i)Δt,(r=1,2,…,N).
Taking the derivative of (72) with respect to time and substituting t=ti into it yields
(76)η̇r(ti)=1ωdr(i)∫ti-1tiNr(τ)[-ςrωr(i)e-ςrωr(i)(ti-τ)sinωdr(i)(ti-τ)+ωdr(i)e-ςrωr(i)(ti-τ)cos⁡ωdr(i)(ti-τ)]dτ-ηr(ti-1)[ωdr(i)(1-ςr2)1/2e-ςrωr(i)Δtsin(ωdr(i)Δt-ψr)+ςrωr(i)(1-ςr2)1/2e-ςrωr(i)Δtcos⁡(ωdr(i)Δt-ψr)]-η̇r(ti-1)[ςrωr(i)ωdr(i)e-ςrωr(i)Δtsinωdr(i)Δt-e-ςrωr(i)Δtcos⁡ωdr(i)Δt](r=1,2,…,N).
It is shown from (75) and (76) that ηr(ti) and η̇r(ti) can be achieved when ηr(ti-1) and η̇r(ti-1) are given. As for ti=0,(77)ηr(0)=(φ(r))TMU(0),η̇r(0)=(φ(r))TMU̇(0).
So the total displacement response can be achieved by combining these modal contributions
(78)U(ti)=∑r=1Nηr(ti)φr(ti),
It is the sum of the steady-state response and the transient state response.Assuming that the investigated trajectory of the moving platform used in the simulation is expressed as(79)x=-0.1+amax⁡T22π(τ-12πsin(2πτ)),y=-0.1+amax⁡T22π(τ-12πsin(2πτ)),z=1.644+amax⁡T22π(τ-12πsin(2πτ)),ϕx=-0.1+amax⁡T22π(τ-12πsin(2πτ)),ϕy=-0.1+amax⁡T22π(τ-12πsin(2πτ)),ϕz=-0.1+amax⁡T22π(τ-12πsin(2πτ)),
where amax=9.8m/s2, τ=t/T, T=2πS/amax s in seconds, and S=0.2 m(rad). The motion period is divided into 512 intervals. The displacement response of the moving platform is shown in Figure 8. It is shown that the displacement response of the moving platform along the x direction is smaller than these displacement responses along the y direction and along the z direction. The angular displacement response of the moving platform rotating about z axis is slightly larger than those angular displacement responses rotating about the x axis and about the y axis.Displacement responses of the moving platform.(a)x direction,  (b) y direction,  (c) z direction, (d) ϕx direction, (e) ϕy direction, and  (f) ϕz direction.
(a)(b)(c)(d)(e)(f)


## 5. Conclusion

Based on the Kineto-elastodynamic assumption, the modeling and the Kineto-elastodynamic characteristics of the 6-PSS parallel structure seismic simulator have been systematically investigated through simulation. The conclusions are drawn from the simulation as follows.(1)
The maps of the natural frequencies with respect to the manipulator configuration have been achieved. It is shown that the second-order natural frequency is much higher than the first-order natural frequency.(2)
From the sensitivity analysis, the first-order natural frequency is sensitive to the radius of the strut and the radius of the lead screw.(3)
The mass of the moving platform should be reduced or the stiffness of the strut should be increased in order to improve the dynamic characteristic of the manipulator, and the stiffness of the sixth strut must be increased from the energy ratios computation.(4)
For the investigated trajectory, the displacement response of the moving platform along thex-direction is smaller than these displacement responses along the y direction and along the z direction. The angular displacement response of the moving platform rotating about z-axis is slightly larger than those angular displacement responses rotating about the x-axis and about the y-axis.



---

*Source: 489695-2011-09-07.xml*