# An Inhomogeneous Cell-Based Smoothed Finite Element Method for Free Vibration Calculation of Functionally Graded Magnetoelectroelastic Structures

Authors: Yan Cai; Guangwei Meng; Liming Zhou
Journal: Shock and Vibration
 (2018)

---

## Abstract

To overcome the overstiffness and imprecise magnetoelectroelastic coupling effects of finite element method (FEM), we present an inhomogeneous cell-based smoothed FEM (ICS-FEM) of functionally graded magnetoelectroelastic (FGMEE) structures. Then the ICS-FEM formulations for free vibration calculation of FGMEE structures were deduced. In FGMEE structures, the true parameters at the Gaussian integration point were adopted directly to replace the homogenization in an element. The ICS-FEM provides a continuous system with a close-to-exact stiffness, which could be automatically and more easily generated for complicated domains, thus significantly decreasing the numerical error. To verify the accuracy and trustworthiness of ICS-FEM, we investigated several numerical examples and found that ICS-FEM simulated more accurately than the standard FEM. Also the effects of various equivalent stiffness matrices and the gradient function on the inherent frequency of FGMEE beams were studied.


---

## Body

## 1. Introduction

Functionally graded magnetoelectroelastic (FGMEE) materials are generally multiphase composites with continuously varying mechanical properties. FGMEE materials can convert magnetic, electric, and mechanical energy from one type into another and have received wide attention recently due to their electroelastic, magnetoelastic, and electromagnetic coupling effects [1, 2]. Therefore, FGMEE materials have been adopted in various smart structures, such as magnetic field probes, smart vibration sensors, optoelectronic devices, and medical ultrasonic transducers [3, 4]. The smart FGMEE structures are commonly fabricated in the beam pattern. However, to better apply FGMEE beams, researchers must analyse the statics and free vibration, and to predict the coupled response of FGMEE beams for practical applications, they should accurately calculate the properties of free vibrations.Several computational techniques were proposed to investigate the electroelastic, magnetoelastic, and electromagnetic coupling effects of smart structures, such as finite element method (FEM), mesh-free method, and scaled boundary FEM [5–10]. Bhangale and Ganesan analyzed the static behaviors of linear anisotropic FGMEE plates using semianalytical FEM and investigated the free vibration of FGMEE plates and cylindrical shells [11, 12]. A layerwise partial mixed FEM was proposed to model MEE plates [13]. Phoenix et al. analysed the static and dynamic behaviors of coupled MEE plates using FEM with the Reissner mixed variational theorem [14]. Buchanan used FEM to study the free vibrations of infinite magnetoelectroelastic cylinders [15]. However, these FEMs overestimated the stiffness of solid structures and were limited by low accuracy. Therefore, Sladek et al. proposed a mesh-free method to more accurately study the static behavior of a circular FGMEE plate [16], but this method reduced the computational effectiveness. Recently, Liu et al. solved the deformations of a nonuniform MEE plate using scaled boundary FEM [17]. By incorporating the nonlocal theory into scaled boundary FEM, Ke and Wang more accurately and effectively studied the free vibrations of MEE beams [18]. However, the effectiveness of scaled boundary FEM is still low and should be improved.FEM, as a powerful computational tool to investigate MEE coupling behaviors, yet overestimates the stiffness of FGMEE structures, which may result in locking behavior and inaccurate eigenvalue solutions [19]. To overcome these limitations, a series of cell-based smoothed FEM (CS-FEM) [20–26] and node-based, edge-based, or face-based smoothed FEMs [27–32] were proposed. In recent years, many CS-FEM-based formulations were proposed. Moreover, CS-FEM does not require the shape function derivatives or high generosity of program and is insensitive to mesh distortion because of the absence of isoparametric mapping.CS-FEM has been successfully extended into dynamical control of piezoelectric sensors and actuators, topological optimization of linear piezoelectric micromotor, and analysis of static behaviors, frequency, and defects of piezoelectric structures [33–43]. Due to its versatility, CS-FEM becomes a simple and effective numerical tool to solve numerous electric and mechanical physical problems. However, the application of CS-FEM to investigate MEE properties is still a challenge.In this work, free vibrations of FGMEE structures were studied. Inhomogeneous CS-FEM (ICS-FEM) for FGMEE materials was formulated by incorporating gradient smoothing into the standard FEM for the multi-physics field of FGMEE. Then the equations of free vibration computation were deduced under the multi-physics coupling field for FGMEE materials. Finally, FGMEE beams in the functional gradient exponential form or power law form were calculated under different boundary conditions. ICS-FEM outperformed FEM when compared with the reference solution.This paper is organized as follows: Section2 introduces the basic formulations for FGMEE materials. Section 3 briefly describes the properties of the FGMEE materials. Section 4 contains the detailed formulation of ICS-FEM. In Section 5, two numerical examples and the model of typical MEMS-based FGMEE energy harvester are investigated in detail. Final conclusions from the numerical results are drawn in Section 6.


