## Abstract

A major advantage of the screw theory is that translations and rotations are treated simultaneously, which can provide greater insight into the vibration phenomena, such as vibration centers and axes. The present study describes how these concepts are extended into beam theory. The stiffness matrix of a beam was derived by incorporating different types of vibrations, such as extension, compression, torsion, and bending, at the same time, which was then used to obtain the equations of motion, including nonstandard forms of boundary conditions. We then presented an analytical method to solve these equations by focusing on two distinct examples, namely the cantilever and robot link. In the first numerical example, the mode shapes of the beam could be regarded as rotations about the vibration centers or axes of the rigid bodies in a discrete system. In the second example, the analytical solutions of mode shapes and natural frequencies of a robot link, for which the revolute joints at both sides are not parallel, were presented to demonstrate the utility of the screw theory. We demonstrated that the screw approach could accurately describe the vibrations of both discrete and continuous systems and that the geometric meaning of the vibration modes of discrete systems can be extended into continuous systems.

## 1 Introduction

The screw theory is often utilized to interpret spatial mechanisms and robotics issues at a high level. The major advantage of this theory is that translations and rotations are considered simultaneously, thus alleviating the complexity of these systems when implementing algebraic and numerical options. From this point of view, the screw theory may be useful for describing vibrations of beams with nonstandard boundary conditions. Previously, Selig and Ding [1,2] applied this theory to derive the equations of static and dynamic beams and identified the natural frequencies and mode shapes for Timoshenko beams with revolute joints at both ends that were not parallel. However, their findings were obtained by using numerical solutions as opposed to implementing an analytical method.

Several research studies have been performed using the screw theory. Ball [3] was the first to describe the oscillation of a rigid body as a repetitive small screw motion. The screw representation has also been used to analyze and design vibration systems. For example, Blanchet [4] introduced the concept of vibration centers and interpreted the mode shape of a planar three-degree-of-freedom (DOF) system as a pure repetitive rotation about the vibration center. That study also showed that the orthocenter of a triangle formed by three vibration centers coincides with its mass center. Some research studies [5–7] utilized transparent geometric constraints to develop design methods for energy harvesters, which focused on uniform power peaks at equally spaced resonant frequencies to widen the effective bandwidth. However, the nature of spatial vibration systems is so complicated that this approach is often impractical for these systems. One of the ways to alleviate the complexity of spatial vibration systems is to employ the decoupling technique. For instance, decoupling with respect to the plane of symmetry was utilized to design engine mounts, absorbers, and energy harvesters [8–10]. In addition to decoupling methods, Lee et al. [11] found three geometric constraints that could be applied to practical design methods. One of the constraints was used to design a spatial vibration system with prespecified ratios of peaks and target frequencies [12]. Furthermore, Ryu et al. [13] investigated the geometrical properties of a dual-body system and developed a design method for prespecified ratios among six energy peaks using a planar symmetric dual-body system.

Various analytical methods for nonstandard beams were previously proposed. A few studies have investigated the dynamic behavior of a damaged beam [14–18]. For example, Čas et al. [19] considered a three-dimensional two-layer composite beam with interlayer slips. Moreover, Liu and Yang [20] presented an analytical method to describe elastically connected double-beam systems, while Rezaiee-Pajand et al. [21] derived the exact elemental stiffness matrix and calculated an accurate solution for a nonprismatic beam–column.

In the present work, we used the screw approach to describe the vibrations of discrete and continuous systems. Using our approach, analytical solutions for nonstandard beams can be obtained, and our findings reveal that the geometric meaning of the vibration modes of discrete systems can be extended into continuous systems. This paper is organized as follows: Sec. 2 provides a geometrical interpretation of the vibrations of a rigid body. In Sec. 3, the stiffness matrix of the beam expressed in terms of the screw theory is introduced, and the equation of motion (EOM) governing the vibrations of the beams is derived, including the boundary conditions. Analytical solutions of the differential equations for a cantilever are then presented. Section 4 describes the discretization of a continuous system into an equivalent multibody system with nonstandard boundary conditions. In Sec. 5, we demonstrate the geometric meaning of vibration modes by comparing the results of continuous and discrete analysis, and the robot link is selected as the second numerical example to demonstrate the utility of the screw theory. Finally, Sec. 6 presents our conclusions and future research directions.

## 2 Screw Theory of Vibrations

In this section, we describe the screw theory of vibrations. Using the expression, the geometrical properties of single and multibody systems, such as vibration centers and axes, are described via the screw theory.

### 2.1 Single Body.

*k*

_{i}is the

*i*th spring constant, and $s^i\u2208R6\xd71$ is a line vector indicating the axis of the

*i*th spring. The line vector is represented by $s^i=[siT;riT\xd7siT]T$, where $si$ is a unit direction vector, and $ri$ is a position vector to a point on the axis of the spring.

*ϕ*( ∈

*R*

^{3×1}) is a small angular vector. Substituting Eq. (3) in Eq. (1), we can obtain that

Equation (4) yields six eigenvectors $D^i(i=1,\u2026,6)$, which correspond to the vibration modes. The vibration modes are generally determined as screws, meaning that the general mode shape is reflected by screw motion in the spatial system. There are three vibration modes in the *xy* planar system, as shown in Fig. 1: $D^i=[Yi\u2212Xi0001]T$ for *i* = 1, 2, 3. In this case, the operational mode shape is rotation about the vibration axis $D^i$ (or vibration center (*X*_{i}, *Y*_{i}, 0)).

