# Free Vibration Analysis of Curvilinearly Tapered Axially Functionally Graded Material Beams
## 3 Solution Methodology
### Finite Element Modeling  
To ensure the accuracy and validation of the proposed semi-analytical methodology for analyzing the free vibration parameters of curvilinearly tapered AFGM beams, a finite element model has been developed and utilized for comparison purposes. It is well-known that the optimal modeling strategy involves using graded finite elements (GFEs), which allow for the variation of material properties at the element level. Additionally, advanced three-dimensional models are essential to precisely capture the non-uniform geometry of axially inhomogeneous beams [31].  
The free vibration analysis of AFGM beams is conducted using ABAQUS [47], incorporating a user-defined material subroutine (UMAT) to implement material gradients into a convenient three-dimensional element, C3D20, as detailed in [31; 33] and illustrated in Figure 2b. The material gradient is assigned within the element volume during the computation of the element stiffness matrix, \(\mathbf{K}^{c}\) using Gauss quadrature, as follows:  
\[\mathbf{K}^{c}=\sum\limits_{p=1}^{n_{p}}\sum\limits_{q=1}^{n_{q}}\sum\limits_{r= 1}^{n_{r}}\Big{(}\mathbf{B}^{T}\mathbf{D}(x)\mathbf{B}|\mathbf{J}|\Big{)}_{p,q,r}w_{p}w_{q}w_{r}, \tag{18}\]  
where the indices \(p\), \(q\), and \(r\) refer to the integration points over the element volume; \(|\mathbf{J}|\) is the determinant of the geometric Jacobian matrix; \(w_{p}\), \(w_{q}\), and \(w_{r}\) are the Gaussian weights; \(\mathbf{B}\) stands for a matrix of gradients of the shape functions for the 3-D element; and the material matrix \(\mathbf{D}(x)\) of the 3-D GFE relevant to this study can be found in [33].  
ABAQUS currently does not support calculating the mass matrix based on a user-defined mass density in frequency analysis [47]. The elements available in ABAQUS use a constant mass density, such as that noted in the following equation:  
\[\mathbf{M}^{e}=\int\limits_{V^{e}}\rho_{0}\mathbf{N}^{T}\cdot\mathbf{N}dV^{e}, \tag{19}\]  
where \(\mathbf{N}\) is an element interpolation matrix.  
To address this limitation, we assigned a constant effective mass density within the GFE by using an average value of the actual mass density distribution over the beam's volume in the following form:  
\[\rho_{avr}=\frac{1}{V}\int\limits_{V}\rho(x)dx \tag{20}\]  
In the case of the free vibration analysis without damping, the global discrete system of governing equations of motion is reduced to the problem of eigenvalues and associated eigenvectors, which is written as follows:  
\[\Big{(}\mathbf{K}-\omega^{2}\mathbf{M}\Big{)}\mathbf{X}=0, \tag{21}\]  
where \(\mathbf{M}\) and \(\mathbf{K}\) are global mass and stiffness matrices, respectively, obtained by assembling the formulated element mass and stiffness matrices; \(\omega\) is the undamped circular frequency; and \(\mathbf{X}\) is the vector of nodal displacements associated with oscillations at each \(\omega\). The free vibration analysis is carried out using linear perturbation load step, obtaining the eigenvalues by means of Lanczos or subspace iteration methods [47]. Figure 2a illustrates a three-dimensional finite element beam model created using the ABAQUS pre-processor, utilizing the 'loft' feature. To provide accurate results, the element aspect ratio in the mesh was consistently maintained close to one throughout all subsequent calculations.