Dynamics analysis of spatial parallel robot with rigid and flexible links

: To analyze the rigid–flexible coupling effects on the dynamic performance of a robot system, a dynamic model of a parallel robot with flexible spatial links is derived in detail using a floating frame of reference (FFR) formulation. Compared to the previous rigid–flexible coupling model where the kinematic chains are all flexible links or where the joints are all flexible components, the inertia matrix and the stiffness matrix are not constant matrix which leading to the differences in respect of dynamic performance in model. To verify the correctness of the derived dynamics equations, the dynamics solutions of the spatial parallel robot from an ideal rigid–body model and the FFR model containing rigid and flexible coordinates were established by an FFR formulation. Furthermore, a finite element analysis (FEA) model, which included rigid links and flexible spatial links, was constructed for comparison. The comparison of the three models showed that the trajectory trends were the same, but the motion trajectories of the end-effector obtained by the FFR and FEA models varied within a certain range, and the maximum variations occurred at the peaks of the trajectories. However, since the FFR model considered the coupling effects of rigid and flexible links and the micro-displacement of the end-effector, the amount of deformation was the largest.


Introduction
In multibody systems, the demand for lightweight designs requires researchers and designers to account for the impact of the link flexibility of the mechanisms. However, research on link flexibility has mainly focused on planar mechanisms and robotic arms, whereas there have been fewer studies focused on rigid flexible spatial parallel mechanisms with non-linear coupling effects. The parallel mechanism with rigid links and flexible spatial links is widely used in aerospace, industrial production, and other practical engineering fields. In particular, in the field of industrial production, robots with fewer degrees of freedom have been widely used [1][2][3][4][5]. Ideally, all the components should be set as rigid links for dynamics analysis, but at high speeds, the slender links will be elastically deformed, and the flexible deformation of the links and the coupling effects with the rigid links will have an important impact on the dynamic performance of the system. Therefore, it is of great significance to derive a reasonable dynamics model [6][7][8][9] of the rigid-flexible spatial parallel mechanisms with nonlinear coupling effects for analyzing the forces acting on a system and ensuring the trajectory accuracy of the end-effector.
The absolute nodal coordinate formulation (ANCF) [10][11][12][13] and the floating frame of reference (FFR) formulation [14,15] are widely used methods for establishing dynamics models of flexible multibody systems. The advantages of the former include the mass matrix being constant and the centrifugal and Coriolis forces being zero. However, the stiffness matrix is highly complicated and inefficient for use in calculating the dynamics solutions. In the FFR formulation, the coupling effects between the rigid and flexible coordinates are included in the dynamics [1,[16][17][18], which can be applied with modal reduction techniques. There are two sets of coordinate systems in the FFR formulation. The first set is used to define the configuration of the local coordinate system, whereas the second set is applied to describe the configuration of the deformed body. Lugrís et al. [19][20][21][22] provide detailed derivations for the inertia terms, shape integrals, and elastic forces of planar mechanisms based on the FFR formulation. For instance, the dynamic performance of the rigid-flexible coupling of a hub-beam was considered in Liu and Liu [23]. Long et al. [24] presented a new method for determining the dynamics of a parallel robot with a flexible platform. In Zhang et al. [25,26], flexible multibody dynamics models of planar parallel robots were established, and the coupling effects were analyzed. Liu et al. [27] verified the importance of thermoplastic coupling using a rigid-flexible-thermal coupling model and a simulation verification. In Liu et al. [28], Lagrange's equations were applied to analyze the dynamic characteristics of a spatial parallel manipulator with a rigid end-effector and three flexible kinematic chains. Han et al. [29] implemented dynamics analysis and a simulation of a flexible beam element using the ANCF and FFR methods. This paper makes six new contributions: 1) By setting the local coordinate system of the flexible spatial links as a non-centroid coordinate system, the influence of the coupling term in the inertia tensor on the dynamic performance of the system is considered.
2) Unlike previous models, this model uses the FFR method to determine the dynamics while accounting for the coupling effects between the rigid and flexible coordinates and the micro-displacement of the end effector. This allows the elastic vibrations and motion errors of the system to be more accurate.
3) Since the coupling term in the inertia tensor has an important influence on the dynamic performance of the system, it is considered in this analysis by setting the local coordinate system of the flexible spatial link as a non-centroid coordinate system. This coupling term was not taken into account in previous studies. 4) Three typical models were compared. The FFR model was established using the FFR formulation and Lagrange's equations, and simulation code was written in MATLAB. The FEA model was implemented using a co-simulation of SOLIDWORKS, ANSYS, and ADAMS. The ideal rigid-body model, which was established based on Lagrange's equations of the first kind, and an analytical solution were obtained.
5) The dynamics equation of the FFR model is a highly coupled, highly nonlinear, high-index differential-algebraic equation. The mass matrix is no longer a symmetric matrix, and the stiffness matrix is no longer a constant matrix. Converting the dynamics equation into a pure differential equation to obtain a numerical solution avoids the problem of numerical divergence caused by the inaccurate estimation of the initial value of the Lagrange multiplier and the new solution problem caused by the expansion of the model. 6) The results indicate that the FFR model is more accurate for calculating the influence of the flexible deformation on the dynamic performance of the spatial parallel robot. The recursive formulas are also applicable to other flexible spatial links.