### 2.2 Multibody System.

*i*th body, and $K12$ is one of the springs connected to each other. $D(=[D1TD2T])$ is the 12 × 1 displacement vector, which can be rewritten in terms of a screw chain as [13]

From Eq. (7), we can generally obtain 12 eigenvectors (or mode shapes) $D~i(i=1,2,\u2026,12)$. Considering that $D~$ is the screw chain, it can be said that there are 24 screw axes of vibration, and the mode shape can be described by a general vibrational screw motion. In a planar system, there exist 12 vibration centers, and the mode shapes are rotational motions [13].

*n*-body vibration system. To obtain the EOM for undamped free vibration, $M$ and $K$ in Eq. (5) become

*n*screws as $=D~ej\omega t$. Using the screw chain, the same form of Eq. (7) is obtained, and thus 6

*n*mode shapes are determined, each consisting of

*n*axes of the screw. In terms of engineering, a system comprising stacked rigid bodies, such as the beam shown in Fig. 2, may be more meaningful than that where all bodies are directly connected to each other. Since only the adjacent bodies are connected by springs, the stiffness matrix can be given by

If we consider *n* → ∞, then there are infinitely many natural frequencies and mode shapes in a continuous system such as the cantilever. If the elastic properties of the beam are reflected in the stiffness matrix described by Eq. (9), we may acquire mode shapes similar to those of the beam. Our goal is to show that the dynamic characteristics of this system become equivalent to the ones of the beam as *n* → ∞ (i.e., each body becomes a plate since the length of the beam remains constant). In the following section, we derive the differential equations of the beam using the screw and solve the equations analytically.

## 3 Vibration of a Continuous System

This section initially derives the stiffness matrix of a beam. Using this stiffness matrix, the EOM for the vibration of a continuous system is obtained and solved analytically.

### 3.1 Stiffness Matrix of a Beam.

A uniform prismatic beam with length *L*, area of cross section *A*, Young's modulus *E*, shear modulus *G*, and density *ρ* is modeled by a stack of *n* identical pieces assumed to be rigid. To track these pieces, we set global and local coordinate frames as follows (Fig. 3):

The global frame (

*O*−*XYZ*) is placed at the foot of the beam, the origin*O*coincides with the centroid of the cross section, and the*X*- and*Y*-axes of the frame are aligned with the principal directions of the cross section.The local frame (

*A*_{i}−*XYZ*) is placed in each rigid body, which is aligned with the principal axes of the rigid body.The local frame (

*B*_{i}−*XYZ*) is located in each cut plane, the origin*B*_{i}coincides with the centroid of the cut plane, and the*X*- and*Y*-axes are aligned with the principal directions.

Under the assumption that the beam is uniform, it is clear that the *Z*-axes of all frames are aligned as shown in Fig. 3.

*p*sets of three orthogonal springs (Fig. 4). First, in order to determine the spring constants, we suppose that the virtual springs parallel to the

*Z*-axis are stretched by Δ

*L*/

*n*, so the

*n*th rigid body is displaced by Δ

*L*when an axial load

*F*

_{Z}is applied, as shown in Fig. 5. In accordance with the Hooke's law, we have

*L*can be expressed as

*F*

_{X}. Similarly, the shear force can be obtained from the shear stress–strain relationship as follows:

*δ*is the sliding displacement, and thus,

*δ*can be given by

*k*

_{Y}can be determined as

*i*th cut plane can be obtained as follows:

*B*

_{i}indicates that the matrix is expressed at

*B*

_{i}−

*XYZ*. Referring to Fig. 4, $s^X,q$, $s^Y,q$, and $s^Z,q$, which are orthogonal three springs illustrated enlarged, are defined as follows:

*p*→ ∞ (and thus Δ

*A*→ 0), $KBi$ becomes

$KBi$ is clearly diagonal since $\u222bAXdA=\u222bAYdA=0$.

*B*

_{i}−

*XYZ*to

*O*−

*XYZ*is given by

*Z*is the distance from

*O*−

*XYZ*to

*B*

_{i}−

*XYZ*. Substituting Eqs. (20) and (22) into Eq. (21), the stiffness matrix can be obtained as

*Z*. It can be said that $K(Z)$ is the stiffness matrix of the virtual springs on the cross section at a distance

*Z*from the origin

*O*.

### 3.2 EOM for Undamped Free Vibration.

*i*th piece, respectively. The mass matrix $MiAi$ can be given by $MAi=\rho A\Delta Z\xd7diag(1,1,1,rgX2,rgY2,rgZ2)$, where

*ρA*Δ

*Z*is the mass of the rigid body, and

*r*

_{gX},

*r*

_{gY}, and

*r*

_{gZ}denote the radius of gyration with respect to the axes of

*A*

_{i}−

*XYZ*. The radius of gyration can be expressed as follows:

*Z*→ 0, we can have the following second-order partial differential equation:

*Z*throughout this paper). It is noted here that the elements of $D^$ are functions of

*Z*. It is also noted that

Since we assumed that the beam is prismatic, i.e., *dA*/*dZ* = 0, the term *A* of Eqs. (31) and (32) can be canceled in Eq. (30). However, if *dA*/*dZ* ≠ 0, e.g., a tapered beam, then the area *A* must be considered.

