Vibration Analysis of Randomly Oriented Carbon Nanotube Based on FGM Beam Using Timoshenko Theory

Abstract

The carbon nanotube (CNT) reinforced functionally graded materials (FGM) are expected to be the new generation materials having wide range of unexplored potential applications in various technological areas such as aerospace and structural and chemical industry. The present work deals with the finite element modeling and free vibration analysis of CNT based functionally graded beam using three-dimensional Timoshenko beam theory. It has been assumed that the material properties of CNT based FG beam vary only along the thickness and these properties are evaluated by rule of mixture. The extended Hamilton principle has been applied to find out the governing equations of CNT based FG beam. Finite element method is used to solve governing equation with the exact shape functions. Initial analysis deals with CNTs assumed to be oriented along the length direction only. But practically it is not possible. So, further work deals with the free vibration analysis of functionally graded nanocomposite beams reinforced by randomly oriented straight single walled carbon nanotubes (SWCNTs). The Eshelby-Mori-Tanaka approach based on an equivalent fiber is used to investigate the material properties of the beam. Results are presented in tabular and graphical forms to show the effects of carbon nanotube orientations, slenderness ratios, and boundary conditions on the dynamic behavior of the beam.

1. Introduction

Functionally graded material (FGM) belongs to a class of advanced material characterized by variation in properties as the dimension varies. FGM was invented with the prime requirement of thermal stability at high temperatures along with better mechanical properties. FGM concept was originated in Japan in 1984 during the space plane project.

Metallic part of FGM will take care of better mechanical properties and ceramic part will take care of thermal stability at high temperatures. The overall properties of FGM are unique and different from any of the individual materials that form it. Due to these advantages FGM has a wide range of applications such as in

• aerospace to make space plane bodies, components of the rocket engines,

• medicine to replace the living tissues like bones and teeth,

• defense as a penetration resistant material for bullet proof vests as well as armor plates,

• energy as thermal barrier, protective blade in gas turbine,

• optoelectronics as graded refractive index material,

• other applications consist of automobile engine component and nuclear reactor components, and so forth.

These applications explicitly require good thermal properties. So, it is important that analysis of thermal behavior of FGM should be carried out which has initiated the work related to this area. Lot of researchers have presented their work so far.

CNTs are known to be discovered by Iijima in 1991. CNT exhibits extraordinary mechanical properties: the Young's modulus is over 1 Tpa. The estimated tensile strength is 200 GPa. It is due to carbon-carbon sp² bonding. These properties are ideal for reinforced composites, nanoelectromechanical systems (NEMS). Carbon nanotubes (CNTs) are also known for their good mechanical and electrical properties. Thus the use of CNT with functionally graded material (FGM) provides definitely improved mechanical, electrical, and thermal properties. Important advantages of CNT based on FGM are listed below:

1. higher strength to density ratio,

2. higher stiffness to density ratio,

3. better fatigue and wear resistance,

4. better elevated temperature properties (higher strength-lower creep rate),

5. ability to fabricate directional mechanical properties,

6. provide multifunctionality,

7. provide ability to control the deformation, dynamic response of the system, wear and corrosion of parts, and so forth,

8. provide ability to remove stress concentrations,

9. provide opportunities to take the benefits of different material systems.

The key of using CNT based on FGM is that one can obtain these properties as per the requirement just by varying the distribution and composition of CNT. That is how one can get directional properties and can control other parameters. Another advantage stated above is the stress concentration free material. It is because the cross section shows that there are no layers inside the material and instead there is a continuous gradation of materials from top to bottom. So, there is no stress concentration and delamination of layers. As it can be manufactured from CNT, polymers, and ceramics, one can take the advantages of each constituent material.

CNT is a new form of carbon, configurationally equivalent to 2D graphene sheet rolled into a tubular form. It is grown by several techniques in the laboratory which is just a few nanometers in diameter and several microns in length. Main problems facing the manufacturers are difficulty in homogeneous dispersion of CNT into the matrix material which results in low interfacial bond strength between these constituents and it affects the quality of the product. There is also a problem of mass production with low cost. To overcome these problems up to certain extent, there are different techniques to manufacture the CNT based on FGM. Most common method is powder metallurgy. Besides this, there are also many methods such as melting and solidification, thermal spray including plasma and cold spraying, electrochemical deposition, vapor deposition, and nanoscale dispersion. The CNT based on functionally graded materials are used in wind turbines, tissue engineering, thin films of shape memory alloys, nanoelectromechanical systems such as micro sensors used in robotics, micro actuators, telecommunications, and transport industry.

All the work that has been carried out up till is presented here regarding the properties of FGM, vibration analysis of a nanocomposed functionally graded material shaft using finite element modeling. So, majorly, it deals with material sciences, rotor dynamics, FEM, and theory of vibrations. The history of this research consists of 2 parts as discussed above, that is, finding out material properties and vibration analysis.