FFR formulation of flexible spatial link
As depicted in Figure 1, the O-XYZ system is a global coordinate system that defines the large rigid-body displacement, and the o-xyz system is a local coordinate system that defines the small deformation body displacement. By using the FFR formulation, the absolute position vector of point k on the flexible link can be written as system is given by Eq (2). In particular, when the flexible link is a spatial link, matrix R includes more than one transformation angle ( 1, ) The finite element method was employed to discretize the deformation field of the flexible spatial link. The displacement vector of any point k on the flexible spatial link caused by the deformation can be expressed as where f q represents the node displacement vector of the spatial beam element in the local coordinate system, and N represents the shape function. The specific expression of f q is follows where o x , o y , o z , b x , b y , and b z represent deformation displacement vectors along the three coordinate axes at the nodes. ox  , oy  , oz  , bx  , by  , and bz  represent the rotational displacement vectors around the three coordinate axes at the nodes. It is assumed that the deformation displacement along the axis is described using a first-order polynomial, and a cubic polynomial interpolation function is used to describe the elastic angular deformation.
1) Axial interpolation function Assuming that the axial interpolation function is a linear function, it can be written as According to the boundary conditions 2) Horizontal interpolation function Assuming that the axial interpolation function is a linear function, it can be written as According to the boundary conditions The axial and horizontal displacements corresponding to nodes o and b can be obtained by solving Eq (6) and Eq (8), therefore, the shape function N can be expressed as

Spatial rigid-flexible coupling dynamics modeling
Lagrange's equations are typically used to derive the dynamics equations as follows [30] where the T is the kinetic energy, V represents the potential energy, C represents constraint equations, vectors q and Q are the generalized coordinates of the generalized external forces, and λ is a Lagrange multiplier vector.
where  is the density of the material. The global position vector k r is equal to The velocity vector of the flexible spatial link is computed from Eq (12), yielding the following equation The matrix H and generalized velocity vector q  can be expressed as where 0 I is a 3 × 3 identity matrix, and 1 2 1 , ,..., n I I I can be expressed as follows: By means of Eq (11), the translational kinetic energy of the flexible spatial link may be written as where A represents the cross-sectional area of the flexible spatial link. The mass matrix When the flexible spatial link degenerates into a rigid spatial link, the deformation coordinates of the flexible link are zero ( 0 f  q ). Hence, the translational kinetic energy of the rigid spatial link can be further written as 2) Rotational kinetic energy of flexible spatial link The rotational kinetic energy of the flexible spatial link is defined as where   , x t θ is the absolute rotation angle of the micro-segment with a distance x and a width dx from the local coordinate system origin in actual motion, and J , x I , y I , and z I represent the respective area moments of inertia of the flexible spatial link with respect to the local coordinate axes x , y and z the rotational kinetic energy of the flexible spatial link can be written as where f R  M is the mass matrix corresponding to the rotational kinetic energy and is expressed as The matrices , , , , When the flexible spatial link degenerates into a rigid spatial link, the deformation coordinates of the flexible link are zero ( 0 f  q ). Hence, the rotational kinetic energy of the rigid spatial link can be further written as 3) Total kinetic energy of flexible spatial link The total kinetic energy of the flexible spatial link is calculated as follows: The total kinetic energy of the rigid spatial link is calculated as follows: However, the generalized coordinates of the flexible spatial link are expressed in the o xyz  system, and they must be expressed in the O XYZ  system. Therefore, the transformation relationships from f q to m u are expressed as follows where B is given by Substituting Eqs (28) and (29) into Eq (26), the total kinetic energy of the flexible spatial link in The specific forms of i q  and O f  M are shown in Appendix 1. Since the mass matrix in the o xyz  system is a symmetric matrix, but the total mass matrix of the system must be converted to the O XYZ  system for representation, so it needs to be multiplied by the transform matrix B , and because the matrix B is a asymmetry matrix. Therefore, the system mass matrix is no longer a symmetric matrix.