As is well known, Eq. (30) is similar to the Sturm–Liouville problem (SLP). Considering that $M$ and $K$ are positive definite, this may be assumed to be a regular problem (although this is not proved in the present paper). One of the important properties of the regular SLP is that there is an infinite but countable number of *ω*^{2}, suggesting the presence of infinitely many natural frequencies and corresponding mode shapes.

### 3.3 Boundary Conditions.

*Z*=

*L*because the beam does not experience any force and moment at its end. Then, it is clear that $D^\u2032(L)=0$ since $K$ is positive definite. Lastly, we consider a beam supported by revolute joints, as shown in Fig. 7. The rigid element attached to the joint can be freely rotated about the joint axis, but it is clamped in all other directions. Also, there will be no torque reaction about the axis. For the revolute joints that are parallel or aligned with the

*X*-axis, the boundary conditions can be written as

Z | $D^$ | $D^\u2032$ | |
---|---|---|---|

Fixed | foot (Z = 0) | $D^=0$ | — |

end (Z = L) | $D^=0$ | — | |

Free | foot (Z = 0) | — | $D^\u2032=0$ |

end (Z = L) | — | $D^\u2032=0$ | |

Revolute joint (X-axis) | foot (Z = 0) | $D^={000;100}T$ | $\varphi X\u2032=0$ |

end (Z = L) | $D^={0L0;100}T$ | $\varphi X\u2032=0$ | |

Revolute joint (Y-axis) | foot (Z = 0) | $D^(Z=0)={000;010}T$ | $\varphi Y\u2032=0$ |

end (Z = L) | $D^(Z=L)={\u2212L00;010}T$ | $\varphi Y\u2032=0$ |

Z | $D^$ | $D^\u2032$ | |
---|---|---|---|

Fixed | foot (Z = 0) | $D^=0$ | — |

end (Z = L) | $D^=0$ | — | |

Free | foot (Z = 0) | — | $D^\u2032=0$ |

end (Z = L) | — | $D^\u2032=0$ | |

Revolute joint (X-axis) | foot (Z = 0) | $D^={000;100}T$ | $\varphi X\u2032=0$ |

end (Z = L) | $D^={0L0;100}T$ | $\varphi X\u2032=0$ | |

Revolute joint (Y-axis) | foot (Z = 0) | $D^(Z=0)={000;010}T$ | $\varphi Y\u2032=0$ |

end (Z = L) | $D^(Z=L)={\u2212L00;010}T$ | $\varphi Y\u2032=0$ |

*Y*-axis can be given by

Besides the aforementioned supports, we can also obtain the boundary conditions for any support using the reciprocal relation between wrenches and instantaneous velocity [2].

### 3.4 Analytical Solutions.

In this subsection, Eq. (30) is solved for the cantilever. We again assume a uniform and prismatic beam with length *L*, area of cross section *A*, Young's modulus *E*, shear modulus *G*, and density *ρ*. The global coordinate frame *O* − *XYZ* is set as illustrated in Fig. 3.

From these equations, we can see that they are decoupled into axial vibrations, rotational vibrations with respect to the *Z*-axis, and bending vibration about the principal directions of the cross section. Consequently, this indicates that the natural frequencies (*ω*′s) are determined from each decoupled equation, and the corresponding mode shapes are axial and rotational vibrations about the *Z*-axis and bending in the principal directions of the cross section.

It should be noted here that *ω*_{p} can be determined by the boundary conditions: *δ*_{X}(0) = 0, $\delta X\u2032(L)=0$, *ϕ*_{Z}(0) = 0, and $\varphi Z\u2032(L)=0$.

*λ*

_{i}and $ui$ for

*i*= 1, …, 4 are the eigenvalues and eigenvectors of $P$, respectively. Then, the solution of $u$ can be expressed as

*a*

_{i}are the coefficients, and

*u*

_{ji}are the

*j*th elements of $ui$. Application of the boundary conditions (

*δ*

_{X}(0) = 0, $\delta X\u2032(L)=0$,

*ϕ*

_{Y}(0) = 0, $\delta Y\u2032(L)=0$) to Eq. (45) yields

*δ*

_{X}+

*Zϕ*

_{Y}) at

*Z*=

*L*is

For a nontrivial solution, the following characteristic equation must be zero: $detU=0$. Recalling that *λ*_{i} and $ui$ are the eigenvalues and the corresponding eigenvectors of $P$, respectively, it can be said that $detU$ depends on *ω*. The determinant of $U$ is not always zero, but there may exist *ω* so that $detU=0$, suggesting that *ω* is the natural frequency for $detU=0$. As previously mentioned, when considering the regular SLP, there are infinitely many values of *ω* that satisfy $detU=0$. If the above is iterated, then we can find the natural frequencies *ω*_{p} for *p* = 1, 2, …, *n*, …, ∞. Once *ω*_{p} is found, the mode shapes depending on *a*_{i}, *λ*_{i}, and $vi$ can be determined.

*μ*

_{i}and

*v*_{i}for

*i*= 1, …, 4 are the eigenvalues and eigenvectors of $Q$, respectively, then $v$ can be written as

*b*

_{i}are coefficients. The boundary conditions (

*δ*

_{Y}(0) = 0, $\delta Y\u2032(L)=0$,

*ϕ*

_{X}(0) = 0, $\delta X\u2032(L)=0$) provide the characteristic equation: $detV=0$, where

*v*

_{ji}is the

*j*th element of $vi$. Again, we need to specify

*ω*and compute $detV$. For $detV=0$, the natural frequencies and the corresponding mode shapes can be determined.