## 2. Basic Formulations

The material properties of a functionally graded material (FGM) plate vary continuously, which is an advantage over the discontinuity across adjoining layers in a laminated plate. The wide range of engineering applications of FGM has attracted many scientists to investigate the behaviors of FGM.Considering the transverse isotropy of the FGMEE medium [9, 44] and for the plane stress problem, we set stress components σy = σyz = τxy = 0, electric displacement component Dy = 0, and magnetic induction component By = 0. The geometric parameters and the chosen Cartesian coordinate system (x,y,z) are illustrated in Figure 1.Figure 1
FGMEE beams: Cartesian coordinate system and geometric parameters.The basic formulations for MEE materials include equilibrium equations, geometric equations, and constitutive equations. The equilibrium equations are(1)∂σx∂x+∂τxz∂z=0,∂τxz∂x+∂σz∂z=0,∂Dx∂x+∂Dz∂z=0,∂Bx∂x+∂Bz∂z=0,where σx,σz, and τxz denote stress components; Dx and Dz are electric displacement components; Bx and Bz are magnetic induction components.The geometric equations are(2)Sx=∂u∂x,Sz=∂w∂z,Sxz=∂w∂x+∂u∂z,Ex=-∂Φ∂x,Ez=-∂Φ∂z,Hx=-∂Ψ∂x,Hz=-∂Ψ∂z,where Sx, Sz, and Sxz denote strain components; u and w are displacement components; Ex and Ez are electric field components; Φ is electrical potential; Hx and Hz are magnetic field components; Ψ is magnetic potential.The constitutive equations are(3)σxσzτxz=C11C130C13C33000C44SxSzSxz+0e310e33e150ExEz+0q310q33q150HxHz,DxDz=00e15e31e330SxSzSxz+ε1100ε33ExEz+m1100m33HxHz,BxBz=00q15q31q330SxSzSxz+m1100m33ExEz+μ1100μ33HxHz,where Cij, εij, and μij are the elastic, dielectric, and magnetic permeability coefficients, respectively; eij, qij, and mij are piezo-electric, piezomagnetic, and magnetoelectric coefficients, respectively. For FGMEE materials, we have(4)Cij=Cij0fz,εij=εij0fz,μij=μij0fz,eij=eij0fz,qij=qij0fz,mij=mij0fz,where f(z) is an arbitrary function; Cij0 and εij0, μij0, eij0, qij0, and mij0 are the values on the plane z = 0.


## 3. FGMEE Materials

An FGMEE structure is characterized by the high heterogeneity of material properties with a distribution prescribing the volume fractions of constituent phases. For particular analysis, it is functional to idealize them as continua with smooth gradual variation of material properties in the spatial coordinates. Hence, the proper micromechanical model should be able to characterize the material property distribution of a system in accurate sense.Previous literatures focus on two types of gradation methods widely applied to solve many problems. Among various methods for composites, some are also used for FGMEE materials, including the exponential and Voigt rule of mixture scheme.For FGMEE materials with exponential variation in the thickness direction (z-direction), (4) can be rewritten as(5)Cij=Cij0eηz/h,εij=εij0eηz/h,μij=μij0eηz/h,eij=eij0eηz/h,qij=qij0eηz/h,mij=mij0eηz/h,where η is the exponential factor governing the degree of z-direction gradient, h is the thickness, the superscript 0 indicates the z-independent coefficients, and η = 0 in homogeneous MEE materials.The volume fraction of an FGMEE structure across the thickness direction is assumed as a simple power law type as follows:(6)VB=2z+h2hn,where -h/2≤z≤h/2 and n is the power law index. The bottom surface of the material (z=-h/2) is VC whereas the top surface (z=h/2) is VB. The total volume of the constituents should be (7)VC+VB=1.Based on (6) and (7), the effective material property is defined as follows: (8)MCeff=MCtopVB+MCbottomVC,where “MC” is general notation for material property. With (3), the effective coefficients can be written as(9)Ceff=(_CB-CC)_VB+CC,εeff=(_εB-εC)_VB+εC,μeff=(_μB-μC)_VB+μC,eeff=(_eB-eC)_VB+eC,qeff=(_qB-qC)_VB+qC,meff=(_mB-mC)_VB+mC,where “eff” stands for effective properties corresponding to a specific value of n. Material coefficients of piezoelectric BaTiO3, magnetostrictive CoFe2O4, and MEE BatiO3-CoFe2O4 are given in Table 1. Figure 2 depicts the through-the-thickness distribution of the volume fraction changing with different values of n. For n = 1.0, the variation of effective material property is linear.Table 1
Magnetoelectroelastic coefficients of material properties [6].
Material constants CoFe2O4 BatiO3-CoFe2O4 BatiO3 C 11 109 N/m2 286 200 166 C 12 109 N/m2 173 110 77 C 13 109 N/m2 170 110 78 C 33 109 N/m2 269.5 190 162 C 44 109 N/m2 45.3 45 43 e 31 C/m2 0 −3.5 −4.4 e 33 C/m2 0 11 18.6 e 15 C/m2 0 0 11.6 ε 11 10−9 C/Vm 0.08 0.9 11.2 ε 33 10−9 C/Vm 0.093 7.5 12.6 μ 11 10−4 Ns2/C2 −5.9 −1.5 0.05 μ 33 10−4 Ns2/C2 1.57 0.75 0.1 q 31 N/Am 580 200 0 q 33 N/Am 700 260 0 q 15 N/Am 560 180 0 m 11 10−12 Ns/VC 0 6.0 0 m 33 10−12 Ns/VC 0 2500 0 ρ kgm−3 5730 5730 5730Figure 2
Variation of the volume fraction function versus the nondimensional thicknessz/h with varying n.


