DOI: 10.11817/j.issn.1672-7207.2016.01.008

基于非均匀光滑有限元法的功能梯度压电梁自由振动分析

蔡斌¹，周立明²

(1. 吉林建筑大学 土木工程学院, 吉林 长春，130118；

2. 吉林大学 机械科学与工程学院, 吉林 长春，130022)

摘要：为了提高求解功能梯度压电材料动响应的精度，克服有限元系统刚度偏硬的缺点，提出非均匀Cell-based 光滑有限元法。基于单元的梯度光滑操作，考虑材料物性沿宽度方向呈梯度连续变化，推导非均匀Cell-based光滑有限元法的基本公式，分析功能梯度压电悬臂梁的材料物性参数遵循不同梯度分布规律时结构自由振动的固有频率与振型，并与FEM求解结果进行对比。研究结果表明：光滑梯度操作可降低有限元系统的刚度，非均匀Cell-based 光滑有限元法的数值解更加接近真实解，从而可为功能梯度压电材料的进一步应用提供参考。

关键词：非均匀光滑有限元法；功能梯度压电材料；自由振动；梯度光滑技术

中图分类号：TB115       文献标志码：A         文章编号：1672-7207(2016)01-0048-06

Inhomogeneous smoothed finite element method for free vibration analysis of functionally gradient piezoelectric cantilevers

CAI Bin¹, ZHOU Liming²

(1. School of Civil Engineering, Jilin Jianzhu University, Changchun 130118, China;

2. School of Mechanical Science and Engineering, Jilin University, Changchun 130022, China)

Abstract: In order to improve the solving precision of dynamic response problems of functionally gradient piezoelectric materials and soften the “over-stiffness”, an inhomogeneous Cell-based smoothed finite element method was presented. Based on gradient smoothing technique of elements and considering the change of gradient of material properties only in the breadth direction, the basic formula of inhomogeneous Cell-based smoothed finite element (ICS-FEM) was derived. Using this proposed method, free vibration characteristics of the functionally graded piezoelectric cantilever were analyzed when the material parameters followed different gradient distributions. The accuracy of the proposed method was confirmed by comparing the FEM results. The results show that the gradient smoothing technique can reduce the stiffness of the finite element system and improve the accuracy of the solution for electromechanical coupling systems by using ICS-FEM, which provides reference for further application of Cell-based smoothed finite element method.

Key words: inhomogeneous smoothed finite element method; functionally graded piezoelectric materials; free vibration; gradient smoothing technique

为了克服有限元系统刚度偏大的缺点，LIU等^[1-2]提出了一些基于梯度光滑的有限元算法。NGUYEN- VAN等^[3]将这种基于梯度光滑的有限元算法拓展到压电领域，随后，NGUYEN-XUAN^[4]采用Edge-Based光滑有限元法对二维压电结构力学问题进行了研究。通过数值分析，发现光滑有限元能够提高求解精度，具有广阔的应用前景。工程中目前流行的压电元件多为多层结构，但元件中某些材料组分和物性的突然变化往往会导致器件在层间界面处存在明显的局部应力失配现象，导致黏结层在高温易蠕变，在低温易开裂，大大缩短了元件的寿命。为解决这类问题，ZHU等^[5]将功能梯度的概念引入到压电智能材料中，制备了功能梯度压电材料，其兼具了压电和梯度二者的优点。此后，一些研究者^[6-7]在功能梯度压电材料这一领域开展研究。基于Euler-Bernoulli梁理论，FU^[8]求解了功能梯度压电梁热-电屈曲的精确解；CHEN等^[9]基于Euler梁理论，获得了功能梯度压电梁的自由振动频率；ZHAO等^[10]采用无网格法求解了不同载荷作用下功梯度压电板的静态弯曲问题；KOMIJANI等^[11]研究了功能梯度压电材料执行器非线性热-机-电多场耦合响应问题；DAI等^[12]分析了功能梯度压电材料的反平面裂纹问题；JODAEI等^[13]采用微分求积法求解了功梯度压电板在不同边界条件下静态分析的三维弹性解。由于材料的非均匀性和多场耦合特性，使得求解功能梯度压电元材料压电方程^[14-15]的难度大大增加，从而导致研究工作大多局限于静力问题，对动力学问题的研究仍然很少。本文作者针对功能梯度压电悬臂梁，考虑材料特性沿宽度方向呈梯度连续变化，基于单元的梯度光滑操作，推导非均匀Cell-Based光滑有限元(ICS-FEM)的基本公式，采用ICS-FEM分析功能梯度压电悬臂梁自由振动固有频率与振型，并与FEM的计算结果进行对比。