There may be a degenerate case that both $detU$ and $detV$ become zero for the same *ω*. Since this is a linear vibration system, the solutions can be expressed by linear combinations of vibrational modes with the same natural frequency. In other words, the solutions can be a mixture of some modes with the same *ω*. The degeneracy sometimes occurs in symmetrical beams such as a cylindrical beam.

Until now, we have described how to solve the differential equations. As mentioned above, $U$ and $V$ are dependent on the boundary conditions, as opposed to $P$ and $Q$ which are not. This means that the natural frequencies and the corresponding mode shapes are affected by the boundary conditions. Therefore, it is important to construct $U$ and $V$ using the boundary conditions and find *ω* through the iterative method presented. However, it should be mentioned that the natural frequencies for bending vibrations cannot be explicitly expressed.

## 4 Discrete Analysis of Beam

*n*rigid bodies, and then the EOM is given by Eq. (5). As already mentioned, in order for a multibody system to have the same dynamic characters of the beam, the mass properties and elasticity are reflected into the mass and stiffness matrices $M$ and $K$. The mass matrix of the

*i*th rigid body $Mi$ can be easily obtained by the coordinate transformation in Eq. (26):

*Z*=

*L*/

*n*. However, we should consider not only elasticity but also boundary conditions when determining the stiffness matrix. We assume that the boundary conditions are given by constraints on both ends of the beam. Then, they affect obviously only $K1$ and $Kn$. So, we can use Eqs. (20) and (21) to determine the stiffness matrices except for $K1$ and $Kn$ as follows:

For $K1$ and $Kn$, it is necessary to modify the EOM according to the type of supports. The following subsections explain how to apply boundary conditions for each case.

### 4.1 Fixed and Free Ends.

For a beam fixed at end, the rows and columns including $Mn$ and $Kn$ are vanished. When a beam is fixed at both ends, we should erase the rows and columns containing all of $M1$, $K1$, $Mn$ and $Kn$.

When the foot of a beam is free, it is clear that $K1=K1,2$ since the first rigid body is connected with the second one. For free end, $Kn=Kn\u22121,n$.

Now that we have determined both $M$ and $K$ for fixed and free ends, eigenanalysis can be performed to obtain natural frequencies and vibration modes after assuming small harmonic displacements.

### 4.2 General Boundary Conditions.

In this subsection, we aim to discuss supports with *p*-DOF. While the method in the previous subsection can cover 0- and 6-DOF supports, they are insufficient for describing supports with *n*-DOF such as revolute, cylindrical, and spherical joints. Therefore, we introduce a more general approach that can cover arbitrary supports in this subsection.

*p*- and

*q*-DOF, respectively. Then, the displacement vectors of the first and last rigid bodies can be expressed as

We have derived the EOM governing vibrations of a multibody system equivalent to a beam for each type of support. The EOM is determined by the DOF of the rigid bodies at both ends, and the natural frequencies and vibration modes can be obtained through eigenanalysis. Furthermore, this method can be applied to a wider range of cases, considering its ability to handle boundary conditions of other bodies. This suggests that it can be also effectively utilized for the discrete analysis of plates.

## 5 Numerical Examples

This section introduces two numerical examples, namely the cantilever and the robot link. First, for the cantilever, we explain the approach employed in this paper and compare the accuracy of this method with numerical solutions (or discrete system analysis). Second, we provide an analytic solution for the robot link (similar to the example described in Ref. [2]) where the revolute joint axes are not parallel and show the discretization into an equivalent multibody system.

### 5.1 Cantilever Beam.

In this example, we consider a uniform and prismatic cantilever with a rectangular cross section (Fig. 8). The parameters of the beam are given in Table 2. Using these parameters, Eqs. (40) and (41) can provide the natural frequencies for axial and rotational vibrations about the *Z*-axis, as listed in Table 3.

Material | $\rho =8050(kg/m3)$, $E=200\xd7109(N/m2)$, $G=7.6923\xd71010(N/m2)$ |

Geometry | h = 0.2 (m) (height), w = 0.1 (m) (width), L = 1.5 (m), R_{X} = 0.0289 (m^{2}), R_{Y} = 0.0577 (m^{2}) |

Material | $\rho =8050(kg/m3)$, $E=200\xd7109(N/m2)$, $G=7.6923\xd71010(N/m2)$ |

Geometry | h = 0.2 (m) (height), w = 0.1 (m) (width), L = 1.5 (m), R_{X} = 0.0289 (m^{2}), R_{Y} = 0.0577 (m^{2}) |

Types | ω_{1} | ω_{2} | ω_{3} | ω_{4} | ω_{5} |
---|---|---|---|---|---|

δ_{Z} | 5219.70 | 15,659.10 | 26,098.51 | 36,537.91 | 46,977.31 |

ϕ_{Z} | 3237.12 | 9711.36 | 16,185.61 | 22,659.85 | 29,134.09 |

Types | ω_{1} | ω_{2} | ω_{3} | ω_{4} | ω_{5} |
---|---|---|---|---|---|

δ_{Z} | 5219.70 | 15,659.10 | 26,098.51 | 36,537.91 | 46,977.31 |

ϕ_{Z} | 3237.12 | 9711.36 | 16,185.61 | 22,659.85 | 29,134.09 |