## 4. ICS-FEM

The solution domainΩ is discretized into np elements containing Nn nodes, the approximation displacement u-, the approximation electrical potential Φ-, and the approximation magnetic potential Ψ-. For FGMEE materials, we have(10)u-=∑i=1npNiuui=Nuu,Φ-=∑i=1npNiΦΦi=NΦΦ,Ψ-=∑i=1npNiΨΨi=NΨΨ,whereu, Φ, and Ψ are the vectors of node displacement, node electrical potential, and node magnetic potential, respectively; Nu, NΦ, and NΨ are displacement shape, electrical potential shape, and magnetic potential shape functions of ICS-FEM, respectively. Nu, NΦ, and NΨ were expressed in similar shape functions. Four-node element divided into four smoothing subdomains [27], field nodes, edge smoothing nodes, center smoothing nodes, edge Gaussian point, outer normal vector distribution, and shape function values are shown in Figure 3.Figure 3
Smoothing subdomains and the values of shape functions.At any pointxk in the smoothing subdomain Ωik, the smoothed strain ε-(xk), smoothed electric field E-(xk), and smoothed magnetic field H-(xk) are(11)S-xk=∫ΩikSxκx-xkdΩ,E-xk=∫ΩikExκx-xkdΩ,H-xk=∫ΩikHxκx-xkdΩ,where Sx, Ex, and Hx are the strain, electric field, and magnetic field in FEM, respectively; κx-xk is the constant function:(12)κx-xk=1Aikx∈Ωik0x∉Ωik,where(13)Aik=∫ΩikdΩ.Substituting (12) into (11), we get(14)S-xk=1Aik∫ΓiknukudΓ,E-xk=1Aik∫ΓiknΦkΦdΓ,H-xk=1Aik∫ΓiknΨkΨdΓ,where Γik is the boundary of Ωik; nuk, nΦk, and nΨk are the outer normal vector matrices of the boundary:(15)nuk=nxk00nzknzknxk,nΦk=nxknzk,nΨk=nxknzk.Eqs. (14) can be rewritten as(16)S-xk=∑i=1neB-uixkui,E-xk=-∑i=1neB-ΦixkΦi,H-xk=-∑i=1neB-ΨixkΨi,where ne is the number of smoothing elements(17)B-uixk=1Aik∫ΓkNiunxk00NiunzkNiunzkNiunxkdΓ,B-Φixk=1Aik∫ΓkNiΦnxkNiΦnzkdΓ,B-Ψixk=1Aik∫ΓkNiΨnxkNiΨnzkdΓ.At the Gaussian point of the smoothing boundaryxbG, (17) are rewritten as(18)B-uixk=1Aik∑b=1nbNiuxbGnxk00NiuxbGnzkNiuxbGnzkNiuxbGnxklbk,B-Φixk=1Aik∑b=1nbNiΦxbGnxkNiΦxbGnzklbk,B-Ψixk=1Aik∑b=1nbNiΨxbGnxkNiΨxbGnzklbk,where lbk is the length of the smoothing boundary; nb is the total number of boundaries for each smoothing subdomain.As for the essential difference, FEM has to derive the shape function matrix of the element, but ICS-FEM avoids this step and simply uses the shape function atxbG, which reduces the requirement for continuity of the shape function and improves the accuracy and convergence by the use of gradient smoothing.The thermodynamic potential of a 2D FGMEE problem is given as(19)G=GS,E,H,whereS,E, andH are independent variables of strain, electric field, and magnetic field, respectively.By applying (1) into (19), we get the variational expression of MEE plane:(20)G=12STCS-12ETεE-12HTμH-SeE-SqH-EmH.By minimizing (20) for nodal variables of shape functions for strain-displacement, electric field–electric potential, and magnetic field–magnetic potential, we get the ICS-FEM equations for MEE plane:(21)Kuu-ω2Mu+KuΦΦ+KuΨΨ=0,KuΦTu-KΦΦΦ-KΦΨΨ=0,KuΨTu-KΦΨTΦ-KΨΨΨ=0,where ω is the eigenvalues.Different elemental stiffness matrices used for FGMEE beams are expressed as follows:(22)Kuu=∑i=1ncB-uiTCB-uiAik,KuΦ=∑i=1ncB-uiTeB-ΦiAik,KuΨ=-∑k=1ncB-uiTqB-ΨiAik,KΦΨ=∑i=1ncB-ΦiTmB-ΨiAik,KΦΦ=∑i=1ncB-ΦiTεB-ΦiAik,KΨΨ=∑i=1ncB-ΨiTμB-ΨiAik,M=∑eMe,Me=diagm1,m1,m2,m2,m3,m3,m4,m4,where nc=np×ne; mi=ρitAik(i=1,2,3,4) is the mass of smoothing element i; t is the smoothing element thickness; ρi is density of Gaussian integration point in smoothing subdomain i; C, ε, μ, e, q, and m are the matrices of elastic constant, dielectric coefficient, magnetic permeability, piezomagnetic coefficient, piezomagnetic coefficient, and magnetoelectric coefficient, respectively.The inhomogeneous smoothing element was adopted to calculate its stiffness matrix. Because the parameters of four smoothing subdomainsAik(i=1,2,3,4) differed in element s, the actual parameters at the Gaussian integration point were taken directly to simulate the changes of material property in each element.By eliminating the terms of electric and magnetic potentials using a condensation technique, we get the equivalent stiffness matrixKeq:(23)Kequ+Mu¨=0,where (24)Keq=Kuu+KuΦKII-1KI+KuΨKIV-1KIII,KI=KuΦT-KΦΨKΨΨ-1KuΨT,KII=KΦΦ-KΦΨKΨΨ-1KΦΨT,KIII=KuΨT-KΦΨTKΦΦ-1KuΦT,KIV=KΨΨ-KΦΨTKΦΦ-1KΦΨ.The eigenvectors corresponding toΦ and Ψ are given as(25)Φ=KII-1KIu,Ψ=KIV-1KIIIu.To study the effect of magnetoelectric constant on system frequencies, we derivedKeq_reduced by neglecting the magnetoelectric coupling effect.By makingKΦΨ=0, we get the reduced cell-based finite element equations:(26)Kuu-ω2Mu+KuΦΦ+KuΨΨ=0,KuΦTu-KΦΦΦ=0,KuΨTu-KΨΨΨ=0.The reduced stiffness matrixKeq_reduced is (27)Keq_reduced=Kuu+KuΦKΦΦ-1KuΦT+KuΨKΨΨ-1KuΨT.To evaluate the effect of PE phase on beam frequency, we derived the stiffness matrixKeq_ΦΦ by setting the magnetic potential = 0:(28)Keq_ΦΦ=Kuu+KuΦKΦΦ-1KuΦT.To study the magnetic effect of PM phase on system frequency, we obtainedKeq_ΨΨ by plugging electric potential to zero in (26):(29)Keq_ΨΨ=Kuu+KuΨKΨΨ-1KuΨT.