In field properties of FGM metals are good conductors of heat and ceramics acts as good insulators. So, they are used individually according to the requirement, but there are certain areas where both are needed simultaneously. In these applications metals and ceramics are bonded directly. At elevated temperatures it is found that thermal stresses are produced which is responsible for debonding of metal and ceramics as well; the delimitation in ceramics takes place. To avoid this, scientist Kawasaki and Watanabe proposed graded interlayer concept for metals and ceramics. Such material possesses super heat resistant property and sufficient toughness which also reduces the thermal stresses at elevated temperatures. This concept is known as functionally graded material (FGM) [1]. So, one can say FGM is a material in which mass and weight proportion is varied along the cross section. Generally, ceramics have high percentage than metals on the outside so that it can be used at elevated temperatures. The concept of FGM is discovered after the concept of composite materials. In case of composite materials the metal fiber is surrounded by the matrix, which gives high specific strength and stiffness with minimum specific density. In materials science “functionally graded material (FGM)” may be characterized by the variation in composition and structure gradually over volume resulting in corresponding changes in the properties of the material. The materials can be designed according to specific function and applications. Various approaches based on the bulk (particulate processing) perform processing layer processing and melt processing are used to fabricate the functionally graded materials. The basic unit for FGM representation is maxel; the term maxel was introduced in 2005 by Rajeev Dwivedi and Radovan Kovacevic at Research Centre for Advanced Manufacturing (RCAM); the attributes of maxel include the location and volume fraction of individual material components.

An overview of FGM was taken by Mahamood et al. [2]. The fabrication processes, area of application, some recent research studies, and the most promising FGM fabrication method (solid freeform SFF) were presented by them. The processing technique of FGM consists of vapour deposition, powder metallurgy, centrifugal, and SFF. For future research, SFF method needs modification through extensive characterization so that the cost of FGM will be reduced and reliability will be increased. Durodola and Attia [3] focused mainly on FGM materials for rotating hollow and solid disks. In their paper, FGM was considered as a nonhomogeneous orthotropic material. The property variation is considered as exponential. The hoop stresses and radial stresses were calculated after finding out Young's modulus for FGM layer wise which vary along the radial direction only. Different values were obtained for different values of exponent. Aboudi et al. [4] considered the higher order micromechanical theory. It was more advantageous than micromechanics approach which was based on representative volume element. RVE do not deal with continuously changing properties due to nonuniform inclusion spacing. In their paper they considered only unidirectional grading of properties. It was applied for thermal and mechanical fields in FGM. Giunta et al. [5] considered the distribution of material properties according to power law in terms of the volume fraction of the material constituents. The properties such as Young's modulus, Poisson's ratio, and density vary along 1 direction and two perpendicular directions together or independently.

In field vibration analysis of CNT based on FGM the later inventions in the field of material science have proved that CNT with functionally graded material (FGM) provides improved mechanical, electrical, and thermal properties. Andrews et al. [6] fabricated the multiwalled carbon nanotube/polymer composites by shear mixing. Results showed that as concentration of nanotubes increases both strength and stiffness increase. Yas and Heshmati [7] worked on dynamics of functionally graded nanocomposite beams with randomly oriented single-walled carbon nanotubes. They found that under the action of moving load, CNT-FGM beam with symmetrical distribution gives superior properties than that of unsymmetrical distribution. Alshorbagy et al. [8] presented the dynamic characteristics of functionally graded beam with gradation of material axially or transversally through the thickness based on the power law. They converted geometrically nonuniform beam into axially or transversely uniform geometrical beam. The work carried out by Ke et al. [9] gave the information about the nonlinear free vibration of functionally graded nanocomposite beams with single-walled carbon nanotubes (SWCNTs). It was based on Timoshenko beam theory and von Karman geometric nonlinearity. They also considered the gradation of material properties in thickness direction only. Khosrozadeh and Hajabasi [10] worked on free vibration of embedded double-walled carbon nanotubes considering nonlinear interlayer van der Waals forces. It was based on Euler-Bernoulli theory. The variation of the interlayer distance along the circumference of DWCNTs causes the inner and outer tubes deflections which were modeled as the interlayer van der Waal force. Natsuki et al. [11] also worked on double-walled carbon nanotubes. Ray and Batra [12] demonstrated the active control of smart structures using single-walled carbon nanotube reinforced 1–3 piezoelectric composite. Micromechanical analysis was carried out to find out the values of effective piezoelectric moduli. Some more about the nonlinear dynamic response of nanotube-reinforced composite plates was carried out by Wang and shen [13]. The material properties were considered as function of temperature. The higher order shear deformation theory was considered along with von Karman type kinematic nonlinearity to find the equation of motion. Work carried out up till was not restricted to only vibration analysis; the vibration control was also included. For safe and efficient functioning of all rotating machines, reduction of rotor vibration is very important. Das et al. [14] worked on the vibration control and stability of the system using electromagnetic exciters. They provided suitable force of actuation which was obtained by varying the control current. The rotor shaft was modeled according to Rayleigh beam theory with viscous damping. To visualise the effect of damping, the unbalance response was plotted for controlled and uncontrolled situations. The active vibration control using fuzzy logic controllers was demonstrated by Hossain Nezhad Shiraji et al. [15]. They worked on both PID controller and FLC to dampen the plate vibration. For the fuzzy inputs and outputs the designed FLC used 5 membership functions. The other method of vibration analysis is modal analysis which was given by Kim and Lee [16]. The natural characteristics of the system such as frequency, damping, and mode shapes were used to describe a structure. The modal parameters may be determined by FEA or other such methods. Mostly modal analysis is used in dynamic problem, vibration, or acoustics.

Regarding the above-mentioned histories, it seems that there is a need for proposing a method for vibration analysis of randomly orientation CNT based on FG materials is introduced. Results are presented in tabular and graphical forms to show the effects of carbon nanotube orientations, slenderness ratios, and boundary conditions on the dynamic behavior of the beam. So the objectives of the present work are as stated below:

1. finding out the material properties of CNT based on FGM,

2. free vibration analysis of Timoshenko beam,

3. finding out the effect of CNT orientation.