Total kinetic energy of each kinematic chain
By assembling the differential equations of the rigid and flexible links, the kinetic energy of each kinematic chain can be expressed as where q represents the generalized coordinates of the kinematic chain, and i M represents the mass matrix of the kinematic chain.

Generalized elastic forces of each kinematic chain
According to the virtual work principle, the elastic force virtual work of any beam element in the flexible spatial link can be written as where  and  represent stress and strain vectors, respectively. For linear isotropic materials, into Eq (32), the beam element stiffness matrix can be written as where E and D represent the partial differential operator of the strain-displacement and the shear modulus of the flexible spatial link, respectively. The specific expression of f k is shown in Appendix 1.
Because the rigid links cannot produce generalized elastic forces, the stiffness matrix of the kinematic chain in the o xyz  system can be expressed as Since f k is a constant matrix, i k is also a constant matrix, but the total stiffness matrix of the system must be converted to the O XYZ  system for representation, so it needs to be multiplied by the transform matrix B , and because the matrix B is a time-varying related to generalized coordinates. Therefore, the system stiffness matrix is no longer a constant matrix, and its specific form i K is shown in Appendix 1.
The generalized elastic forces of the kinematic chain may be expressed by deriving the elastic potential energy relative to the elastic coordinates where i e V represent the elastic potential energy of the kinematic chain i.

Generalized external forces of each kinematic chain
Based upon the aforementioned generalized forces, the generalized forces of the kinematic chain due to gravity of the spatial parallel robot can be written as The virtual work of the robotic actuated force (torque) can be written as The robotic actuated force (torque) of the kinematic chain can be derived by differentiation as By means of Eq (10), the matrix form of the dynamics equations of the kinematic chain can be obtained as follows where i V Q includes the Coriolis and centrifugal force matrices of the kinematic chain, which can be written as

Dynamics equations of spatial rigid flexible coupling model
The dynamics equations of the spatial rigid-flexible coupling model can be expressed by assembling the dynamics equations of the kinematic chain and the constraint equations as follows where M and K are the mass and stiffness matrices of the spatial rigid-flexible coupling model, respectively, ˆg F and ˆF Q are the generalized gravitational and robotic actuated force (torque), respectively, ˆV Q represents the Coriolis and centrifugal forces, q and T q C are generalized coordinates and constraint Jacobian matrix of the spatial rigid-flexible coupling model, respectively. Among them, the constrained Jacobian matrix T q C can be obtained by deriving the generalized coordinates from the constrained equations. The matrices M , K and q can be expressed as where n and m represent the number of the rigid and flexible coordinates of the spatial parallel robot, respectively.

Physical description
A typical 3-RRRU spatial parallel robot with rigid links and flexible spatial links is shown in Figure 2. This robot has three triangular symmetric kinematic chains. It is composed of ten parts, a driving link A i B i (i = 1,2,3), an intermediate link B i C i (i = 1,2,3), a passive link C i P i (i = 1,2,3), and the end-effector P 1 P 2 P 3 . The fixed base and driving links are connected by a revolute pair A i (i = 1,2,3), whereas joints B i (i = 1,2,3) and C i (i = 1,2,3) are passive revolute pairs. The end-effector and passive links are connected by universal joints. Among these ten parts, the passive link C i P i (i = 1,2,3) is a slender link that will be elastically deformed under high-speed motion. Therefore, it is treated as a flexible-body, and the other parts are treated as a rigid body. The link C i P i (i = 1,2,3) can perform arbitrary movements in space. However, the component connected to it is a rigid component, so the end-effector can only perform translational movement along three coordinate axes. End-effector The joint variables of the robotic kinematic chain with rigid links and flexible spatial links is illustrated in Figure 3. system is the Y-axis, but the rotation axis of passive link CiPi (i = 1,2,3) is the Z-axis. Therefore, the passive link CiPi (i = 1,2,3) is a spatial link.
in the Z direction, and Ai Z is parallel to Z .