## 5. Results and Discussion

### 5.1. C-F Beam

The free vibrations on FGMEE beams were calculated by changing the exponential factor (Figure4). The material properties of FGMEE beams were governed by the z-direction exponential variation. The following geometrical parameters were considered: length L = 0.3 m and width H = 0.02 m with the assumption of plane stress. Boundary conditions were u = w = Φ = Ψ = 0 at the clamped end. Table 2 gives the various structural frequencies in the study.Table 2
List offused in the study.
Structural frequency Matrix used to compute the structural frequency f u u K u u f e q K e q f e q _ r e K e q _ r e d u c e d f e q _ Φ Φ K e q _ Φ Φ f e q _ Ψ Ψ K e q _ Ψ ΨFigure 4
Geometry of an FGMEE beam and the coordinates.Firstly, the convergence of ICS-FEM was verified by using BaTiO3-CoFe2O4 FGMEE beams, with properties listed in Table 1. The natural frequencies of these beams were calculated using ICS-FEM with different meshes (30 × 2, 60 × 4, 120 × 8, 150 × 10) (Figure 5). The simulation results with different meshes agree well, which prove the good convergence of ICS-FEM. The first- to third-order modes of clamp-free BaTiO3–CoFe2O4 FGMEE beams using different elements with exponential factor η = 1.0 and 5.0 were calculated by ICS-FEM, and the results were summarized in Figures 6 and 7, respectively. It was found that the first- to third-order modes of BaTiO3–CoFe2O4 FGMEE beams in the same gradient distribution were basically not affected by equivalent stiffness matrix or mesh number. the first- to third-order modes of BaTiO3–CoFe2O4 FGMEE beams were basically not different between η = 1.0 and 5.0.Figure 5
Structural frequencies (a)fuu; (b) feq; (c) feq-re; (d) feq_ΦΦ; (e) feq_ΨΨ of clamp-free BaTiO3–CoFe2O4 FGMEE beams.
(a) (b) (c) (d) (e)Figure 6
The first- to third-order modes of clamp-free BaTiO3–CoFe2O4 FGMEE beams using different elements with exponential factor η = 1.0 by ICS-FEM.Figure 7
The first- to third-order modes of clamp-free BaTiO3–CoFe2O4 FGMEE beams using different elements with exponential factor η = 5.0 by ICS-FEM.Secondly, the free vibration frequencies of CoFe2O4 FGMEE beams were studied by both ICS-FEM and FEM using extremely irregular elements, with domain discretization shown in Figure 8. The frequencies of clamp-free CoFe2O4 FGMEE beams with different values of exponential factor are shown in Figure 9; the first eleven natural frequencies calculated by ICS-FEM are smaller than those calculated by FEM. The validity of ICS-FEM is verified by the agreements between the calculations and the reference solutions. The shape of quadrilateral element in FEM cannot be severely distorted but was eliminated in ICS-FEM. ICS-FEM abstains from calculating the derivative of the shape functions of an element, and the area integral of the solution domain is converted to the boundary integral. The stiffness of FGMEE structures is improved because ICS-FEM does not require continuity of the shape function. The ICS-FEM provides a continuous system with a close-to-exact stiffness, which could be automatically and more easily generated for complicated domains, thus significantly decreasing the numerical error. The free vibration of CoFe2O4 FGMEE beam, a pure CoFe2O4 material without piezoelectric or magnetoelectric material coefficients, influences structural frequency feq because the magnetic effect is marginally higher compared with fuu. fΦΦ coincides with fuu since piezoelectric phase is absent in CoFe2O4. Similarly, feq_re coincides with feq of the CoFe2O4 FGMEE beam as the magnetoelectric effect is absent in pure CoFe2O4 FGMEE beam. The natural frequencies of CoFe2O4 FGMEE beams increase with the rise of the exponential factor.Figure 8
Domain discretization using four-node extremely irregular elements.Figure 9
Structural frequenciesf of clamp-free CoFe2O4 FGMEE beams with exponential factor η = (a) 0, (b) 0.5, (c) 1.0, (d) 2.0, and (e) 5.0.
(a) (b) (c) (d) (e)Finally, the comparison of calculation time between ICS-FEM and FEM at Intel (R) Xeon (R) CPU E3-1220 v3 @ 3.10 GHz, 16 G RAM is shown in Figure10, with the element number of 60, 240, 960, and 1500. As showed in Figure 10, the time required to solve algebraic equations by ICS-FEM is similar to that of FEM. Because the stiffness construction of ICS-FEM is based on smoothing cells inside each element, no coupling occurs between nodal degrees-of-freedom that are the distance of up to two elements. In other words, the bandwidth of ICS-FEM stiffness matrix is the same as that of FEM. Nevertheless, ICS-FEM is more effective in terms of generalized displacement (including displacement, electrical potential and magnetic potential) and computational efficiency (computation time for the same accuracy).Figure 10
Comparison of computational efficiency.