2. Modelling of Beam

2.1. Calculation of Mechanical Properties of Beam

The modeling of beam starts with the calculation of material properties. As CNT volume is assumed to vary along the thickness only, material properties for each layer are calculated first and finally effective values for the entire beam are calculated. Here, linear variation of volume fraction of CNT (V_cnt) is considered. It is calculated by following formula:

$V_{\text{cnt}}=\left(\mathrm{1}-\frac{\mathrm{2}z}{h}\right)V_{\text{cnt}}^{*}.$

(1)

Here V_cnt^* depends on mass fraction and density of CNT and density of matrix and it is given by

$V_{\text{cnt}}^{*}={\mathrm{\Lambda }}_{\text{cnt}}+\left(\frac{{\rho }^{\text{cnt}}}{{\rho }^m}\right)-\left(\frac{{\rho }^{\text{cnt}}}{{\rho }^m}\right){\mathrm{\Lambda }}_{\text{cnt}},$

(2)

where Λ_cnt called as is mass fraction of CNT. For uniform CNT distribution V_cnt equals the V_cnt^*. After calculating volume fraction Young' modulus, shear modulus, Poisson's ratio, and density of a beam can be calculated as follows:

$\begin{array}{c}E={\eta }_j{V}_{\text{cnt}}{E}^{\text{cnt}}+V_m{E}^m,\\ \frac{{\eta }_{\mathrm{2}}}{G}=\frac{V_{\text{cnt}}}{G^{\text{cnt}}}+\frac{V_m}{G^m},\\ v=V_{\text{cnt}}{v}^{\text{cnt}}+V_m{v}^m,\\ \rho =V_{\text{cnt}}{\rho }^{\text{cnt}}+V_m{\rho }^m,\end{array}$

(3)

where E^cnt and G^cnt are Young's modulus in longitudinal, transverse direction and shear modulus for carbon nanotube. E ^m , G ^m are Young's modulus & shear modulus for matrix (ceramic). η _j (j = 1,2, 3) is CNT efficiency parameter.

2.2. Problem Formulation

Here first order shear deformation theory is considered for Timoshenko beam. The axial and transverse displacements for this beam are as follows,

$\begin{array}{c}u\left(x,y,z,t\right)=u_{\mathrm{0}}\left(x,t\right)-z\phi \left(x,t\right),\\ w\left(x,y,z,t\right)=w_{\mathrm{0}}\left(x,t\right).\end{array}$

(4)

In this case three degrees of freedom per node have been considered. They are longitudinal displacement, transverse displacement (along Z axis), and angular rotation about Y-axis. This is shown in Figure 1.

Figure 1:

Degree of freedom for single element.

Strain-displacement relation for small displacement can be derived from (4) as follows:

${\epsilon }_{xx}=\frac{{\partial u}_{\mathrm{0}}}{\partial x}-z\frac{\partial \phi }{\partial x},\hspace{1em}\hspace{1em}{\gamma }_{xz}=-\phi +\frac{{\partial w}_{\mathrm{0}}}{\partial x}.$

(5)

From above equation, the constitutive relations can be formed for longitudinal and shear stresses as follows:

$\left\{\begin{array}{c}{\sigma }_{xx}\\ {\tau }_{xz}\end{array}\right\}=\left[\begin{array}{cc}E\left(z\right)& \mathrm{0}\\ \mathrm{0}& G\left(z\right)\end{array}\right]\left\{\begin{array}{c}{\epsilon }_{xx}\\ {\gamma }_{xz}\end{array}\right\}.$

(6)

Constitutive relations are formed to calculate strain energy and kinetic energy of beam. For given case the final expressions for these energies are as follows:

$\begin{array}{c}{u}_b=\frac{\mathrm{1}}{\mathrm{2}}{{\displaystyle \int }}_v\left({\sigma }_{xx}{\epsilon }_{xx}+{\tau }_{xz}{\gamma }_{xz}\right)dA\hspace{0.17em}dx,\\ T_b=\frac{\mathrm{1}}{\mathrm{2}}\rho {{\displaystyle \int }}_v\left({\dot{u}}^{\mathrm{2}}+{\dot{w}}^{\mathrm{2}}\right)dA\hspace{0.17em}dx,\end{array}$

(7)

where the symbol $\dot{\left(\right)}$ shows the first differentiation of a property with respect to time. After this, Hamilton's principle is applied to get the differential equations in terms of two translational and one rotational degree of freedom, that is, u, w, and ϕ. According to Hamilton's principle

${{\displaystyle \int }}_{T_1}^{T_2}[\partial T_b-\partial u_b]=0.$

(8)

Using (7) in above equation, the governing equations can be

$\begin{array}{c}\delta u:I_{\mathrm{0}}\left(x\right){\ddot{u}}_{\mathrm{0}}-I_{\mathrm{1}}\left(x\right)\ddot{\phi }-\frac{\partial }{\partial x}\left(A_{\mathrm{11}}\left(x\right)\frac{{\partial u}_{\mathrm{0}}}{\partial x}\right)+\frac{\partial }{\partial x}\left(B_{\mathrm{11}}\left(x\right)\frac{\partial \phi }{\partial x}\right)=\mathrm{0},\\ \delta w:I_{\mathrm{0}}\left(x\right){\ddot{w}}_{\mathrm{0}}-\frac{\partial }{\partial x}\left(A_{\mathrm{55}}\left(x\right)(\frac{\partial w}{\partial x}-\phi )\right)=\mathrm{0},\\ \delta \phi :I_{\mathrm{2}}\left(x\right)\ddot{\phi }-I_{\mathrm{1}}\left(x\right){\ddot{u}}_{\mathrm{0}}+\frac{\partial }{\partial x}\left(B_{\mathrm{11}}\left(x\right)\frac{{\partial u}_{\mathrm{0}}}{\partial x}\right)-\frac{\partial }{\partial x}\left(D_{\mathrm{11}}\left(x\right)\frac{\partial \phi }{\partial x}\right)-A_{\mathrm{55}}\left(x\right)\left(\frac{\partial w}{\partial x}-\phi \right)=\mathrm{0}.\end{array}$