The dynamic analysis and comparison of the three robotic models
To ensure the correctness of the derived dynamic equations, it is necessary to verify the influence of the flexible deformation on the system by comparing the end-effector motion trajectories in the FFR and the ideal rigid-body models and to verify the accuracy of the derived formula by comparing the motion trajectory deviations of the FFR and FEA models relative to that of the ideal rigid-body model. The dynamics equations of the three models were established using Lagrange's equations of the first kind and solved by numerical integration.

FFR model of robotic dynamics
1) Dynamic equations of the spatial flexible link Dynamics models can be established using the dynamics equations derived in section 2 and 3. First, by means of Eq (12), the vector k r is described as follows: where The absolute rotation angle vector can be written as  Finally, the translational and rotational kinetic energies of the flexible spatial link of the 3-RRRU spatial parallel robot can be obtained by substituting Eqs (43) and (44) into Eqs (11) and (20), respectively, and the total kinetic energy of the flexible spatial link in the global coordinate system can be obtained according to Eq (30), as follows: where   2) Dynamic equations of rigid components Because the driving link A i B i (i = 1,2,3) is purely rotated, its rotational kinetic energy is The intermediate link B i C i (i = 1,2,3) has both translational and rotational kinetic energy, and its kinetic energy is expressed as follows 2 2 2 1 1 2 1 2 1 2 2 2 1 1 The kinetic energy of each kinematic chain can be obtained by substituting Eqs (45)-(47) into Eq (31).
The potential energy of rigid component of the kinematic chain can be written as The generalized gravitational force i g F of the kinematic chain can be obtained by substituting Eq (48) into Eq (36). The generalized coordinates and the driving force of the C i P i (i = 1,2,3) can be expressed as The actuated torque i F Q of the robotic system can be obtained by substituting Eq (49) into Eq (38). However, the elastic deformation of the flexible links will cause positional errors in the kinetic trajectory of the end-effector, Therefore, the micro-displacement and micro-rotation generalized coordinates of the end-effector can be assumed to be . P U represents geometrical constraint relationship [31] between the actual kinetic trajectory and ideal kinetic trajectory of the end-effector. Thus, the kinetic and potential energy of the end-effector can be written as

3) Constraint relationship
The constraint equations are defined as follows.
(1) Constraint equation of the length of link B i C i (i = 1,2,3) (2) Constraint equation for the distance between point P i (i = 1,2,3) and P (3) Constraint equation for joint C i (i = 1,2,3) can be expressed by the position vector of the rigid joint C i (i = 1,2,3) equal to the position vector of the flexible joint C i (i = 1,2,3) (4) Constraint equation for joint P i (i = 1,2,3) can be expressed by the position vector of the rigid joint P i (i = 1,2,3) equal to the position vector of the flexible joint P i (i = 1,2,3) The constraint equations can be written as The velocity constraint matrix and the angular velocity constraint matrix can be obtained by calculating the first and second derivative of the position constraint matrix Ĉ with respect to time t, the specific expression are as follows The Eq (56), , , and Eq (51) are substituted into Eq (41) to obtain the matrix form of the system dynamics equations of the spatial rigid-flexible coupling model.

4) Numerical solution
The FFR model of the spatial parallel robot with rigid links and flexible spatial links was established using the FFR formulation. The number of generalized coordinates for the FFR model was 42, which included 12 rigid and 30 flexible coordinates. By means of Eq (41), the dynamics equations can be rewritten as Since the Eq (59) in this study were differential algebraic equations with a high index of 3, and numerical divergence occurred readily. Therefore, to obtain the solution, they were combined with the acceleration Eq (58) to solve.
Because the dynamics solution was obtained at the acceleration level, as the number of iterations increased, it was possible to calculate the position and velocity constraint equation errors. Therefore, positional and velocity violations were eliminated using the Baumgarte stabilization method, which can be expressed as follows where  and  are system feedback control parameters, 1 ε and 2 ε represent the displacement and velocity constraint stabilization parameters, respectively, expressed as follows It is necessary to estimate the initial values of the generalized position coordinate vector q and the Lagrange multiplier vector λ . q can be estimated based on the initial position of the rigid-body dynamics, but it is much more difficult for λ . Therefore, the numerical calculation was performed by reducing the index to 1. Equation (59) can be rewritten in the following form Multiplying both sides of Eq (62) by q C yields the following Vector λ can be obtained by combining Eqs (60) where h is the step size. Computational programs were written in MATLAB to perform dynamic simulation of the FFR model. 5) Frequency solution By means of Eq (41), the free vibration equations of the 3-RRRU spatial parallel robot with rigid links and flexible spatial links can be written as ˆ  Mq Kq 0  (65) According to the modal analysis theory [32], the natural frequency and modal expressions of the system can be obtained as When the mass matrix is a non-singular matrix, let  -1 S M K , then Eq (66) can be rewritten as where m  is the natural angular frequency of the system, and m f is the natural frequency of the system.