### 5.2. S-S Beam

The free vibrations on BaTiO3–CoFe2O4 FGMEE beams were studied by changing the gradient function form in Figure 4. The geometrical parameters were the same as the C-F beam. The simply supported boundary conditions were used: u=w=Ψ=Φ=0 at (x=0, z=h/2) and w=Ψ=Φ=0 at (x=L, z=h/2).Firstly, the free vibration frequencies for BaTiO3–CoFe2O4 FGMEE beams were calculated by both ICS-FEM with 60 × 4 meshes and FEM with 240 × 16 meshes (Figure 11). Results show the first eleven natural frequencies calculated by ICS-FEM are closer to the reference solutions than those calculated by FEM, indicating ICS-FEM is more efficient than FEM due to the reduced number of meshes. The differences in natural frequencies feq, feq_re, and feq_ΦΦ are marginal, so the magnetic effect only slightly impacts the natural frequencies of FGMEE beams. The natural frequencies increase with the rise of the exponential factor.Figure 11
Structural frequenciesf for simply supported BaTiO3-CoFe2O4 FGMEE beams with exponential factor η = (a) 0, (b) 0.5, (c) 1.0, (d) 2.0, and (e) 5.0.
(a) (b) (c) (d) (e)Secondly, as for BaTiO3–CoFe2O4 FGMEE beams with power law type (n=0, 1.0, 5.0, 10.0, 20.0), the bottom surface (z=-h/2) is BaTiO3–CoFe2O4 whereas the top surface (z=h/2) is CoFe2O4. The free vibration frequencies for FGMEE beams calculated by ICS-FEM with 60 × 4 meshes and FEM with 240 × 16 meshes are shown in Figure 12. Results show the first eleven natural frequencies calculated by ICS-FEM are closer to the reference solutions than those calculated by FEM. Also the differences in feq, feq_re, and feq_ΦΦ are marginal, so the magnetic effect does not largely impact the natural frequencies of FGMEE beams. Meanwhile, the natural frequency of BaTiO3–CoFe2O4 FGMEE beams is between those of BaTiO3–CoFe2O4 MEE beams and CoFe2O4 MEE beams.Figure 12
Structural frequenciesf of simply supported BaTiO3-CoFe2O4 FGMEE beams with power law index n = (a) 0, (b) 1.0, (c) 5.0, (d) 10.0, and (e) 20.0.
(a) (b) (c) (d) (e)