(9)

The values of stiffness coefficients and mass moments from (5) are given below. K _s is shear correction factor

$\begin{array}{c}\left[\begin{array}{ccc}{A}_{\mathrm{11}}\left(x\right)& B_{\mathrm{11}}\left(x\right)& D_{\mathrm{11}}\left(x\right)\end{array}\right]={{\displaystyle \int }}_{A}E\left(z\right)\left[\begin{array}{ccc}\mathrm{1}& z& z^{\mathrm{2}}\end{array}\right]dA,\\ A_{\mathrm{55}}\left(x\right)={{\displaystyle \int }}_A{K}_{s}G\left(z\right)dA,\\ \left[\begin{array}{ccc}{I}_{\mathrm{0}}\left(x\right)& I_{\mathrm{1}}\left(x\right)& I_{\mathrm{2}}\left(x\right)\end{array}\right]={{\displaystyle \int }}_A\rho \left(z\right)\left[\begin{array}{ccc}\mathrm{1}& z& z^{\mathrm{2}}\end{array}\right]dA.\end{array}$

(10)

2.3. Finite Element Analysis

To solve the above governing equation finite element analysis is implemented. Finite element analysis aim is to find out the field variables (displacement) at nodal points by approximate analysis. These field variables are related to the field variables of nodal points inside the element. Here employing the approximate solution method the governing equations are approximated by a system of ordinary differential equations. The solution of the above partial differential equations is assumed in the following form:

$\begin{array}{c}{u}_{\mathrm{0}}=c_{\mathrm{1}}+c_{\mathrm{2}}x+c_{\mathrm{3}}{x}^{\mathrm{2}},\\ w_{\mathrm{0}}=c_{\mathrm{4}}+c_{\mathrm{5}}x+c_{\mathrm{6}}{x}^{\mathrm{2}}+c_{\mathrm{7}}{x}^{\mathrm{3}},\\ \phi =c_{\mathrm{8}}+c_{\mathrm{9}}x+c_{\mathrm{10}}{x}^{\mathrm{2}}.\end{array}$

(11)

In the above equations, the order of interpolation of slope is less than that of order of transverse displacement. This is to avoid the shear locking of the element. Some constants from above equation can be defined as

$\begin{array}{c}{c}_{\mathrm{3}}=\frac{\mu \left(c_{\mathrm{8}}-c_{\mathrm{5}}\right)}{\mathrm{2}},\hspace{1em}\hspace{1em}{c}_{\mathrm{7}}=\frac{\eta \left(c_{\mathrm{8}}-c_{\mathrm{5}}\right)}{\mathrm{6}},\\ c_{\mathrm{6}}=\frac{c_{\mathrm{9}}}{\mathrm{2}},\hspace{1em}\hspace{1em}{c}_{\mathrm{10}}=\frac{\delta (c_{\mathrm{8}}-c_{\mathrm{5}})}{\mathrm{2}},\end{array}$

(12)

where

$\begin{array}{c}\mu =\frac{B_{\mathrm{11}}{A}_{\mathrm{55}}}{D_{\mathrm{11}}{A}_{\mathrm{11}}-B_{\mathrm{11}}^{\mathrm{2}}},\hspace{1em}\hspace{1em}\eta =\frac{A_{\mathrm{11}}{A}_{\mathrm{55}}}{D_{\mathrm{11}}{A}_{\mathrm{11}}-B_{\mathrm{11}}^{\mathrm{2}}},\\ \delta =\frac{\mathrm{1}}{\mathrm{12}+\eta l_e^{\mathrm{2}}}.\end{array}$

(13)

So, the longitudinal, transverse displacement and angular rotation in terms of field variables can be written as

$\left\{u\right\}=\left[N\left(x\right)\right]\left\{a\right\}.$

(14)

In the above equation [N(x)] is matrix having functions of x and {a} is vector of independent constants. Applying the boundary conditions at nodes 1 and 2, respectively, depending upon the end supports, shape functions for each case can be obtained as follows:

${\left[G\right]}^{-\mathrm{1}}=\left[\begin{array}{c}N\left(\mathrm{1}\right)\\ N\left(\mathrm{2}\right)\end{array}\right],\hspace{1em}\hspace{1em}\left\{a\right\}=\left[G\right]\left\{\widehat{u}\right\}.$

(15)

Here, $\left\{\widehat{u}\right\}$ is the nodal displacement vector. So, further

$\left\{u\right\}=\left[N\left(x\right)\right]\left\{a\right\}=\left[N\left(x\right)\right]\left[G\right]\left\{\widehat{u}\right\}=\left[\psi \left(x\right)\right]\left\{\widehat{u}\right\}.$

(16)

The ψ(x) represents matrix of exact shape functions. For this case, ψ ^u (x), ψ ^w (x), and ψ^ϕ(x) give the exact shape functions for axial displacement, transverse displacement, and angular rotation, respectively. The values of these shape functions are given below:

$\begin{array}{c}{\psi }_{\mathrm{11}}^u=\left(\mathrm{1}-\frac{x}{l_e}\right),\hspace{1em}\hspace{1em}{\psi }_{\mathrm{12}}^u=\mathrm{6}\mu \delta \left(\frac{x}{l_e}-\mathrm{1}\right),\\ {\psi }_{\mathrm{13}}^u=\mathrm{3}\mu \delta x\left(x-l_e\right),\hspace{1em}\hspace{1em}{\psi }_{\mathrm{14}}^u=\frac{x}{l_e},\\ {\psi }_{\mathrm{15}}^u=-{\psi }_{\mathrm{12}}^u,\hspace{1em}\hspace{1em}{\psi }_{\mathrm{16}}^u={\psi }_{\mathrm{13}}^u,\hspace{1em}\hspace{1em}{\psi }_{\mathrm{11}}^w=\mathrm{0},\\ {\psi }_{\mathrm{12}}^w=\frac{\delta }{l_e}\left(\eta l_e^{\mathrm{3}}+\mathrm{12}{l}_e-\mathrm{12}x+\mathrm{2}\eta x^{\mathrm{3}}-\mathrm{3}\eta x^{\mathrm{2}}{l}_e\right),\\ {\psi }_{\mathrm{13}}^w=\frac{x\delta }{l_e}(\eta l_e^{\mathrm{3}}+\mathrm{6}{l}_e-\mathrm{6}x+\eta x^{\mathrm{2}}{l}_e),\\ {\psi }_{\mathrm{14}}^w=\mathrm{0},\hspace{1em}\hspace{1em}{\psi }_{\mathrm{15}}^w=\frac{x\delta }{l_e}\left(-\mathrm{12}+\mathrm{2}\eta x^{\mathrm{2}}-\mathrm{2}\eta x{l}_e^{\mathrm{2}}\right),\\ {\psi }_{\mathrm{16}}^w=\frac{x\delta }{l_e}(-\mathrm{6}{l}_e+\mathrm{6}\eta x^{\mathrm{2}}{l}_e-\eta x{l}_e^{\mathrm{2}}),\\ {\psi }_{\mathrm{11}}^{\phi }=\mathrm{0},\hspace{1em}\hspace{1em}{\psi }_{\mathrm{12}}^{\phi }=\frac{\mathrm{6}x\delta \eta }{l_e}\left(-l_e+x\right),\\ {\psi }_{\mathrm{13}}^{\phi }=\frac{\delta }{l_e}\left(\eta l_e^{\mathrm{3}}+\mathrm{12}{l}_e-\mathrm{12}x+\mathrm{3}\eta x^{\mathrm{2}}{l}_e-\mathrm{4}\eta x{l}_e^{\mathrm{2}}\right),\\ {\psi }_{\mathrm{14}}^{\phi }=\mathrm{0},\hspace{1em}\hspace{1em}{\psi }_{\mathrm{15}}^{\phi }=-{\psi }_{\mathrm{12}}^{\phi },\\ {\psi }_{\mathrm{16}}^{\phi }=\frac{x\delta }{l_e}\left(\mathrm{12}+\mathrm{3}\eta x{l}_e-\mathrm{2}\eta l_e^{\mathrm{2}}\right).\end{array}$

(17)

Finally, the equation of motion is derived as follows:

$\left[M\right]\left\{\ddot{q}\right\}+\left[K\right]\left\{q\right\}=\left\{F\right\},$

(18)

where [M] denotes the global mass matrix and [K] denotes the global stiffness matrix. Initially, elemental stiffness and mass matrices need to find out.

3. Alidation of Present Theory

The elemental stiffness matrix can be obtained as given

$K^e={{\displaystyle \int }}_{\mathrm{0}}^{l_e}{\left[B\right]}^T\left[D\right]\left[B\right]dx.$

(19)

In this equation [B] is strain displacement matrix and [D] is constitutive matrix. Integrating the above equation the elemental stiffness matrix can be obtained. The stiffness matrix is symmetric, so the terms in upper triangular matrix are given below and the remaining terms are zero:

$\begin{array}{c}{K}_{\mathrm{11}}^e=-K_{\mathrm{14}}^e=\frac{A_{\mathrm{11}}^e}{l_e},\hspace{1em}\hspace{1em}{K}_{\mathrm{13}}^e=-K_{\mathrm{16}}^e=-\frac{B_{\mathrm{11}}^e}{l_e},\\ K_{\mathrm{22}}^e=-K_{\mathrm{25}}^e=\frac{A_{\mathrm{55}}^e{\psi }^e}{l_e},\hspace{1em}\hspace{1em}{K}_{\mathrm{23}}^e=K_{\mathrm{26}}^e=\frac{A_{\mathrm{55}}^e{\psi }^e}{\mathrm{2}},\\ K_{\mathrm{33}}^e=K_{\mathrm{66}}^e=\frac{D_{\mathrm{11}}^e}{l_e}+\frac{A_{\mathrm{55}}^e{\psi }^e{l}_e}{\mathrm{4}},\hspace{1em}\hspace{1em}{K}_{\mathrm{44}}^e={-K}_{\mathrm{11}}^e,\\ K_{\mathrm{46}}^e={-K}_{\mathrm{13}}^e,\hspace{1em}\hspace{1em}{K}_{\mathrm{55}}^e={-K}_{\mathrm{22}}^e,\hspace{1em}\hspace{1em}{K}_{\mathrm{56}}^e={-K}_{\mathrm{23}}^e.\end{array}$