1  基本方程

功能梯度压电材料问题的基本方程如下：

                (1)

             (2)

                       (3)

                   (4)

式中：；i，j和k分别为x，y和

z方向的单位向量；σ为应力张量；f为体力向量；ε为应变张量；u为位移向量；为u对时间t的二次导数；D为电位移向量；ρ为质量密度；p^e为自由电荷密度；E为电场强度；f为电势。在边界上满足如下条件：

σ·n=，在力边界Γ_t上             (5)

u=，在位移边界Γ_u上          (6)

D·n=-，在电位移边界Γ_D上   (7)

f=，在电势边界Γ_f上           (8)

式中：和分别为表面张力和表面电荷密度；n为介质的外法向矢量；和分别为给定位移和电势。功能梯度压电材料的本构方程为

σ=C_Eε-e^TE                 (9)

D=eε+gE                  (10)

式中：C_E为弹性模量张量；e为压电常数张量；g为介电常数张量。对于横观各向同性功能梯度压电材料，x-z平面为各向同性面，研究x-z平面内力电耦合问题，满足

ε_yy=ε_zy=ε_xy=E_y=0            (11)

则式(9)和式(10)的矩阵形式为

  (12)

   (13)

材料物性参数遵循如下梯度分布规律：

M(z)=M_bf_b(z)                (14)

式中：M(z)为c_ij，e_ij，g_ii和ρ等在结构中的实际物理参数；M_b为功能梯度压电材料结构底部的物理参数；f_b(z)为底部材料的体积分数分布函数，见图1。f_b(z)一般有3种函数形式：

f_b(z)=e^α(z/h)                (15)

f_b(z)=(1+z/h)^b                 (16)

f_b(z)=[1+gsin(pz/(2h))]                 (17)

式中：α，b和g分别为指数函数、幂函数和正弦函数分布规律的形状因子；z为沿z轴坐标；h为梁宽。

图1  功能梯度压电材料物参梯度分布图

Fig. 1  Gradient distribution of functionally graded piezoelectric materials

2  非均匀Cell-based光滑有限元法

将求解域离散成np个单元，包含N_n个节点，每个单元内的广义位移u和广义电势表示为

            (18)

            (19)

式中：N_u和N_f分别为ICS-FEM位移形函数和电势形函数；q和j分别为节点位移向量和节点电势向量。

图2所示为四节点单元划分为4个光滑子元的场

节点、边中间光滑节点、中心光滑节点、边高斯点、外法向向量分布情况及形函数值。

在光滑子元内任意点x^k处，光滑应变和光滑电场强度分别为

          (20)

         (21)

式中：和分别为FEM中的应变和电场强度；Φ(x-x^k)为光滑函数，

          (22)

                (23)

将式(22)代入式(20)和(21)，得

            (24)

             (25)

式中：为光滑域的边界；和为光滑域边界的外法线向量矩阵，

，             (26)

图2  光滑子元及形函数值

Fig. 2  Smoothing subcells and values of shape functions

将式(24)和式(25)改写为

              (27)

              (28)

式中：n_e为光滑子元个数；

      (29)

         (30)

N_i为相关节点i的形函数。

在高斯点处，式(29)和式(30)可改写为

    (31)

        (32)

式中：和分别为光滑边界的高斯点和长度，n_b为每个光滑节点域的边界总数。

