## 3 Finite Element Analysis Considering Support Stiffness and Axial Force  
### Stiffness and Mass Matrices of the Shaft Element  
FEA, which is widely used to solve structural vibration response [17; 18], was used in this study to calculate the natural frequencies of the turning workpiece. In order to solve the modal response of the finite element model, we used MATLAB to build finite elements. In this model, the effects of shear deformation and moment of inertia were ignored because the ratio of the diameter to the length was very small; thus, an Euler-Bernoulli beam element was used [23]. The finite elements are uniformly distributed, and the model converges rapidly through iterative solution. Therefore, only 6 units were set up in this model, each unit had 4 nodes and 4 degrees of freedom.  
A typical transversely vibrating element of the shaft is shown in Figure 3. The transverse displacement \(w(x,t)\) of each point in the element is assumed to be consistent with the static deformation curve of the beam along the \(x\) direction, obeying the cubic equation [31]  
\[w(x,t)=a(t)+b(t)x+c(t)x^{2}+d(t)x^{3}, \tag{4}\]  
where \(a(t),b(t),c(t)\text{and}d(t)\) are constants.  
The displacement of each node must satisfy the boundary conditions  
\[\begin{array}{l}w(0,t)=w_{1}(t),\frac{\partial w}{\partial x}(0,t)=\theta_{1 }(t)\ ;\\ w(l,t)=w_{2}(t),\frac{\partial w}{\partial  
From Equations (6) and Figure 3, \(w(x,t)\) can be expressed as  
\[\begin{array}{c}w(x,t)=\Big{(}1-3\frac{x^{2}}{l^{2}}+2\frac{x^{3}}{l^{3}}\Big{)} w_{1}(t)+\Big{(}\frac{x}{l}-2\frac{x^{2}}{l^{2}}+\frac{x^{3}}{l^{3}}\Big{)}l \theta_{1}(t)\\ +\Big{(}3\frac{x^{2}}{l^{2}}-2\frac{x^{3}}{l^{3}}\Big{)}w_{2}(t)+\Big{(}-\frac {x^{2}}{l^{2}}+\frac{x^{3}}{l^{3}}\Big{)}l\theta_{2}(t),\end{array} \tag{7}\]  
or more conveniently as  
\[w(x,t)=\sum_{i=1}^{4}N_{i}(x)\mathbf{w}(t), \tag{8}\]  
where \(\mathbf{w}(t)=\big{[}w_{1}(t)\ \ \theta_{1}(t)\ \ \ w_{2}(t)\ \ \theta_{2}(t)\big{]}^{T}\). Each element is determined by two nodes, each of which introduces two constraints (for a total of four). Therefore, the shape function \(N_{i}(x)\) is  
\[\begin{array}{c}N_{1}(x)=1-3\big{(}\frac{x}{l}\big{)}^{2}+2\big{(}\frac{x}{l }\big{)}^{3},\\ N_{2}(x)=x-2l\big{(}\frac{x}{l}\big{)}^{2}+l\big{(}\frac{x}{l}\big{)}^{3},\\ N_{3}(x)=3\big{(}\frac{x}{l}\big{)}^{2}-2\big{(}\frac{x}{l}\big{)}^{3},\\ N_{4}(x)=-l\big{(}\frac{x}{l}\big{)}^{2}+l\big{(}\frac{x}{l}\big{)}^{3}.\end{array} \tag{9}\]  
The kinetic and deformation energies of each element can be expressed as  
\[\begin{array}{c}T(t)=\frac{1}{2}\int_{0}^{l}\rho A\Big{[}\frac{\partial w(x, t)}{\partial t}\Big{]}^{2}dx=\frac{1}{2}\dot{\mathbf{w}}(t)^{T}\mathbf{M}\dot{\mathbf{w}}(t), \\ V(t)=\frac{1}{2}\int_{0}^{l}EI\Big{[}\frac{\partial^{2}w(x,t)}{\partial x^{2} }\Big{]}^{2}dx=\frac{1}{2}\mathbf{w}(t)^{T}\mathbf{K}\mathbf{w}(t),\end{array} \tag{10}\]  
where \(\rho\) is the mass density, \(E\) is the Young's modulus of the material, \(I\) is the moment of inertia of the cross section, \(A\) is the cross-sectional area, and \(\dot{\mathbf{w}}(t)=\big{[}dw_{1}(t)/dt\ \ d\theta_{1}(t)/dt\ \ dw_{2}(t)/dt\ \ d\theta_{2}(t)/dt\big{]}^{T}\).  
After substituting Equation (7) into Equation (10), the mass and stiffness matrices of the transverse vibration of the shaft are obtained by the Rayleigh-Lietz energy method [32]  
\[\mathbf{M}\mathbf{s}=\frac{\rho Al}{420}\begin{bmatrix}156&22l&54&-13l\\ 22l&4l^{2}&13l&-3l^{2}\\ 54&13l&156&-22l\\ -13l&-3l^{2}&-22l&4l^{2}\end{bmatrix}, \tag{11}\]  
\[\mathbf{K}\mathbf{s}=\frac{EI}{l^{3}}\begin{bmatrix}12&6l&-12&6l\\ 6l&4l^{2}&-6l&2l^{2}\\ -12&-6l&12&-6l\\ 6l&2l^{2}&-6l&4l^{2}\end{bmatrix}. \tag{12}\]  
Figure 3: Transverse vibration of shaft element.