Ideal rigid-body model of robotic dynamics
Ideally, the spatial links C i P i (i = 1,2,3) of the spatial parallel robot can be treated as a rigid body, and the kinetic energy i i C P T can be directly calculated using Eq (27). The specific expression is 3 1 The potential energy of the spatial link C i P i (i = 1,2,3) can be obtained by the work-energy theorem, and its expression is as follows: The kinetic and potential energies of the end-effector can be calculated as Therefore, the dynamic model of the ideal rigid-body can be established by introducing the kinetic and potential energy of each rigid link into the Lagrange equations of the first kind, which can be solved by numerical methods.

FEA model of robotic dynamics
To compare and verify with the FFR model, the structural features and numerical solution methods of the FEA model are the same as the FFR model. When the spatial parallel mechanism is simulated by finite element software, there are many factors affecting its dynamic characteristics. Therefore, the following assumptions must be made before establishing a FEA model of a spatial rigid-flexible coupled parallel robot.
1) The elastic deformation of the spatial slender link during system operation will have an important impact on the dynamic performance of the system. Therefore, in the FEA model, the spatial link is considered to be a flexible body, and the other links are considered to be rigid bodies for analysis.
2) Each pair of components are connected by a rigid motion pair.
3) There is negligible friction between pairs. The FEA model was implemented using a co-simulation of three software packages: ANSYS, SOLIDWORKS, and ADAMS. The workflow of the co-simulation using these three software packages is shown in Figure 6.  Table 1. The dynamics equations were deduced using Lagrange's equations, and an integral solver (SI2, GSTIFF) was used for the simulation analysis in ADAMS, and the solution accuracy is 10 -4 . Table 1. Material properties of robot.
The simulation time was set to 5 s, and the step size was 0.001 s. To verify the correctness of the rigid-flexible coupled model, the results of the fully rigid inverse dynamics model were loaded into the rigid-flexible coupled model based on the same driver. The validity of the model was determined by examining the movement trends and trajectory errors of the two types of models.