ICS-FEM与FEM的本质差异为：FEM需要对单元形函数矩阵求导，则ICS-FEM只需采用光滑元边界高斯点处的形函数，不涉及形函数求导，降低了形函数的连续性要求，从而提高了计算方法的精度和收敛性。功梯度压电耦合系统的动力学模型可由Hamilton原理导出，形式如下：

    (33)

式中：

               (34)

   (35)

m_i=ρ_it(i=1, 2, 3, 4)为对应节点i的第i个光滑元的质量；t为光滑元厚度；ρ_i为第i个光滑子元的质心密度；diag表示该矩阵的对角线方阵。

               (36)

                (37)

            (38)

          (39)

              (40)

式中：n_c=n_p×n_e。采用非均匀光滑单元计算单元刚度矩阵，在单元s中4个光滑子元(i=1, 2, 3, 4)的物理参数不相同，直接取高斯积分点处的实际物理参数进行计算，这样可以使每个单元内部也能够体现材料属性的变化。

由式(33)可知位移和电势是耦合的。在分析特征值时，F=Q=0，先将和电势相关的自由度进行凝聚，

            (41)

将式(41)代入式(33)，有

               (42)

式中：

                (43)

        (44)

3  数值算例

3.1  算例1

如图3所示，一功能梯度压电悬臂梁，长度L=20 cm，宽度h=5 cm，厚度b=1 cm。假设材料成分沿厚度方向呈梯度变化，考虑3种梯度分布形式：(Ⅰ) 底部材料为PZT-4，弹性系数、密度、压电系数和介电系数均按指数函数变化，=0.5；(Ⅱ) 底部材料为PZT-4，弹性系数、压电系数和介电系数均按指数函数模式变化，=0.5，密度按幂函数变化，=0.5；(Ⅲ) 底部材料为PZT-4，弹性系数、压电系数和介电系数均按指数函数变化，=0.5，密度无变化。光滑子元个数取为4，材料参数如表1所示，求解其固有频率。

图3  悬臂梁模型

Fig. 3  Geometry of cantilever beam

表1  材料常数

Table 1  Material constants

表2所示为功能梯度压电悬臂梁在I，II和III这3种材料梯度变化模式下，采用40×10均匀分布单元，由ICS-FEM和FEM求解得到的前5阶固有频率。从表2可以看出：ICS-FEM计算得到的固有频低于FEM求解的固有频率，两者的相对误差最大为2.85%，从而验证了ICS-FEM的正确性。CS-FEM不需对单元形函数矩阵求导，将求解域内面积分改变为边界积分，形函数的连续性要求低。

表2  不同材料梯度分布形式下的前5阶固有频率

Table 2  First five natural frequency values of materials with different gradient distributions

3.2  算例2

功能梯度压电悬臂梁，长度L=40 cm，其他参数不变，底部材料为PZT-4，顶部材料为C-91，弹性系数、密度、压电系数和介电系数均按指数函数变化。材料参数如表1所示，求解其固有频率和振型。

图4所示为功能梯度压电悬臂梁采用80×10均匀分布单元，由ICS-FEM和FEM求解得到的前5阶固有振型。从图4可看出：两者振型相吻合。这进一步验证了ICS-FEM的正确性和有效性。

图4  功能梯度悬臂梁前5阶振型

Fig. 4  First five order modes of vibration of functionally gradient piezoelectric cantilevers

表3所示为功能梯度压电悬臂梁采用80×10均匀分布单元，由ICS-FEM和FEM求解得到的前5阶固有频率，以Ansys采用320×40单元所得结果作为参考解。从表3可以看出：在同样单元条件下，ICS-FEM计算得到的固有频率低于FEM求解的固有频率，其数值解更加接近真实解，从而验证了ICS-FEM可降低有限元系统刚度，提高求解精度，有效地求解功能梯度压电结构的动响应问题。

表3  功能梯度压电悬臂梁的前5阶固有频率

Table 3  First five natural frequency of functionally gradient piezoelectric cantilevers              kHz

4  结论

1) ICS-FEM求解功能梯度压电结构的动响应问题是正确的、有效的。