(20)

Similarly, elemental stiffness matrix can be obtained by using following formula:

$M^e={{\displaystyle \int }}_{\mathrm{0}}^{l_e}\left[\rho \right]{\left[N\right]}^T\left[N\right]dx.$

(21)

Based upon the present theory a Matlab code has been developed and implemented to a beam with h = 0.1 m, L/h = 10 to study the convergence of nondimensional fundamental frequency for V_cnt^* = 0.12. It is calculated by the following formula:

$\Omega =\frac{\omega l^{\mathrm{2}}}{h}\sqrt{\frac{{\rho }_m}{E_m}}.$

(22)

The material properties of the beam are as given in Table 1.

Table 1:

Material properties of beam.

Material property	CNT	Matrix
Young's modulus (GPa)	600	2.5
Poisson's ratio	0.19	0.3
Density (Kg/m3)	1130	1190

Table 2 shows the nondimensional frequency for different numbers of elements.

Table 2:

Variation of nondimensional frequency for different numbers of elements.

Number of elements	Nondimensional frequency (Ω)
2	1.8193
5	1.6748
10	1.6540
15	1.6501
20	1.6488
21	1.6486

The graph is plotted between number of elements along X-axis and nondimensional frequency along Y-axis.

From the Table 2 and Figure 2 it is found that the result converges for 21 numbers of elements. So, for further analysis 21 elements were used. Now, present theory is validated for V_cnt^* = 0.12, 0.17, and 0.28. Results are obtained with both ends of the beam clamped and it is compared with the results given by Ke at al. [9]. Table 3 compares present nondimensional fundamental frequencies with the results available in the literature. From the table it is found that the results agree well with reference.

Table 3:

Comparison of nondimensional fundamental frequencies.

S. number	V cnt *	Ωref [9]	Ωpresent
1	0.12	1.6063	1.6490
2	0.17	2.0557	2.0898
3	0.28	2.3325	2.2256

Figure 2:

Convergence study.

Proposed theory is also validated by comparing static deflection of free end of cantilever beam with that presented by Rokni et al. [17] (Table 5). The beam dimensions are L = 0.14 m, b = 0.02 m, and h = 0.01 m. For this case the assumed material properties are given in Table 4.

Table 4:

Material properties of beam.

Material property	CNT	Matrix
Young's modulus (GPa)	900	1.9
Poisson's ratio	0.28	0.34
Density (Kg/m3)	2100	1050

Table 5:

Comparison of results for cantilever beam.

S. Number	Property	Reference	Present
1	Effective Young's modulus E (Gpa)	3.5979	3.4347
2	Deflection at free end (mm)	2.9897	3.1316

The deflection at the free end of the beam is calculated by following formula:

$\delta =\frac{F{L}^{\mathrm{3}}}{\mathrm{3}EI}.$

(23)

4. Results

Now, the results are obtained for above beam with simply supported end conditions for varying volume fraction of CNT. The effect of variation of volume fraction of CNT on nondimensional fundamental frequency is shown in Table 6 for slenderness ratio that is equal to 10 in all cases.

Table 6:

Effect of variation of volume fraction of CNT on nondimensional fundamental frequency.

S. number	CNT volume fraction Vcnt*	Nondimensional frequency (Ω)
1	0.1	0.8118
2	0.2	0.9293
3	0.3	1.0307

From Table 6 it is found that nondimensional frequency increases as volume fraction of CNT increases. Now, the results are obtained to observe the effect of slenderness ratio on dimensional frequency for constant CNT volume fraction. Table 7 shows effect of slenderness ratio (L/h) on nondimensional fundamental frequency with V_cnt = 0.1. It shows that as L/h ratio increases frequency also increases.

Table 7:

Effect of slenderness ratio on nondimensional fundamental frequency.

S. number	(L/h) ratio	Nondimensional frequency (Ω)
1	10	0.5805
2	20	1.1632
3	30	1.7506
4	40	2.3450
5	50	2.949

Finally, mode shapes are plotted for first case of Table 7. First five modes are shown in Figure 3. Mode shapes are the vibration patterns for the particular frequency and they are obtained from eigenvectors.

Figure 3:

Five modes shapes corresponding to the first five natural frequencies.

Now, the dynamic analysis is carried out for given beam. Displacement versus time diagram is plotted to understand the time response with varying CNT volume fraction. For simply supported case, the maximum displacement obtained for V_cnt = 0.1, 0.2, and 0.3 are 0.029 m, 0.022 m, and 0.016 m, respectively. Figure 4 shows the response of sinusoidal force of magnitude 100sin ⁡ (80t) for the same case of mode shapes for 0.45 seconds in the tip of beam. As Figure 4 shows, if V_cnt increases, maximum displacement decreases.

Figure 4:

Displacement versus time graph for V_cnt = 0.1, 0.2, and 0.3.

5. Effect of Carbon Nanotube Orientation

In the previous analysis the CNTs are assumed to be unidirectional along the length. But practically, it is not possible. The effect of CNT orientation on natural frequency is discussed in this section.

5.1. Material Properties of FG-CNTRC (Carbon Nanotube Reinforced Composite)

An embedded carbon nanotube in a polymer matrix is considered to be replaced with an equivalent long fiber for predicting the mechanical properties of the carbon nanotube/polymer composite. The inverse rule of mixture is used for calculating material properties of equivalent fiber [7]:

$\begin{array}{c}{E}_{LEF}=\frac{E_{LC}}{V_{EF}}-\frac{E_m{V}_m}{V_{EF}},\\ \frac{\mathrm{1}}{G_{EF}}=\frac{\mathrm{1}}{G_C{V}_{EF}}-\frac{V_m}{G_m{V}_{EF}},\\ {\nu }_{EF}=\frac{V_C}{V_{EF}}-\frac{{{\nu }_{m}V}_m}{V_{EF}},\end{array}$

(24)

where E, G, ν, and V are the longitudinal Young's modulus, shear modulus, Poisson's ratio, and volume fraction, respectively. Suffixes EF, C, and m show properties for equivalent fiber, CNT, and matrix material, respectively. The equivalent fiber for SWCNT with chiral index of (10, 10) is a solid cylinder with diameter of 1.424 nm. Material properties for this CNT taken from reference are mentioned in Table 8.

Table 8:

Material properties of equivalent fiber.

Mechanical property	Equivalent fiber [9]
Longitudinal Young's modulus (E LEF )	649.12 (Gpa)
Transverse Young's modulus (E TEF )	11.27 (Gpa)
Longitudinal shear modulus (G EF )	5.13 (Gpa)
Poisson's ratio (υ EF )	0.284

5.2. Modelling of Beam

The effect of randomly oriented, straight CNTs is formulated here. The orientation of a straight CNT is characterized by two Euler angles α and β, as shown in Figure 5. Thus giving rise to the base vectors $\overrightarrow{e_i}$ of the global $\left(\begin{array}{cccc}0& -x_1& x_2& x_3\end{array}\right)$ and $\overrightarrow{e_i^{\prime }}$ of the local coordinate systems $\left(\begin{array}{cccc}0& -x_1^{\prime }& x_2^{\prime }& x_3^{\prime }\end{array}\right)$ which are related through the transformation matrix g:

$\overrightarrow{e_i}=g\overrightarrow{e_i^{\prime }},$

(25)

where g is given by

$g=\left[\begin{array}{ccc}\cos \beta & -\cos \alpha \sin \beta & \sin \alpha \sin \beta \\ \sin \beta & \cos \alpha \cos \beta & -\sin \alpha \cos \beta \\ \mathrm{0}& \sin \alpha & \cos \alpha \end{array}\right].$

(26)

The orientation distribution of CNTs in a composite is characterized by a probability density function p(α, β) satisfying the normalizing condition:

${{\displaystyle \int }}_{\mathrm{0}}^{\mathrm{2}\pi }{{\displaystyle \int }}_{\mathrm{0}}^{\pi /\mathrm{2}}p\left(\alpha ,\beta \right)\sin \alpha \hspace{0.17em}d\alpha \hspace{0.17em}d\beta =\mathrm{1}.$

(27)

Consider CNTs to be completely randomly oriented, the density function for this case is

$p\left(\alpha ,\beta \right)=\frac{\mathrm{1}}{\mathrm{2}\pi }.$

(28)

The Hill's elastic module of the reinforcing phase can be calculated from the equality of following two matrices [7]:

$C_r\hspace{0.5em}=\left[\begin{array}{cccccc}{n}_r& l_r& l_r& \mathrm{0}& \mathrm{0}& \mathrm{0}\\ l_r& k_r+m_r& k_r-m_r& \mathrm{0}& \mathrm{0}& \mathrm{0}\\ l_r& k_r-m_r& k_r+m_r& \mathrm{0}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& \mathrm{0}& \mathrm{0}& p_r& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& \mathrm{0}& \mathrm{0}& \mathrm{0}& m_r& \mathrm{0}\\ \mathrm{0}& \mathrm{0}& \mathrm{0}& \mathrm{0}& \mathrm{0}& p_r\end{array}\right]\hspace{0.5em}={\left[\begin{array}{cccccc}\frac{\mathrm{1}}{E_{ELF}}& \frac{-{\nu }_{EF}}{E_{ELF}}& \frac{-{\nu }_{EF}}{E_{ELF}}& \mathrm{0}& \mathrm{0}& \mathrm{0}\\ \frac{-{\nu }_{EF}}{E_{ELF}}& \frac{\mathrm{1}}{E_{ELF}}& \frac{-{\nu }_{EF}}{E_{ELF}}& \mathrm{0}& \mathrm{0}& \mathrm{0}\\ \frac{-{\nu }_{EF}}{E_{ELF}}& \frac{-{\nu }_{EF}}{E_{ELF}}& \frac{\mathrm{1}}{E_{ELF}}& \mathrm{0}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& \mathrm{0}& \mathrm{0}& \frac{\mathrm{1}}{G_{EF}}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& \mathrm{0}& \mathrm{0}& \mathrm{0}& \frac{\mathrm{1}}{G_{EF}}& \mathrm{0}\\ \mathrm{0}& \mathrm{0}& \mathrm{0}& \mathrm{0}& \mathrm{0}& \frac{\mathrm{1}}{G_{EF}}\end{array}\right]}^{-\mathrm{1}},$

(29)

where C _r is the stiffness tensor of the equivalent fibers and when CNTs are assumed randomly oriented in the matrix, the composite is assumed to be isotropic then bulk modulus K and shear modulus G of the beam are derived as:

$\begin{array}{c}K=K_m+\frac{V_{\text{cnt}}({\delta }_r-\mathrm{3}{K}_m{\alpha }_r)}{\mathrm{3}\left(V_m{+V}_{\text{cnt}}{\alpha }_r\right)},\\ G=G_m+\frac{V_{\text{cnt}}({\eta }_r-\mathrm{2}{G}_m{\beta }_r)}{\mathrm{2}\left(V_m{+V}_{\text{cnt}}{\beta }_r\right)}.\end{array}$