Simulation results and analysis
To validate the presented method, a circular trajectory of the geometric center point of the end-effector was used to carry out the simulation. The circle defined in the plane 0.7 m z P   with a radius of 0.1 m can be written as where  represents the angular velocity of the end-effector, which was set to 2 rad/s. The shearing modulus G = 8 × 10 4 MPa. The geometric properties of all the parts are shown in Table 2, and the solution accuracy is 10 -4 . When the dynamic analysis was performed on an ideal rigid-body model, the material properties of all parts were steel. Table 2. Geometric properties of robot.
The relationship of the natural frequency of the spatial parallel robot with rigid links and flexible spatial links with time can be solved by MATLAB software programming, and the result is shown in Figure 7. In order to verify the correctness of natural frequency by using the derived dynamics model, based on the FEA simulation model, the natural frequency of the system is analyzed through the ADAMS/Vibration module. The results are shown in Table 3.
It can be seen from Table 3 that the first 6 natural frequencies in the vibration simulation model are all 0, that is, the rigid body mode corresponding to the system. Therefore, the minimum natural frequency of the system appears at the 7th order, and the corresponding natural frequency is 2.21E-005Hz. The minimum natural frequency of the derived dynamics model calculated by Eq (67) is 3.72E-005Hz. Due to the large number of elements of the flexible spatial link and the system's assembly errors, the minimum natural frequency of the simulation model is smaller than the derived dynamics model, but the two values are basically the same, thus verifying the correctness of the establishment of the dynamics model. A comparison of the displacement of the end-effector between the ideal rigid-body and FFR models is shown in the figures below. The displacement trajectories of the geometric center point of the end-effector in the X and Y directions with time are shown in Figures 8 and 9, respectively.  As shown in Figures 7 and 8, the trajectory trends of the end-effector of the ideal rigid-body and FFR models were basically the same, but the FFR model predictions yielded certain range of fluctuations relative to the ideal rigid-body model trajectory. The maximum deviations of the fluctuations are shown in Table 4. The fluctuations reached a maximum at the peak point of the motion trajectory. Table 4. Maximum deviations of the end-effector between the Ideal rigid-body and FFR models.
The displacement of the geometric center of the end-effector in the Z direction is shown in Figure 9. As shown in Figure 10, when the end-effector was moving in a plane, the Z-direction displacement of the rigid system did not change, and the rigid flexible coupling system floated between −0.70083 and −0.6993 m. The circular motion trajectories of the end-effector of the ideal rigid-body and FFR models on the three planes are shown in Figure 11. Figure 11 shows that the motion trajectory of the ideal rigid-body model was an ideal circular curve, but the FFR model was no longer an ideal circular curve. This was because the elastic deformation of the flexible spatial links in the FFR model affected the dynamic performance of the system. Therefore, the elastic deformation of the flexible spatial links cannot be ignored, and it is necessary to establish a dynamic model that follows rigid-flexible coupling dynamics.  As shown in Figures 12 and 13, the FFR and FEA models exhibited the same trends and fluctuated within a certain range. This was because the FFR and FEA models considered the elastic deformation of the flexible spatial links. The maximum fluctuations of both models in the X and Y directions mainly appeared at the peaks of the trajectories. However, compared to the FEA model, the ranges of the FFR model's values in the X and Y directions were −0.008942 to 0.0022 m and −0.00471 to 0.000805 m, respectively.
The displacements of the end-effector for the FFR and FEA models in the Z direction are shown in Figure 14.
As shown in Figure 14, the motion trajectories of the FFR and the ideal rigid-body models both fluctuated within a certain range, and the range of the FFR model was larger than that of the ideal rigid-body model. The floating values in the Z direction were −0.00099 to 0.00125 m. As shown in Figure 15, the deviations of the end-effector in the X, Y, and Z directions showed periodic variations. The error fluctuations in the X direction were the largest, followed by those in the Y direction, and the error changes in the Z axis were the smallest. This was related to the relative position of the geometric center of the end-effector and the global coordinate system. The comparison of the maximum deviations of the motion trajectories of the end-effector of the FFR and FEA models in the X, Y, and Z directions relative to the ideal rigid-body model are shown in Figure 16.
As shown in Figure 16, the maximum deviations of the motion trajectory of the end-effector of the FFR model relative to the ideal rigid-body in the X, Y, and Z direction were larger than those of the FEA model. The deviations of the trajectories between the FFR and FEA models relative to the ideal rigid-body in the X, Y, and Z directions were 0.00372, 0.00015, and 0.000066 m, respectively. This was because the FFR model accounted for the micro-displacement of the end-effector in the dynamic modeling process and strengthened the influence of the flexible spatial link deformation on the dynamic performance of the system by ignoring damping. In addition, the FFR model further considered the nonlinear coupling effects of the rigid links and flexible spatial links. Therefore, the FFR model could more accurately describe the nonlinear deformation of the flexible links and the robotic rigid-flexible coupling effects during the motion, resulting in fluctuations in the trajectory compared to the trajectories of the other two methods.

Conclusions
Compared with an ideal rigid-body model, the FFR and FEA models exhibited a range of fluctuations in the X, Y, and Z directions because they considered flexible deformation of the spatial links, and the results were consistent with actual observations. The peak values of the end-effector were larger for the FFR than those obtained using the FEA model. This occurred because the FFR model considered not only the deformation of the spatial flexible links, but also the nonlinear coupling effects between the rigid and flexible spatial links and the micro-displacement of the end-effector. The FFR model could more accurately describe the influence of the flexible links on the spatial parallel robot. Because the effectiveness of the dynamics equations was verified, this method can be used for the dynamics modeling of other flexible spatial links where the nonlinear coupling effects cannot be neglected.
where G represents the shear modulus of the flexible spatial link, ©2020 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)