We then compute the eigenvalues and the corresponding eigenvectors of $P$ and $Q$, respectively, iteratively for 0 < *ω* < 10,000[rad/s] with a step size of 0.01. For $detU=0$ or $detV=0$, the solutions are listed in Table 4. It is noted that $U$ and $V$ are described by Eqs. (48) and (52), respectively. Now, we can line up the natural frequencies and the corresponding mode shapes for all vibrations, as illustrated in Fig. 9. As shown, our results are similar to those of the classical beam theory.

ω_{1} | ω_{2} | ω_{3} | |
---|---|---|---|

$\delta X$, $\varphi Y(P)$ | 444.25 | 2603.72 | 6676.63 |

$\delta Y$, $\varphi X(Q)$ | 224.16 | 1379.68 | 3758.78 |

ω_{1} | ω_{2} | ω_{3} | |
---|---|---|---|

$\delta X$, $\varphi Y(P)$ | 444.25 | 2603.72 | 6676.63 |

$\delta Y$, $\varphi X(Q)$ | 224.16 | 1379.68 | 3758.78 |

$Mi$, $Ki$, and $Ki\u22121,i$ were calculated by using Eqs. (53)–(55). To perform eigenanalysis, we again assume small vibrations, i.e., $D=D~ej\omega t$, where $D~={D^2TD^3T\cdots D^nT}T$. $D^i$ is the vibrational axis of the *i*th rigid body, and the eigenequation is $\u2212\omega 2MD~+KD~=0$, which was solved in matlab for *n* = 5, 10, 20, and 30.

Our results for the natural frequencies are listed in Table 5. It can be seen that the numerical solutions approach exact solutions as *n* increases. For bending vibration (except *ω*_{5}), the error rates between the exact and numerical results decrease from approximately 27% (*n* = 5) to 3% (*n* = 30). In addition, the error rate for the axial vibration (*ω*_{5}) decreases from 11% to 2%. The error rates are calculated by $(E\u2212N)/E\xd7100(%)$, where *E* and *N* denote the exact and numerical solutions, respectively. These results suggest that the error rates depend on the type of vibration and the number of rigid bodies *n*. In contrast, the error rates are independent of high and low frequencies. It is generally known that a larger error occurs at a higher frequency. In this sense, a beam modeled as a stack of rigid bodies and the concept of virtual springs may be utilized for the analysis of beams with more complex geometry.

n | ω_{1} | ω_{2} | ω_{3} | ω_{4} | ω_{5} | ω_{6} | |
---|---|---|---|---|---|---|---|

Exact solutions | — | 224.16 | 444.25 | 1379.68 | 2603.72 | 3237.12 | 3758.78 |

Numerical solutions (error rate) | 5 | 279.68 (24.8%) | 552.98 (24.4%) | 1759.19 (27.5%) | 3262.97 (25.3%) | 3578.57 (10.5%) | 4802.90 (27.8%) |

10 | 248.92 (11.0%) | 492.83 (10.9%) | 1538.20 (11.5%) | 2883.30 (10.7%) | 3403.62 (5.1%) | 4195.63 (11.6%) | |

20 | 235.90 (5.2%) | 467.31 (5.2%) | 1452.68 (5.3%) | 2733.36 (5.0%) | 3319.22 (2.5%) | 3956.33 (5.3%) | |

30 | 231.86 (3.4%) | 459.37 (3.4%) | 1427.09 (3.4%) | 2688.08 (3.2%) | 3291.60 (1.7%) | 3886.24 (3.4%) |

n | ω_{1} | ω_{2} | ω_{3} | ω_{4} | ω_{5} | ω_{6} | |
---|---|---|---|---|---|---|---|

Exact solutions | — | 224.16 | 444.25 | 1379.68 | 2603.72 | 3237.12 | 3758.78 |

Numerical solutions (error rate) | 5 | 279.68 (24.8%) | 552.98 (24.4%) | 1759.19 (27.5%) | 3262.97 (25.3%) | 3578.57 (10.5%) | 4802.90 (27.8%) |

10 | 248.92 (11.0%) | 492.83 (10.9%) | 1538.20 (11.5%) | 2883.30 (10.7%) | 3403.62 (5.1%) | 4195.63 (11.6%) | |

20 | 235.90 (5.2%) | 467.31 (5.2%) | 1452.68 (5.3%) | 2733.36 (5.0%) | 3319.22 (2.5%) | 3956.33 (5.3%) | |

30 | 231.86 (3.4%) | 459.37 (3.4%) | 1427.09 (3.4%) | 2688.08 (3.2%) | 3291.60 (1.7%) | 3886.24 (3.4%) |

As listed in Table 6, the vibration axes of the first, third, and sixth modes lie on the *ZX* plane, whereas those of the second and fourth modes lie on the *YZ* plane, and the vibration axes of the fifth modes are aligned with the *Z*-axis. All mode shapes are shown in Fig. 10, and our results are almost the same as the ones obtained by the screw approach. Referring to Fig. 11, the mode shapes can be thought of as rotations of rigid bodies about the respective vibration axes, and this thought is the most important concept in the screw theory of a linear vibration, such as the vibration center and axis [4]. Therefore, our results suggest that this approach can be extended to describe the vibration of two connected rigid bodies [13] but also to a *n*-body vibration system.

Form of vector | Results of vibration axes | |
---|---|---|

1st mode | $D^i=\varphi X{0Z0;100}T$(parallel to the X-axis) | $D^5=0.07{00.120;100}T$ |

$D^10=0.14{00.230;100}T$ | ||

$D^20=0.21{00.380;100}T$ | ||

$D^30=0.21{00.430;100}T$ | ||