(30)

Finally, Young's modulus for oriented CNT considering modified shear and bulk modulus will be

$E=\frac{\mathrm{9}KG}{\mathrm{3}K+G}.$

(31)

Keeping all other values constant and using Young's modulus from the above equation, (7), (9), and (24) are modified to obtain the effect of CNT orientation on elemental stiffness and mass matrices.

Figure 5:

Representative volume element (RVE) with randomly oriented, straight CNTs [7].

As CNT volume is assumed to vary along the thickness only, material properties for each layer are calculated first and finally effective values for entire beam are calculated. Here, linear variation of volume fraction of symmetrical CNT V_cnt is considered. It is calculated by following formula:

$V_{\text{cnt}}=\frac{\mathrm{4}\left|z\right|}{h}{V}_{\text{cnt}}^{*}.$

(32)

5.3. Numerical Results

This theory has been implemented to a beam with h = 1 m, L/h = 20 to study the convergence of nondimensional fundamental frequency for V_cnt^* = 0.075. It is calculated by following formula:

${\lambda }^{\mathrm{2}}=\omega l^{\mathrm{2}}\sqrt{\frac{{\rho }_{m}A}{E_{m}I}}.$

(33)

Now, the material properties assumed for the beam are as mentioned in the Table 9.

Table 9:

Material properties of beam.

Material property	CNT	Matrix
Young's modulus (Gpa)	900	10
Poisson's ratio	0.28	0.3
Density (Kg/m3)	2100	1150

With these properties and applying the present theory, results of the first five nondimensional frequencies for clamped-clamped (C-C) are obtained. The effect on SFG-CNTR beam based on Timoshenko beam theory with different number of elements is obtained to study the convergence and is shown in Table 10. The Graph is also plotted to observe the convergence of frequencies as shown in Figure 6.

Table 10:

Effect of number of elements on nondimensional fundamental frequency.

Mode number	20	40	60	80	100
1	5.7046	5.4562	5.4229	5.4089	5.4030
2	7.7032	7.728	7.6686	7.6506	7.6413
3	9.6638	9.4572	9.3968	9.3708	9.3599
4	9.8956	10.9441	10.8519	10.8237	10.809
5	11.7088	12.2242	12.1418	12.1037	12.0874

Figure 6:

Convergence study of frequencies with consideration of CNT orientation.

It is observed from Figure 6 and Table 10 that the convergence of the present results occurs with a number of elements N = 100.

To validate the results beam slenderness ratio L/h = 20 and V_cnt^* = 0.075 are selected with clamped-clamped condition of Timoshenko beam and are verified with the results given by Yas and Heshmati [7].

Table 11 compares the effect of variation of volume fraction of CNT on dimensionless fundamental frequency. The first five frequencies for clamped-clamped, SFG-CNTRC (carbon nanotube reinforced composite) beams are shown in the table.

Table 11:

Effect of variation of volume fraction of CNT on nondimensional fundamental frequency.

Mode number	Ref [7]	Present study (CNT volume fraction (Vcnt))
	.075	0.075	0.1	0.125	0.15	0.175
1	5.031699	5.0403	5.7963	6.1273	6.4122	6.663
2	7.85222	7.6413	8.1981	8.6656	9.0706	9.4244
3	10.3392	9.3599	10.0411	8.6656	11.108	11.5424
4	—	10.809	11.5967	12.2579	12.8307	13.3311
5	—	12.0874	12.9672	13.7077	14.3448	14.9053

From this table it is clear that the values of natural frequencies increase due to the effect of CNT orientation. Table 12 shows the effect of slenderness ratio (L/h) on nondimensional fundamental frequency. Frequencies are calculated for slenderness ratios 20, 40, 60, and 80. It is found that as L/h ratio increases frequency also increases.

Table 12:

Effect of slenderness ratio on nondimensional fundamental frequency.

Mode number	L/h = 20	L/h = 40	L/h = 60	L/h = 80
1	5.403	7.6883	9.5166	21.2886
2	7.6413	7.6883	13.4554	15.7780
3	9.3599	13.3185	16.4855	19.3357
4	10.809	15.379	19.0325	22.3165
5	12.0874	17.1994	21.2886	24.9678

Displacement versus time diagram is plotted to understand the time response with varying CNT volume fraction. Figure 7 shows the response of sinusoidal force of magnitude 100 N which is acting at the centre of the beam with clamped-clamped boundary condition of FG beam and exciting frequency of 80 Hz. The maximum displacement in this case is found to be 0.09 micrometers.

Figure 7:

Plot of displacement versus time of CNT based on FG beam under sinusoidal loading.

6. Conclusions

In the present work, CNT based on FG Timoshenko beam has been modeled using finite element method. Firstly convergence of results has been studied and then different analyses such as static and dynamic analyses have been carried out. Nondimensional fundamental frequencies are calculated for different volume fractions of CNT and slenderness ratio. Mode shapes and displacement diagrams are also plotted for this beam. From this analysis it can be concluded that

• as volume fraction of CNT increases, fundamental frequency also increases. The reason behind this is more volume of CNT provides more stiffness to the beam which results in higher frequencies.

• Nondimensional fundamental frequency increases as slenderness ratio increases.

• The response of displacement versus time shows that for increasing values of volume fraction displacement reduces.

Conflict of Interests

The authors declare that there is no conflict of interests regarding the publication of this paper.