### 5.3. Typical MEMS-Based FGMEE Energy Harvester

The model of FGMEE energy harvester developed by ICS-FEM is shown in Figure13. The free vibrations on the FGMEE energy harvester were studied by changing the exponential factor. The geometrical parameters were L=30 mm, H=2 mm, and its structure was fabricated with BaTiO3 FGMEE.Figure 13
Typical MEMS-based energy harvester fabricated with BaTiO3 FGMEE. (a) Model of the energy harvester; (b) simplified model of ICS-FEM.
(a) (b)The free vibration frequencies for the FGMEE energy harvester calculated by ICS-FEM with 60 × 4 meshes and FEM with 240 × 16 meshes are shown in Figure14. The first eleven natural frequencies calculated by ICS-FEM are closer to the reference solutions than those calculated by FEM, indicating ICS-FEM is more efficient than FEM owing to the reduced number of meshes. The ICS-FEM does not take the derivative of the shape functions of the element and can be much easily generated automatically for complicated domains, thus significantly decreasing the numerical errors. The natural frequencies feq and feq_ΦΦ agree well with each other since the piezomagnetic phase is absent from the BaTiO3 FGMEE energy harvester. Moreover, feq_re coincides with feq of the BaTiO3 FGMEE energy harvester as the magnetoelectric effect is absent in pure BaTiO3 materials. The fuu and feq are very close, so the piezoelectric effect only slightly affects the natural frequencies of the energy harvester, which increase with the rise of the exponential factor.Figure 14
Structural frequenciesf for BaTiO3 FGMEE energy harvester with exponential factor η = (a) 0, (b) 0.5, (c) 1.0, (d) 2.0, and (e) 5.0.
(a) (b) (c) (d) (e)The Wilson-θ method and the equivalent stiffness matrix [Keq] were employed to solve the dynamic response of the FGMEE energy harvester. The parameters were set as time step = 0.005 s, θ=1.4; without damping; sine-wave transient load with a time period of 2 s; 4 cycles of loading (Figure 15). The dynamic behaviors for the harvester calculated by ICS-FEM with 60 × 4 meshes and FEM with 240 × 16 meshes are shown in Figures 16 and 17, respectively. The temporal variations of displacement uy and electric potential Φ calculated by ICS-FEM are closer to the reference solutions than those by FEM, which validate the accuracy of ICS-FEM. The temporal variations of uy and Φ of the FGMEE energy harvester decrease with the increase of the exponential factor when the material properties of the harvester are in exponential distribution.Figure 15
Sine-wave load at pointA of the BaTiO3 FGMEE energy harvester.Figure 16
Temporal variations of displacementuy at point B of BaTiO3 FGMEE energy harvester with different values of exponential factor calculated by ICS-FEM and FEM.Figure 17
Temporal variations of electric potentialΦ at point A of BaTiO3 FGMEE energy harvester with different values of exponential factor calculated by ICS-FEM and FEM.


## 5.1. C-F Beam