2) ICS-FEM可降低有限元系统刚度，在相同单元条件下，其计算精度高于FEM模拟结果的精度，更接近真实值。

参考文献：

[1] LIU G R, DAI K Y, NGUYEN T T. A smoothed finite element method for mechanics problems[J]. Comput Mech, 2007, 39(6): 859-877.

[2] LI ERIC, HE Z C, CHEN L, et al. An ultra-accurate hybrid smoothed finite element method for piezoelectric problem[J]. Engineering Analysis with Boundary Elements, 2015, 50: 188-197.

[3] NGUYEN-VAN H, MAI-DUY N, TRAN-CONG T. A node- based element for analysis of planar piezoelectric structures[J]. CMES: Computer Modeling in Engineering and Sciences, 2008, 36(1): 65-95.

[4] NGUYEN-XUAN H, LIU G R, NGUYEN-THOI T, et al. An edge-based smoothed finite element method for analysis of two-dimensional piezoelectric structures[J]. Smart Materials and Structures, 2009, 18(6): 065015.

[5] ZHU Xinhua, MENG Zhongyan. Operational principle, fabrication and displacement characteristics of a functionally gradient piezoelectric ceramic actuator[J]. Sensors and Actuators A: Physical, 1995, 48(3): 169-176.

[6] KOMIJANI M, REDDY J N, ESLAMI M R. Nonlinear analysis of microstructure-dependent functionally graded piezoelectric material actuators[J]. Journal of the Mechanics and Physics of Solids, 2014, 63: 214-227.

[7] ZHOU Liming, MENG Guangwei, LI Feng, et al. Cell-based smoothed finite element method-virtual crack closure technique for a piezoelectric material of crack[J]. Mathematical Problems in Engineering, 2015: 371083.

[8] FU Yiming, WANG Jianzhe, MAO Yiqi. Nonlinear analysis of buckling, free vibration and dynamic stability for the piezoelectric functionally graded beams in thermal environment[J]. Applied Mathematical Modelling, 2012, 36(9): 4324-4340.

[9] CHEN Haiyan, JIN Dengren, MENG Zhongyan. Dynamic characteristics of functionally gradient piezoelectric actuators [C]// Properties and Applications of Dielectric Materials. Xi’an, China: IEEE, 2000: 983-987.

[10] ZHAO X, LEE Y Y, LIEW K M. Free vibration analysis of functionally graded plates using the element-free method[J]. Journal of Sound and Vibration, 2009, 319(3): 918-939.

[11] KOMIJANI M, REDDY J N, ESLAMI M R. Nonlinear analysis of microstructure-dependent functionally graded piezoelectric material actuators[J]. Journal of the Mechanics and Physics of Solids, 2014, 63(1): 214-227.

[12] DAI Yao, CHONG Xiao, LI Shimin. The higher order crack tip fields for anti-plane crack in functionally graded piezoelectric materials[J]. Applied Mechanics and Materials, 2014, 472: 617-620.

[13] JODAEI A. 3D elasticity solution for static analysis of functionally graded piezoelectric annular plates on elastic foundations using SSDQM[J]. Mechanica, 2014, 49(1): 215-237.

[14] LEZGY-NAZARGAH M, VIDAL P, POLIT O. An efficient finite element model for static and dynamic analyses of functionally graded piezoelectric beams[J]. Composite Structures, 2013, 104: 71-84.

[15] LI Y S, FENG W J, CAI Z Y. Bending and free vibration of functionally graded piezoelectric beam based on modified strain gradient theory[J]. Composite Structures, 2014, 115(1): 41-50.

(编辑  陈灿华)

收稿日期：2015-03-10；修回日期：2015-05-08

基金项目(Foundation item)：国家自然科学基金资助项目(51305157)；吉林省科技厅基金资助项目(20130305006GX) (Project(51305157) supported by the National Natural Science Foundation of China; Project(20130305006GX) supported by Department of Science and Technology Fund of Jilin Province)

通信作者：周立明，博士，讲师，从事计算固体力学研究；E-mail: lmzhou@jlu.edu.cn