2nd mode | $D^i=\varphi Y{\u2212Z00;010}T$(parallel to the Y-axis) | $D^5=0.07{\u22120.1100;010}T$ |

$D^10=0.14{\u22120.2200;010}T$ | ||

$D^20=0.21{\u22120.3700;010}T$ | ||

$D^30=0.22{\u22120.4100;010}T$ | ||

3rd mode | $D^i=\varphi X{0Z0;100}T$(parallel to the X-axis) | $D^5=\u22120.09{00.100;100}T$ |

$D^10=\u22120.08{00.060;100}T$ | ||

$D^20=0.12{01.250;100}T$ | ||

$D^30=0.21{01.190;100}T$ | ||

4th mode | $D^i=\varphi Y{\u2212Z00;010}T$(parallel to the Y-axis) | $D^5=\u22120.08{\u22120.0700;010}T$ |

$D^10=\u22120.08{0.0200;010}T$ | ||

$D^20=0.12{\u22121.2500;010}T$ | ||

$D^30=0.21{\u22121.1800;010}T$ | ||

5th mode | $D^i=\varphi Z{000;001}T$(aligned with the Z-axis) | $D^5=0.06{000;001}T$ |

$D^10=0.12{000;001}T$ | ||

$D^20=0.22{000;001}T$ | ||

$D^30=0.26{000;001}T$ | ||

6th mode | $D^i=\varphi X{0Z0;100}T$(parallel to the X-axis) | $D^5=0.11{00.070;100}T$ |

$D^10=\u22120.02{02.090;100}T$ | ||

$D^20=\u22120.06{00.520;100}T$ | ||

$D^30=0.22{01.310;100}T$ |

Form of vector | Results of vibration axes | |
---|---|---|

1st mode | $D^i=\varphi X{0Z0;100}T$(parallel to the X-axis) | $D^5=0.07{00.120;100}T$ |

$D^10=0.14{00.230;100}T$ | ||

$D^20=0.21{00.380;100}T$ | ||

$D^30=0.21{00.430;100}T$ | ||

2nd mode | $D^i=\varphi Y{\u2212Z00;010}T$(parallel to the Y-axis) | $D^5=0.07{\u22120.1100;010}T$ |

$D^10=0.14{\u22120.2200;010}T$ | ||

$D^20=0.21{\u22120.3700;010}T$ | ||

$D^30=0.22{\u22120.4100;010}T$ | ||

3rd mode | $D^i=\varphi X{0Z0;100}T$(parallel to the X-axis) | $D^5=\u22120.09{00.100;100}T$ |

$D^10=\u22120.08{00.060;100}T$ | ||

$D^20=0.12{01.250;100}T$ | ||

$D^30=0.21{01.190;100}T$ | ||

4th mode | $D^i=\varphi Y{\u2212Z00;010}T$(parallel to the Y-axis) | $D^5=\u22120.08{\u22120.0700;010}T$ |

$D^10=\u22120.08{0.0200;010}T$ | ||

$D^20=0.12{\u22121.2500;010}T$ | ||

$D^30=0.21{\u22121.1800;010}T$ | ||

5th mode | $D^i=\varphi Z{000;001}T$(aligned with the Z-axis) | $D^5=0.06{000;001}T$ |

$D^10=0.12{000;001}T$ | ||

$D^20=0.22{000;001}T$ | ||

$D^30=0.26{000;001}T$ | ||

6th mode | $D^i=\varphi X{0Z0;100}T$(parallel to the X-axis) | $D^5=0.11{00.070;100}T$ |

$D^10=\u22120.02{02.090;100}T$ | ||

$D^20=\u22120.06{00.520;100}T$ | ||

$D^30=0.22{01.310;100}T$ |

### 5.2 A Robot Link.

In this subsection, a robot link is used as an example similar to the one proposed in Ref. [2] to demonstrate the utility of the screw theory. In general, the classical beam theory does not deal with hinges that are not parallel at both ends, such as the robot link. However, the boundary conditions only affect $U$ and $V$, as shown in the methodology presented in this paper, and we can therefore easily obtain analytical solutions for the robot link.

*X*-axis, and the other hinge is parallel to the

*Y*-axis. The system parameters are given in Table 7. Referring to Table 1, the boundary conditions are given by

Material | $\rho =8050(kg/m3)$, $E=200\xd7109(N/m2)$, $G=7.6923\xd71010(N/m2)$ |

Geometry | $h=w=0.1(m)$, $L=1.5(m)$, $RX=0.0289(m2)$, $RY=0.0289(m2)$ |

Material | $\rho =8050(kg/m3)$, $E=200\xd7109(N/m2)$, $G=7.6923\xd71010(N/m2)$ |

Geometry | $h=w=0.1(m)$, $L=1.5(m)$, $RX=0.0289(m2)$, $RY=0.0289(m2)$ |

*ω*

_{a}and

*ω*

_{r}as follows:

To construct $U$ and $V$, the boundary conditions are arranged as

*λ*

_{i}and eigenvectors $ui$ of $P$, the boundary conditions of Eq. (68) can be rewritten as