The free vibrations on FGMEE beams were calculated by changing the exponential factor (Figure4). The material properties of FGMEE beams were governed by the z-direction exponential variation. The following geometrical parameters were considered: length L = 0.3 m and width H = 0.02 m with the assumption of plane stress. Boundary conditions were u = w = Φ = Ψ = 0 at the clamped end. Table 2 gives the various structural frequencies in the study.Table 2
List offused in the study.
Structural frequency Matrix used to compute the structural frequency f u u K u u f e q K e q f e q _ r e K e q _ r e d u c e d f e q _ Φ Φ K e q _ Φ Φ f e q _ Ψ Ψ K e q _ Ψ ΨFigure 4
Geometry of an FGMEE beam and the coordinates.Firstly, the convergence of ICS-FEM was verified by using BaTiO3-CoFe2O4 FGMEE beams, with properties listed in Table 1. The natural frequencies of these beams were calculated using ICS-FEM with different meshes (30 × 2, 60 × 4, 120 × 8, 150 × 10) (Figure 5). The simulation results with different meshes agree well, which prove the good convergence of ICS-FEM. The first- to third-order modes of clamp-free BaTiO3–CoFe2O4 FGMEE beams using different elements with exponential factor η = 1.0 and 5.0 were calculated by ICS-FEM, and the results were summarized in Figures 6 and 7, respectively. It was found that the first- to third-order modes of BaTiO3–CoFe2O4 FGMEE beams in the same gradient distribution were basically not affected by equivalent stiffness matrix or mesh number. the first- to third-order modes of BaTiO3–CoFe2O4 FGMEE beams were basically not different between η = 1.0 and 5.0.Figure 5
Structural frequencies (a)fuu; (b) feq; (c) feq-re; (d) feq_ΦΦ; (e) feq_ΨΨ of clamp-free BaTiO3–CoFe2O4 FGMEE beams.
(a) (b) (c) (d) (e)Figure 6
The first- to third-order modes of clamp-free BaTiO3–CoFe2O4 FGMEE beams using different elements with exponential factor η = 1.0 by ICS-FEM.Figure 7
The first- to third-order modes of clamp-free BaTiO3–CoFe2O4 FGMEE beams using different elements with exponential factor η = 5.0 by ICS-FEM.Secondly, the free vibration frequencies of CoFe2O4 FGMEE beams were studied by both ICS-FEM and FEM using extremely irregular elements, with domain discretization shown in Figure 8. The frequencies of clamp-free CoFe2O4 FGMEE beams with different values of exponential factor are shown in Figure 9; the first eleven natural frequencies calculated by ICS-FEM are smaller than those calculated by FEM. The validity of ICS-FEM is verified by the agreements between the calculations and the reference solutions. The shape of quadrilateral element in FEM cannot be severely distorted but was eliminated in ICS-FEM. ICS-FEM abstains from calculating the derivative of the shape functions of an element, and the area integral of the solution domain is converted to the boundary integral. The stiffness of FGMEE structures is improved because ICS-FEM does not require continuity of the shape function. The ICS-FEM provides a continuous system with a close-to-exact stiffness, which could be automatically and more easily generated for complicated domains, thus significantly decreasing the numerical error. The free vibration of CoFe2O4 FGMEE beam, a pure CoFe2O4 material without piezoelectric or magnetoelectric material coefficients, influences structural frequency feq because the magnetic effect is marginally higher compared with fuu. fΦΦ coincides with fuu since piezoelectric phase is absent in CoFe2O4. Similarly, feq_re coincides with feq of the CoFe2O4 FGMEE beam as the magnetoelectric effect is absent in pure CoFe2O4 FGMEE beam. The natural frequencies of CoFe2O4 FGMEE beams increase with the rise of the exponential factor.Figure 8
Domain discretization using four-node extremely irregular elements.Figure 9
Structural frequenciesf of clamp-free CoFe2O4 FGMEE beams with exponential factor η = (a) 0, (b) 0.5, (c) 1.0, (d) 2.0, and (e) 5.0.
(a) (b) (c) (d) (e)Finally, the comparison of calculation time between ICS-FEM and FEM at Intel (R) Xeon (R) CPU E3-1220 v3 @ 3.10 GHz, 16 G RAM is shown in Figure10, with the element number of 60, 240, 960, and 1500. As showed in Figure 10, the time required to solve algebraic equations by ICS-FEM is similar to that of FEM. Because the stiffness construction of ICS-FEM is based on smoothing cells inside each element, no coupling occurs between nodal degrees-of-freedom that are the distance of up to two elements. In other words, the bandwidth of ICS-FEM stiffness matrix is the same as that of FEM. Nevertheless, ICS-FEM is more effective in terms of generalized displacement (including displacement, electrical potential and magnetic potential) and computational efficiency (computation time for the same accuracy).Figure 10
Comparison of computational efficiency.


## 5.2. S-S Beam

The free vibrations on BaTiO3–CoFe2O4 FGMEE beams were studied by changing the gradient function form in Figure 4. The geometrical parameters were the same as the C-F beam. The simply supported boundary conditions were used: u=w=Ψ=Φ=0 at (x=0, z=h/2) and w=Ψ=Φ=0 at (x=L, z=h/2).Firstly, the free vibration frequencies for BaTiO3–CoFe2O4 FGMEE beams were calculated by both ICS-FEM with 60 × 4 meshes and FEM with 240 × 16 meshes (Figure 11). Results show the first eleven natural frequencies calculated by ICS-FEM are closer to the reference solutions than those calculated by FEM, indicating ICS-FEM is more efficient than FEM due to the reduced number of meshes. The differences in natural frequencies feq, feq_re, and feq_ΦΦ are marginal, so the magnetic effect only slightly impacts the natural frequencies of FGMEE beams. The natural frequencies increase with the rise of the exponential factor.Figure 11
Structural frequenciesf for simply supported BaTiO3-CoFe2O4 FGMEE beams with exponential factor η = (a) 0, (b) 0.5, (c) 1.0, (d) 2.0, and (e) 5.0.
(a) (b) (c) (d) (e)Secondly, as for BaTiO3–CoFe2O4 FGMEE beams with power law type (n=0, 1.0, 5.0, 10.0, 20.0), the bottom surface (z=-h/2) is BaTiO3–CoFe2O4 whereas the top surface (z=h/2) is CoFe2O4. The free vibration frequencies for FGMEE beams calculated by ICS-FEM with 60 × 4 meshes and FEM with 240 × 16 meshes are shown in Figure 12. Results show the first eleven natural frequencies calculated by ICS-FEM are closer to the reference solutions than those calculated by FEM. Also the differences in feq, feq_re, and feq_ΦΦ are marginal, so the magnetic effect does not largely impact the natural frequencies of FGMEE beams. Meanwhile, the natural frequency of BaTiO3–CoFe2O4 FGMEE beams is between those of BaTiO3–CoFe2O4 MEE beams and CoFe2O4 MEE beams.Figure 12
Structural frequenciesf of simply supported BaTiO3-CoFe2O4 FGMEE beams with power law index n = (a) 0, (b) 1.0, (c) 5.0, (d) 10.0, and (e) 20.0.
(a) (b) (c) (d) (e)


