3.1. Dynamics of the end-effector

In global frame OXYZ, the length vector of the ith cable L _i is defined as [Reference Peng, Wu, Lin, Zhou, Liu and Wang44]:

(1) \begin{equation} \boldsymbol{{L}}_{i}=\boldsymbol{{B}}_{i}-\boldsymbol{{X}}_{P}-\boldsymbol{{R}}\boldsymbol{{x}}_{{P_{i}}} \end{equation}

where $\boldsymbol{{B}}_{i}= \overrightarrow{OB_{i}}, \boldsymbol{{X}}_{P}= \overrightarrow{OP}=[X_{P},{Y}_{P},Z_{P}]^{{T}}$ represents the position vector of point P in global frame OXYZ; $\boldsymbol{{x}}_{P_{i}}$ is the position vector of point P _i in local frame Pxyz; and R is the transformation matrix from Pxyz to OXYZ.

As shown in Fig. 1, the end-effector is considered as a rigid body, and the equation of motion for the end-effector can be established based on the Newton–Euler method [Reference Peng, Wu, Lin, Zhou, Liu and Wang44]:

(2) \begin{equation} \boldsymbol{{M}}\ddot{\boldsymbol{{X}}}+\boldsymbol{{W}}_{C}-\boldsymbol{{W}}_{G}=\boldsymbol{{J}}_{A}^{T}\boldsymbol{{T}} \end{equation}

where $\boldsymbol{{M}}=\left[\begin{array}{c@{\quad}c} (m\boldsymbol{{I}})_{3\times 3} & \textbf{0}_{3\times 3}\\[3pt] \textbf{0}_{3\times 3} & \boldsymbol{{A}}_{G}\boldsymbol{{H}} \end{array}\right]$ denotes the mass matrix; $\boldsymbol{{H}}=\left[\begin{array}{c@{\quad}c@{\quad}c} \cos \theta \cos \psi & {-}\sin \psi & 0\\[3pt] \cos \theta \sin \psi & \cos \psi & 0\\[3pt] {-}\sin \theta & 0 & 1 \end{array}\right]$ ; m is the mass of the end-effector; A _G is the inertia matrix of the end-effector about its mass center (point P). W _C denotes the nonlinear Coriolis force, W _G denotes the gravity vector of the end-effector, and $\boldsymbol{{J}}_{A}^{{T}}$ denotes the Jacobian matrix of the CDPR. T is the tension vector, T = [T ₁, T ₂, …, T ₈]^T, and T _i (i = 1, 2, …, 8) represents the tension in the ith cable. The calculation expressions for the above variables can be found in ref. [Reference Peng, Wu, Lin, Zhou, Liu and Wang44].

Figure 3. Force analysis of the differential cable element.

3.2. Modeling of the cable

It is assumed that the cables are homogeneous. The cables are always in a tension state with the diameter d ≤ 3 mm. Therefore, the cables’ gravity is not taken into account. The bending stiffness of the cables is negligible, because the cables can sustain only tensile forces. Considering the elastic effect of the cables, the deformation of the cables conforms to Hooke’s law.

Consider the transverse vibration of the cable, since the sag of the cable is close to zero [Reference Irvine and Caughey45]. Taking one cable as an example, as shown in Fig. 3, force analysis of the differential cable element is carried out. The axial direction is along oy _wi , and the transverse direction is along ox _wi . An approximation ds _i $\approx$ dy _wi is made for the differential cable element. The cable tension is T _i u _i , and u _i is the unit directional vector of cable tension for the ith cable. The transverse component of the cable tension is T _xi , $\left| \boldsymbol{{T}}_{xi}\right| \approx -T_{i}\frac{\partial x_{wi}}{\partial y_{wi}}$ , and $\left| \boldsymbol{{T}}_{xi}+\Delta \boldsymbol{{T}}_{xi}\right| \approx T_{i}\dfrac{\partial x_{wi}}{\partial y_{wi}}+T_{i}\dfrac{\partial }{\partial y_{wi}}\dfrac{\partial x_{wi}}{\partial y_{wi}}dy_{wi}$ . The transverse motion equation of the ith cable can be obtained from the balance of forces:

(3) \begin{equation} \rho _{w}Ads_{i}\dot{v}_{xi}=T_{i}\frac{\partial ^{2}x_{wi}}{\partial y_{wi}^{2}}dy_{wi}+F_{xi} \end{equation}

where $\rho$ _w is the cable’s density, A is the cross-sectional area of the unstrained cable, ds _i is the arc length of the differential cable element, v _xi is the velocity in the direction of ox _wi , $\dot{v}_{xi}$ is the acceleration in the direction of ox _wi , and F _xi is the external force in the direction of ox _wi .

The axial speed of the ith cable (in the direction of oy _wi ), denoted by v _yi (t), is time dependent. At axial coordinate y _wi and time t, the transverse vibration of the cable is specified by the transverse displacement x _wi (y _wi , t) [Reference Chen, Tang and Zu36]. So $v_{xi}=\dfrac{dx_{wi}}{dt}$ . Then

(4) \begin{align} \dot{v}_{xi} & =\frac{\partial }{\partial t}\!\left(\frac{dx_{wi}}{dt}\right)+\frac{\partial }{\partial y_{wi}}\!\left(\frac{dx_{wi}}{dt}\right)\frac{dy_{wi}}{dt}\nonumber\\[3pt] & =\frac{\partial ^{2}x_{wi}}{\partial t^{2}}+2v_{yi}\frac{\partial ^{2}x_{wi}}{\partial y_{wi}\partial t}+\dot{v}_{yi}\frac{\partial x_{wi}}{\partial y_{wi}}+v_{yi}^{2}\frac{\partial ^{2}x_{wi}}{\partial y_{wi}^{2}} \end{align}

Consider the approximation of ds _i $\approx$ dy _wi . Suppose the external force F _xi is zero. By substituting Eq. (4) into Eq. (3), the partial differential equation for the transverse displacement of the ith cable can be derived [Reference Chen, Tang and Zu36]:

(5) \begin{equation} \rho _{w}A\!\left(\frac{\partial ^{2}x_{wi}}{\partial t^{2}}+2v_{yi}\frac{\partial ^{2}x_{wi}}{\partial y_{wi}\partial t}+\dot{v}_{yi}\frac{\partial x_{wi}}{\partial y_{wi}}+v_{yi}^{2}\frac{\partial ^{2}x_{wi}}{\partial y_{wi}^{2}}\right)=T_{i}\frac{\partial ^{2}x_{wi}}{\partial y_{wi}^{2}} \end{equation}

To obtain the numerical solution of Eq. (5), the Galerkin method is applied to discretize the spatial variable. Then the solution of Eq. (5) can be expressed as:

(6) \begin{equation} x_{wi}=X_{wi}\sin \!\left(\frac{\begin{array}{l} \pi \end{array}y_{wi}}{L_{i}}\right) \end{equation}

where X _wi is the function of time t, y _wi ∈[0, L _i ].

Substituting Eq. (6) into Eq. (5) gives

(7) \begin{align} & \rho _{w}A\!\left(\sin \!\left(\frac{\pi y_{wi}}{L_{i}}\right)\ddot{X}_{wi}+2v_{yi}\!\left(\frac{\pi }{L_{i}}\right)\cos \!\left(\frac{\pi y_{wi}}{L_{i}}\right)\dot{X}_{wi}+\dot{v}_{yi}\!\left(\frac{\pi }{L_{i}}\right)\cos \!\left(\frac{\pi y_{wi}}{L_{i}}\right)X_{wi}-v_{yi}^{2}\!\left(\frac{\pi }{L_{i}}\right)^{2}\sin \!\left(\frac{\pi y_{wi}}{L_{i}}\right)X_{wi}\right)\nonumber\\ & \qquad +T\!\left(\frac{\pi }{L_{i}}\right)^{2}\sin \!\left(\frac{\pi y_{wi}}{L_{i}}\right)X_{wi}=0 \end{align}

Setting $y_{wi}=\dfrac{L_{i}}{2}$ , the transverse vibration equation for the mid-span of the ith cable can be obtained:

(8) \begin{equation} \ddot{X}_{wi}-\left(\frac{\pi }{L_{i}}\right)^{2}\left(v_{yi}^{2}-\frac{T_{i}}{\rho A}\right)X_{wi}=0 \end{equation}

where the axial velocity of the cable v _yi is related to the motion of the end-effector, and the cable tension T _i is related to the motion of the end-effector and the cable vibration. The expression of these two variables is needed to be established to get the complete vibration equation of the cable.

3.2.1. Expression of v_yi

Based on Eq. (1), each side of it can be double as:

(9) \begin{equation} L_{i}^{2}=\left(\boldsymbol{{X}}_{P}+\boldsymbol{{R}}\boldsymbol{{x}}_{{P_{i}}}-\boldsymbol{{B}}_{i}\right)^{{T}}\left(\boldsymbol{{X}}_{P}+\boldsymbol{{R}}\boldsymbol{{x}}_{{P_{i}}}-\boldsymbol{{B}}_{i}\right) \end{equation}

Taking derivatives of both sides of Eq. (9) with respect to time t gives

(10) \begin{equation} L_{i}\dot{L}_{i}=\left(\boldsymbol{{L}}_{i}\right)^{{T}}\left(\dot{\boldsymbol{{X}}}_{P}+\dot{\boldsymbol{{R}}}\boldsymbol{{x}}_{{P_{i}}}\right) \end{equation}

Taking derivatives of the transformation matrix R with respect to time t gives

(11) \begin{equation} \dot{\boldsymbol{{R}}}=\boldsymbol{{R}}^{\circ }\boldsymbol{{R}} \end{equation}

where $\boldsymbol{{R}}^{\circ }=\left[\begin{array}{c@{\quad}c@{\quad}c} 0 & -\omega _{Z} & \omega _{Y}\\ \omega _{Z} & 0 & -\omega _{X}\\ -\omega _{Y} & \omega _{X} & 0 \end{array}\right]$ , $\boldsymbol{\omega}$ is the angular velocity vector of the end-effector, and $\boldsymbol{\omega }=\boldsymbol{{H}}[\dot{\varphi },\dot{\theta },\dot{\psi }]^{\textrm{T}}=[\omega _{X},\omega _{Y},\omega _{Z}]^{\textrm{T}}$ .

Substituting Eq. (11) into Eq. (10) gives

(12) \begin{equation} L_{i}\dot{L}_{i}=\left(\boldsymbol{{L}}_{i}\right)^{{\textrm{T}}}\dot{\boldsymbol{{X}}}_{P}+\left(\boldsymbol{{r}}_{i}\times \boldsymbol{{L}}_{i}\right)^{{\textrm{T}}}\boldsymbol{\omega } \end{equation}

The expression of the axial velocity of the ith cable can be derived from Eq. (12):

(13) \begin{equation} v_{yi}=\dot{L}_{i}=\boldsymbol{{u}}_{i}^{{\textrm{T}}}\dot{\boldsymbol{{X}}}_{P}+\left(\boldsymbol{{r}}_{i}\times \boldsymbol{{u}}_{i}\right)^{{\textrm{T}}}\boldsymbol{\omega } \end{equation}

3.2.2. Expression of T_i

The total tension T _i of the ith cable consists of three components: (1) the initial tension T _0i , (2) the variation of the cable tension T _ei caused by the 6-DOF motion of the end-effector, and (3) the variation of the cable tension T _δi caused by the vibration of the cable.

Using the definitions of L _0i and L _i listed in Section 2, the boundary elongation △L _i caused by the motion of the end-effector can be described as:

(14) \begin{equation} \Delta L_{i}=L_{i}-L_{0i}=\sqrt{\left(\boldsymbol{{B}}_{i}-\boldsymbol{{X}}_{P}-\boldsymbol{{R}}\boldsymbol{{x}}_{{P_{i}}}\right)^{{\textrm{T}}}\left(\boldsymbol{{B}}_{i}-\boldsymbol{{X}}_{P}-\boldsymbol{{R}}\boldsymbol{{x}}_{{P_{i}}}\right)}-L_{0i} \end{equation}

Due to $\partial x_{wi}/\partial y_{wi} \ll 1$ in the vibration of the cable that is in tension tightly, from Eq. (6) with an approximation made, the extension δ _i in the elastic deformation caused by the vibration of the ith cable can be derived:

(15) \begin{align} \delta _{i} & =\int \!\left(ds_{i}-dy_{wi}\right)=\int _{0}^{L_{i}}\!\left\{dy_{wi}\!\left[1+\left(\partial x_{wi}/\partial y_{wi}\right)^{2}\right]^{1/2}-dy_{wi}\right\}\nonumber\\ & \approx \int _{0}^{L_{i}}\frac{1}{2}\left(\partial x_{wi}/\partial y_{wi}\right)^{2}dy_{wi} =\frac{\pi ^{2}X_{wi}^{2}}{4L_{i}} \end{align}

Finally, the total tension T _i of the ith cable can be written as:

(16) \begin{equation} T_{i}=T_{0i}+k_{wi}\!\left(\Delta L_{i}+\frac{\pi ^{2}X_{wi}^{2}}{4L_{i}}\right) \end{equation}

where k _wi is the stiffness coefficient of the cable. As can be seen from Eq. (16), cable tension T _i varies with the motion of the end-effector and vibration of the cable.

3.2.3. Vibration equation of the cable

Based on the analysis above, considering the cable damping, the transverse vibration equation of the ith cable can be described as:

(17) \begin{equation} \ddot{X}_{wi}+c_{wi}\dot{X}_{wi}+\left(\frac{\pi }{L_{i}}\right)^{2}\left(\frac{T_{i}}{\rho _{w}A}-v_{yi}^{2}\right)X_{wi}=0 \end{equation}

where cable tension $T_{i}=T_{0i}+k_{wi}\!\left(\Delta L_{i}+\dfrac{\pi ^{2}X_{wi}^{2}}{4L_{i}}\right), v_{yi}=\boldsymbol{{u}}_{i}^{\textrm{T}}\dot{\boldsymbol{{X}}}_{P}+\left(\boldsymbol{{r}}_{i}\times \boldsymbol{{u}}_{i}\right)^{\textrm{T}}\boldsymbol{\omega }$ , c _wi is the damping coefficient.

3.3. Dynamic equation of CDPRs

Combining Eq. (2) and Eq. (17), and writing in matrix form leads to:

(18) \begin{equation} \left\{\begin{array}{l} \ddot{\boldsymbol{{X}}}_{w}+\boldsymbol{{C}}_{w}\dot{\boldsymbol{{X}}}_{w}+\boldsymbol{{K}}_{w}\boldsymbol{{X}}_{w}=0\\[5pt] \boldsymbol{{M}}\ddot{\boldsymbol{{X}}}+\boldsymbol{{W}}_{C}-\boldsymbol{{W}}_{G}=\boldsymbol{{J}}_{A}^{{T}}\boldsymbol{{T}} \end{array}\right. \end{equation}

where

X _w = [X _w1, X _w2,…, X _w8] ^T , X _wi (i = 1, 2, …, 8) represents the transverse deflection at the mid-span of the cable;

C _w = diag(c _w1, c _w2, …, c _w8), c _wi (i = 1, 2, …, 8) represents the damping coefficient of the ith cable;

K _w = diag(K _w1, K _w2, …, K _w8), $K_{wi}=\left(\dfrac{\pi }{L_{i}}\right)^{2}\left(\dfrac{T_{i}}{m_{s}}-v_{yi}^{2}\right)$ (i = 1, 2, …, 8);

T = [T ₁, T ₂,…, T ₈] ^T , T _i (i = 1, 2, …, 8) represents the tension in the ith cable, $T_{i}=T_{0i}+k_{wi}\!\left(\Delta L_{i}+\dfrac{\pi ^{2}X_{wi}^{2}}{4L_{i}}\right)$ ; and

$v_{yi}=\boldsymbol{{u}}_{i}^{\textrm{T}}\dot{\boldsymbol{{X}}}_{P}+(\boldsymbol{{r}}_{i}\times \boldsymbol{{u}}_{i})^{\textrm{T}}\boldsymbol{\omega }$ represents the axially moving velocity of the ith cable.

Table I. The parameters of the CDPR.

Equation (18) is a nonlinear coupled equation for the motion of CDPR, which is a multibody system having multiple DOFs. To obtain accurate approximate solutions at low computational cost, the numerical method is an effective solution. To solve Eq. (18), it can be transformed into a state equation as below:

(19) \begin{equation} \dot{\boldsymbol{{X}}}_{c}=\left[\begin{array}{l@{\quad}l} \begin{array}{l@{\quad}l} 0_{6\times 6} & \boldsymbol{{I}}_{6\times 6}\\[3pt] 0_{6\times 6} & 0_{6\times 6} \end{array} & \ \ \ \ \ 0_{12\times 16}\\[3pt]\ \ \ \ \ 0_{16\times 12} & \begin{array}{l@{\quad}l} 0_{8\times 8} & \boldsymbol{{I}}_{8\times 8}\\[3pt] -\boldsymbol{{K}}_{w} & -\boldsymbol{{C}}_{w} \end{array} \end{array}\right]\boldsymbol{{X}}_{c}+\left[\begin{array}{c} 0_{6\times 1}\\[3pt] \boldsymbol{{D}}^{*}\\[3pt] 0_{16\times 1} \end{array}\right] \end{equation}

where $\boldsymbol{{X}}_{c}=\left[\begin{array}{l} \boldsymbol{{X}}_{6\times 1}\\[3pt] \dot{\boldsymbol{{X}}}_{6\times 1}\\[3pt] \boldsymbol{{X}}_{{w_{8\times 1}}}\\[3pt] \dot{\boldsymbol{{X}}}_{{w_{8\times 1}}} \end{array}\right]$ ; $\boldsymbol{{D}}^{\ast }=\boldsymbol{{M}}^{-1}(\boldsymbol{{J}}_{A}^{T}\boldsymbol{{T}}-\boldsymbol{{W}}_{C}+\boldsymbol{{W}}_{G})$ .

The unknown X _c to be solved is a 28-by-1 vector, and Eq. (19) contains a large sparse matrix and the coupled terms D ^* and K _w , which add more difficulties to the solution of the equation. Therefore, Eq. (19) is solved by the implicit variable-order Runge–Kutta numerical integration method, which can solve the coupled problem effectively.

4.1. Simulation conditions

The simulation has been carried out on a 6-DOF CDPR suspended by 8 cables. The cables are Kevlar cables, which are made of aramid fiber material. The damping of the Kevlar cable with the diameter of 0.5 mm is very small and is ignored here, referring to ref. [Reference Peng, Wu, Lin, Zhou, Liu and Wang44]. Detailed coordinates of B _i and P _i are described in ref. [Reference Ji, Lin, Hu, Peng and Wang46]. A summary of the parameters used in the simulation is shown in Table I.

Note that X _w0i is the initial deflection at the mid-span of the ith cable (i = 1, 2, …, 8) perpendicular to its axis (Fig. 2), and X _P0 is the initial disturbance of the end-effector along the direction of OX. With the initial disturbance of X _P0 = 0.1 mm and X _w0i = 0.1 mm (i = 1, 2, …, 8), Eq. (19) is solved using MATLAB with a time step of 10⁻⁴ s.

One can assume that the motion of the end-effector is a pitching oscillation: $\theta =A^{*}\sin\! (2\pi ft)$ , A ^* = 2°. The variation interval of oscillation frequency f falls in [2.06, 20] Hz. The vibration of the cables excited by the oscillation of the end-effector and the resulting dynamic response of the CDPR are investigated.

Figure 4. Attitude amplitude of the end-effector: (a) displacement amplitude and (b) rotation angle amplitude.

Figure 5. Cable amplitude: (a) cable 1, cable 2, cable 6, and cable 7, (b) cable 3, cable 4, cable 5, and cable 8.

4.2. Simulation results and discussions

The frequency–amplitude relationship of the CDPR is demonstrated in Figs. 4 and 5. As shown in Figs. 4 and 5, the attitude amplitude variations of the end-effector are consistent with the cables. When the pitching oscillation frequency of the end-effector $f \lt 8$ Hz, the displacement amplitudes of the end-effector are less than 0.5 mm, the rotation angle amplitudes are less than 0.3°, and the cable amplitudes are only 0.1 mm. It can be seen that the CDPR is stable in the frequency interval of $f \lt 8$ Hz. When $f \gt 8$ Hz, both the end-effector and the cables begin to occur large vibration. Especially, at the pitching oscillation frequency of f = 10 Hz and f = 14 Hz, the amplitudes of both the end-effector and the cables reach the maximum. At f = 10 Hz, the maximum displacement amplitude of the end-effector is up to 33.10 mm, the maximum rotation angle amplitude is 5.04°, and the maximum cable amplitude is 95.87 mm. At f = 14 Hz, the maximum displacement amplitude of the end-effector is up to 32.19 mm, the maximum rotation angle amplitude is 7.76°, and the maximum cable amplitude is 127.30 mm. Consequently, two peaks can be found at f = 10 Hz and f = 14 Hz, according to Figs. 4 and 5. It can be seen that the CDPR is unstable when $f \gt 8$ Hz.

Figure 6. Pitching angle variation of the end-effector: f = 2.06 Hz.

Figure 7. Pitching angle variation of the end-effector: f = 14 Hz.

According to Figs. 4 and 5, it can be seen that obvious deviations in the attitude of the end-effector are caused by the vibrations of the cables and the end-effector, when the frequency is higher than 8 Hz. Therefore, it can be inferred that the reliable motion frequency for the CDPR is within 8 Hz. When the frequency is higher than 8 Hz, the influence of the flexible cables on the attitude deviation of the end-effector cannot be ignored. The flexibility of the cables needs to be considered in the motion control of the CDPR.

Besides, as shown in Figs. 4 and 5, it turns out that, as the oscillation frequency of the end-effector increases, the end-effector and the cables exhibit the dynamics process from steady state to unstable large-amplitude vibration and finally to stable small-amplitude vibration.

In order to further discuss the significance of the excitation effect caused by the axially moving length-variable cables, the results with and without excitation effect are compared. Equation (18) represents the case with excitation effect. In the case without excitation effect, we set X _w = 0 in Eq. (18), that is, the cable vibration and its excitation effect on the end-effector are not considered.

As shown in Fig. 6, when the pitching oscillation frequency f = 2.06 Hz, the pitching angle amplitude |θ |_max without excitation effect is 2.00°, while |θ |_max with excitation effect is 2.03°. The relative increment is only 1.5% when f = 2.06 Hz. As shown in Fig. 7, when f = 14 Hz, the pitching angle amplitude |θ |_max without excitation effect is 2.00°, while |θ |_max with excitation effect is 7.80°. The relative increment is 290% when f = 14 Hz. The resulting data are shown in Table II.

By comparing the results with and without excitation effect of the axially moving length-variable cables at different frequencies, it can be concluded that the excitation effect can cause large attitude deviations of the end-effector in some specific frequencies f. The large attitude deviation of the end-effector can reduce the positioning accuracy of the CDPR. The proposed dynamic model of the CDPR can be used to analyze the response of the end-effector, so as to evaluate the influence of excitation effect on the attitude accuracy of the CDPR, which is of great importance.

Table II. The results of the CDPR with and without excitation effect.

Moreover, the phase diagrams of vibrations of the end-effector and the cables at the frequency of 10 Hz and 14 Hz are shown in Figs. 8 and 9. The pitching angle of the end-effector and the amplitude of cable 5 are taken as examples. It can be seen that both the end-effector and the cables occur large vibrations at the frequency of 10 Hz and 14 Hz, and the vibrations are non-linear. This indicates that the dynamics of the CDPR exhibit non-linear characteristics, due to the influence of flexible cables.

Figure 8. Phase diagram of pitching angle variation of the end-effector: (a) f = 10 Hz, (b) f = 14 Hz.

Figure 9. Phase diagram of the cable (cable 5): (a) f = 10 Hz, (b) f = 14 Hz.

It can be concluded that attitude deviation of the end-effector is caused by the cable flexibility. This effect is particularly significant when the motion frequency of the end-effector is in some specific intervals. With the proposed nonlinear dynamic model of the CDPRs, these intervals can be quickly identified. It would be very useful for the design of the CDPRs.