By iteration, we can find the six degenerate solutions for 0 < *ω* < 20 000[rad/s] as listed in Table 8. That is, $detU=detV=0$ for the solutions. Considering the symmetry of the robot link, we can say that the degeneracy is natural. Therefore, two vibration modes occur with the same natural frequency. It means that mode shape can be expressed by linear combinations of the two vibration modes, and therefore, there are infinitely many mode shapes for a natural frequency. Figure 12 illustrates mode shapes where |*ϕ*_{X}(*Z* = 0)| = |*ϕ*_{Y}(*Z* = *L*)| among numerous vibration modes for the first and second natural frequencies, respectively. The thick lines depict the axes of revolute joints. This figure indicates that the robot link does not have pure bending modes. In terms of finding the analytic solutions for the robot link, the proposed screw approach can be useful for the vibration analysis of a beam that experiences an arbitrary bending.

ω_{1} | ω_{2} | ω_{3} | ω_{4} | ω_{5} | ω_{6} |
---|---|---|---|---|---|

971.34 | 3069.50 | 6186.08 | 10,146.86 | 14,781.84 | 19,943.06 |

ω_{1} | ω_{2} | ω_{3} | ω_{4} | ω_{5} | ω_{6} |
---|---|---|---|---|---|

971.34 | 3069.50 | 6186.08 | 10,146.86 | 14,781.84 | 19,943.06 |

*n*. The displacement vector is defined as $D={aD2T\cdots Dn\u22121Tb}T$. For small harmonic displacements, the vector can be written as $D=D~ej\omega t$ where $D~={a^D^2T\cdots D^n\u22121Tb^}T$. Once again, the eigenanalysis was conducted for

*n*= 5, 10, 20, and 30.

Table 9 shows the results of the natural frequencies *ω*_{i} with degeneracy as well as the nondegenerate ones. The vibration modes for *n* = 30 are listed in Table 10. Observation of Tables 9 and 10 reveals that the degenerate vibration modes are two linearly independent vectors and that the nondegenerate ones are axial or rotational vibrations when the robot link is fixed at both ends, i.e., *a* = *b* = 0. In addition, the error rates decrease as *n* increases, as shown in Table 10. Consequently, we can say that the numerical solutions converge to exact solutions.

n | ω_{1,2} | ω_{3,4} | ω_{5,6} | ω_{7} | ω_{8} | ω_{9} |
---|---|---|---|---|---|---|

5 | 1235.91 (27.1%) | 3964.45 (29.1%) | 7828.93 (26.6%) | 7886.40 | 12,716.44 | 14,572.17 |

n | $\omega 1,2$ | $\omega 3,4$ | $\omega 5,6$ | $\omega 7$ | $\omega 8,9$ | $\omega 10$ |

10 | 1082.99 (11.5%) | 3435.11 (11.9%) | 6932.27 (12.1%) | 7157.14 (10.5%) | 11,347.25 (11.8%) | 11,540.54 (10.5%) |

20 | 1023.08 (5.3%) | 3235.22 (5.4%) | 6520.89 (5.4%) | 6807.23 (5.14%) | 10,689.36 (5.3%) | 10,976.33 (5.14%) |

30 | 1005.01 (3.5%) | 3176.54 (3.5%) | 6401.46 (3.5%) | 6694.21 (3.4%) | 10,496.16 (3.4%) | 10,794.10 (3.4%) |

n | ω_{1,2} | ω_{3,4} | ω_{5,6} | ω_{7} | ω_{8} | ω_{9} |
---|---|---|---|---|---|---|

5 | 1235.91 (27.1%) | 3964.45 (29.1%) | 7828.93 (26.6%) | 7886.40 | 12,716.44 | 14,572.17 |

n | $\omega 1,2$ | $\omega 3,4$ | $\omega 5,6$ | $\omega 7$ | $\omega 8,9$ | $\omega 10$ |

10 | 1082.99 (11.5%) | 3435.11 (11.9%) | 6932.27 (12.1%) | 7157.14 (10.5%) | 11,347.25 (11.8%) | 11,540.54 (10.5%) |

20 | 1023.08 (5.3%) | 3235.22 (5.4%) | 6520.89 (5.4%) | 6807.23 (5.14%) | 10,689.36 (5.3%) | 10,976.33 (5.14%) |

30 | 1005.01 (3.5%) | 3176.54 (3.5%) | 6401.46 (3.5%) | 6694.21 (3.4%) | 10,496.16 (3.4%) | 10,794.10 (3.4%) |

$a^$ | $D^16$ | $b^$ | ||
---|---|---|---|---|

ω_{1,2} | $D~1$ | $0.23$ | ${0\u22120.130\u22120.05400}T$ | $0$ |

$D~2$ | $0$ | ${0.0240000.070}T$ | $0.20$ | |

ω_{3,4} | $D~3$ | $0.22$ | ${0\u22120.140\u22120.2100}T$ | $0$ |

$D~4$ | $0$ | ${0.16000\u22120.180}T$ | $0.21$ | |

ω_{5,6} | $D~5$ | $0.23$ | ${0\u22120.130\u22120.05400}T$ | $0$ |

$D~6$ | $0$ | ${\u22120.060000.120}T$ | $0.21$ | |

ω_{7} | $D~7$ | $0$ | ${00000\u22120.26}T$ | $0$ |

ω_{8,9} | $D~8$ | $0.23$ | ${0\u22120.130\u22120.05400}T$ | $0$ |

$D~9$ | $0$ | ${\u22120.130000.160}T$ | $0.21$ | |

ω_{10} | $D~10$ | $0$ | ${00\u22120.26000}T$ | $0$ |

$a^$ | $D^16$ | $b^$ | ||
---|---|---|---|---|

ω_{1,2} | $D~1$ | $0.23$ | ${0\u22120.130\u22120.05400}T$ | $0$ |