## 5.3. Typical MEMS-Based FGMEE Energy Harvester

The model of FGMEE energy harvester developed by ICS-FEM is shown in Figure13. The free vibrations on the FGMEE energy harvester were studied by changing the exponential factor. The geometrical parameters were L=30 mm, H=2 mm, and its structure was fabricated with BaTiO3 FGMEE.Figure 13
Typical MEMS-based energy harvester fabricated with BaTiO3 FGMEE. (a) Model of the energy harvester; (b) simplified model of ICS-FEM.
(a) (b)The free vibration frequencies for the FGMEE energy harvester calculated by ICS-FEM with 60 × 4 meshes and FEM with 240 × 16 meshes are shown in Figure14. The first eleven natural frequencies calculated by ICS-FEM are closer to the reference solutions than those calculated by FEM, indicating ICS-FEM is more efficient than FEM owing to the reduced number of meshes. The ICS-FEM does not take the derivative of the shape functions of the element and can be much easily generated automatically for complicated domains, thus significantly decreasing the numerical errors. The natural frequencies feq and feq_ΦΦ agree well with each other since the piezomagnetic phase is absent from the BaTiO3 FGMEE energy harvester. Moreover, feq_re coincides with feq of the BaTiO3 FGMEE energy harvester as the magnetoelectric effect is absent in pure BaTiO3 materials. The fuu and feq are very close, so the piezoelectric effect only slightly affects the natural frequencies of the energy harvester, which increase with the rise of the exponential factor.Figure 14
Structural frequenciesf for BaTiO3 FGMEE energy harvester with exponential factor η = (a) 0, (b) 0.5, (c) 1.0, (d) 2.0, and (e) 5.0.
(a) (b) (c) (d) (e)The Wilson-θ method and the equivalent stiffness matrix [Keq] were employed to solve the dynamic response of the FGMEE energy harvester. The parameters were set as time step = 0.005 s, θ=1.4; without damping; sine-wave transient load with a time period of 2 s; 4 cycles of loading (Figure 15). The dynamic behaviors for the harvester calculated by ICS-FEM with 60 × 4 meshes and FEM with 240 × 16 meshes are shown in Figures 16 and 17, respectively. The temporal variations of displacement uy and electric potential Φ calculated by ICS-FEM are closer to the reference solutions than those by FEM, which validate the accuracy of ICS-FEM. The temporal variations of uy and Φ of the FGMEE energy harvester decrease with the increase of the exponential factor when the material properties of the harvester are in exponential distribution.Figure 15
Sine-wave load at pointA of the BaTiO3 FGMEE energy harvester.Figure 16
Temporal variations of displacementuy at point B of BaTiO3 FGMEE energy harvester with different values of exponential factor calculated by ICS-FEM and FEM.Figure 17
Temporal variations of electric potentialΦ at point A of BaTiO3 FGMEE energy harvester with different values of exponential factor calculated by ICS-FEM and FEM.


## 6. Conclusions

The free vibrations on FGMEE structures were studied. Firstly, ICS-FEM for FGMEE materials was formulated by incorporating gradient smoothing into the FEM-based computation for the FGMEE multi-physics field. Then the equations of free vibration computation were deduced for the multi-physics coupling field of FGMEE materials. Finally, the FGMEE beams were calculated with functional gradient exponential form or power law form under different boundary conditions.(i) ICS-FEM reduced the systematic stiffness of the finite element, which improved the computational accuracy compared with FEM under the same element number. ICS-FEM was more efficient than FEM in terms of computation time for the same accuracy.(ii) Due to the material property changes in each smoothing element, the true parameters at the Gaussian integration point were adopted directly. ICS-FEM avoided the derivative of the shape functions but only transformed the area integral to the boundary integral in the solution domain, which omitted the requirement of continuity of the shape function.(iii) The magnetic effect slightly influenced the natural frequencies of FGMEE beams, which increased with the increase of the exponential factor when the material properties of FGMEE beams were under exponential distribution. The natural frequency of FGMEE beams lied in between those of the FGMEE beams using the bottom surface of materials and the FGMEE beams using the upper surface, when the material properties of FGMEE beams were under the power law distribution.(iv) The natural frequencies and general displacements of the FGMEE energy harvester developed by ICS-FEM were more accurate compared with FEM, owing to the reduced number of meshes.



---

*Source: 5141060-2018-03-07.xml*