$D~2$ | $0$ | ${0.0240000.070}T$ | $0.20$ | |

ω_{3,4} | $D~3$ | $0.22$ | ${0\u22120.140\u22120.2100}T$ | $0$ |

$D~4$ | $0$ | ${0.16000\u22120.180}T$ | $0.21$ | |

ω_{5,6} | $D~5$ | $0.23$ | ${0\u22120.130\u22120.05400}T$ | $0$ |

$D~6$ | $0$ | ${\u22120.060000.120}T$ | $0.21$ | |

ω_{7} | $D~7$ | $0$ | ${00000\u22120.26}T$ | $0$ |

ω_{8,9} | $D~8$ | $0.23$ | ${0\u22120.130\u22120.05400}T$ | $0$ |

$D~9$ | $0$ | ${\u22120.130000.160}T$ | $0.21$ | |

ω_{10} | $D~10$ | $0$ | ${00\u22120.26000}T$ | $0$ |

Figure 13 illustrates the mode shapes of robot links. The vibration modes with degeneracy can be expressed by linear combinations of two vectors: $D~=\alpha D~i+\beta D~i+1$, and thus, there are infinitely many mode shapes with a natural frequency. The mode shapes illustrated in Fig. 13 were selected among the mode shapes so that |*ϕ*_{X}(*Z* = 0)| = |*ϕ*_{Y}(*Z* = *L*)|. In this case, the vibration axis of the *i*th piece $D^i$ is not a line but a screw. That is, the mode shape is a repetitive screw motion. It suggests that the mode shape can be considered as screw motions as well as rotations. This is a key insight provided by the screw theory.

So far, we have interpreted the dynamic characteristics of the robot link both analytically and numerically. Through both approaches, we were able to obtain nondegenerate and degenerate solutions and also confirmed that the numerical solutions converge to the analytical solutions. In the future work, we will apply the approaches to problems with comply geometry and plate theory.

## 6 Conclusions

The present study demonstrates that the screw theory could accurately describe the vibrations of discrete and continuous systems. We initially modeled a beam as a stack of rigid bodies connected by virtual springs to each other, and then derived the EOM governing different types of vibrations, such as extension, compression, torsion, and bending, and introduced a method to solve these equations analytically for each boundary condition.

Our study makes three major contributions. First, the stiffness matrix of the beam was derived by using virtual springs acting as normal and shear stresses. This representation could be very useful for deformation problems of continuous systems, such as in ones evaluating the dynamics and statics of the beam. In addition, this approach also enabled us to easily obtain the stiffness matrix of a multibody system by coordinate transformation. Second, it was shown that a continuous system could be approximated as a stack of rigid bodies connected by virtual springs to each other, and the effectiveness of this approach was confirmed by a comparison between analytical and numerical results. Our findings also revealed that the mode shapes could be thought of as rotations (or generally screw motions) of rigid bodies about the vibration axes. Lastly, we analytically solved the differential equations for nonstandard forms of the classical beam theory. We showed that analytical solutions could be easily obtained by the proposed iterative method, provided that the boundary conditions were clearly specified.

Throughout this paper, we only considered deformations in one direction (*Z*-axis). As a result, the stiffness matrix of the beam was derived, and the mode shapes could be understood as rotations of rigid bodies. However, these concepts could be extended into deformations in two directions, such as the plate theory. It is expected that the stiffness matrix can still be used, and the mode shapes of the plate could be considered as screw motions (or rotations) of small rigid plates. In the future, we intend to investigate the applicability of this approach by considering the screw theory.

## Acknowledgment

This research was supported by a grant from R&D Program (Development of core technology for evacuation control in the accident case and passenger safety, PK2302A2) of the Korea Railroad Research Institute and the Korea Agency for Infrastructure Technology Advancement (KAIA) grant funded by the Ministry of Land, Infrastructure and Transport (Grant RS-2023-00238018, Development of technology for cognizing, predicting and responding to high-risk disasters for deep railway).

## Funding Data

Korea Railroad Research Institute (KRRI) and Korea Agency for Infrastructure Technology Advancement (KAIA).

## Conflict of Interest

There are no conflicts of interest. This article does not include research in which human participants were involved. Informed consent was obtained for all individuals. Documentation provided upon request. This article does not include any research in which animal participants were involved.

## Data Availability Statement

The authors attest that all data for this study are included in the paper.

## Nomenclature

- $s$ =
unit direction vector of $s^i$

*A*=area of cross section

*E*=Young's modulus

*G*=shear modulus

*L*=length of a beam

- $D$ =
time-dependent displacement vector

- $K$ =
stiffness matrix

- $M$ =
mass matrix

- $s^i$ =
line vector expressed in Plücker's axis coordinates

- $D~$ =
screw chain

- $D^$ =
time-independent displacement vector (or vibration axes)

- $ri$ =
position vector of $s^i$

- $Ki$ =
stiffness matrix of all springs attached to

*i*th rigid body- $Ki,j$ =
stiffness matrix of springs connecting

*i*th and*j*th bodies- $Mi$ =
mass matrix of

*i*th rigid body*O*−*XYZ*=global frame

- $D^(Z)$ =
screw deflection

*A*_{i}−*XYZ*=local frame placed in

*i*th rigid body*B*_{i}−*XYZ*=local frame located in

*i*th cut plane- $\delta $ =
small translational vector

*ρ*=density

*ϕ*=small angular vector

*ω*_{p}=*p*th natural frequency