华科大压电材料

1、一、 Y 切石英板的静态分析1.1问题描述计算对象为一块石英板，其尺寸为：长度方向（ x 轴）：L=0.5m厚度方向（ y 轴）：H=2h=0.01m宽度方向（ z 轴）：W=0.05m板厚度方向两个平面有电极，大小为V/2=1000v。用 ANSYS 计算自由和夹支边界两种情况下的电极所在平面单位面积上的电荷 e。该石英板的材料常数为：密度 =2650kg/m3弹性常数矩阵：压电常数矩阵：介电常数矩阵：精选文库1.2理论解边值问题为：电极所在平面单位面积上电荷的理论解为：（ 1）自由表面（ 2）夹支表面1.3 ANSYS 求解(1) 单元选择根据该石英板的材料属性，可以选择 solid226

2、 单元进行求解。 Solid226 单元为耦合场 20 节点六面体单元， KEYOPT(1)=1001 即可激活压电自由度、位移和电压。输入“ /prep7”进入前处理。APDL ： et,1,solid226,1001ET 命令的功能为从单元库中指定一个单元类型。其使用格式为：ET, ITYPE, Ename, KOP1, KOP2, KOP3, KOP4, KOP5, KOP6, INOPR(2) 介电系数矩阵由于该石英板介电常数为不是各项同性的，所以需要用TB ， DPER 和TBDATA 命令指定各向异性介电系数矩阵中的各项， ANSYS 中需要输入的是相对介电常数，因此需要用介电常数

3、除以真空中介电常数0- 8.854 10- 12 。=TBOPT=0 则指定恒应变介电系数矩阵S ； TBOPT=1 则指定恒应力介电系数矩TS阵 ，后者将被转换为恒应变介电系数矩阵 。值得指出的是， solid5 单元只-2精选文库能用 MP 命令制定 PERX， PERY 和 PERZ 三项，所以本次计算不采用。ANSYS 中输入相对介电系数矩阵的各项时也需要注意顺序，其顺序如下：从而该石英的介电系数矩阵可以按如下输入：tb,DPER,1tbdata,1,4.4285,4.4974,4.5642tbdata,5,0.0971TB 和 TBDATA 指令用于输入矩阵中的各项。 其中 DPER

4、 用于指定该矩阵为各向异性相对介电系数矩阵，另外 ANEL 为各项异性弹性系数矩阵， PIEZ 为压电系数矩阵。与 ANSYS 中介电系数矩阵的输入类似，弹性系数矩阵和压电系数矩阵的输入也要遵从 ANSYS 中的规则。(3) 弹性系数矩阵弹性系数矩阵主要包括刚度系数矩阵c 和柔度系数矩阵s 。用户可以用TBOPT 在这两种矩阵中选择。 ANSYS 中刚度系数矩阵 c 中常数的顺序为：从而该石英的弹性系数矩阵可以按如下输入：tb,ANEL,1tbdata,1,86.74E9,-8.25E9,27.15E9,-3.66E9! c11,c21,c31,c61,c41tbdata,7,129.77E9

5、,-7.42E9,5.7E9! c22,c32,c62,c42tbdata,12,102.83E9,9.92E9! c33,c63,c43tbdata,16,29.01E9,2.53E9! c66,c46,c56tbdata,19,38.61E9! c44tbdata,21,68.81E9! c55(4) 压电系数矩阵压电矩阵可以定义压电应力矩阵e 或压电应变矩阵 d 。压电应力矩阵 e 与-3精选文库各项异性弹性刚度矩阵 c 相联系，而压电应变矩阵 d 与弹性柔度矩阵 s相联系。ANSYS 中压电应力矩阵 e 中常数的顺序为：因此命令流可以输入如下：tb,PIEZ,1tbdata,1,0.1

6、71! e11tbdata,4,-0.152! e21tbdata,7,-0.0187! e31tbdata,11,-0.095,0.067! e62,e63tbdata,13,0.067! e41tbdata,17,0.108,-0.0761 ! e52,e53(5) 网格划分由于本次求解使用的 solid226 为高阶单元，因此厚度方向划分两层网格即可，即单元长度可选为 H 的一半 0.005m：LESIZE,ALL,5e-3, , , ,1, , ,1Vmesh,1LESIZE 命令用来对所选择的线设置网格单元大小，格式为：LESIZE, NL1, SIZE, ANGSIZ, NDIV,

7、 SPACE, KFORC, LAYER1, LAYER2,KYNDIV划分后网格如下图：-4精选文库(6) 静态求解输入“ /solu”进入求解部分。通过“ ANTYPE ”命令指定一种分析类型和重启动状态。其格式为：ANTYPE, Antype, Status, LDSTEP, SUBSTEP, Action 本例是静态求解，故该命令可为：/soluantype,static!选取静态求解nsel,s,loc,y,0!定义下表面电压d,all,volt,1000!d,all,ux,0!约束下表面，若自由边界则不需要nsel,s,loc,y,H1d,all,volt,-1000!定义上表面电

8、压!d,all,ux,0!约束上表面，若自由边界则不需要nsel,allsolvefini通过 NSEL 命令来选择一组节点子集。使用格式为：NSEL, Type, Item, Comp, VMIN, VMAX, VINC, KABS通过 D 命令在节点上施加DOF 约束。其格式为：D, NODE, Lab, VALUE, VALUE2, NEND, NINC, Lab2, Lab3, Lab4, Lab5, Lab6(7) 通用后处理-5精选文库输入“ /post1”进入通用后处理部分。本例采用的命令流为：/post1pldisp,2!显示变形前的轮廓和变形后的图形plnsol,d,y! 显

9、示厚度方向的电流密度save,YshearStaticClamped,dbfiniPLDISP 命令用来显示结构的位移，格式为： PLDISP, KUND 。PLNSOL 命令用来显示节点解。其格式为：PLNSOL, Item, Comp, KUND, Fact, FileID后处理结果如下：变形图（位移矢量和）：XZ 平面（厚度方向）显然产生了剪切。厚度方向电位移图：-6精选文库1.4结果对比自由边界，厚度方向电位移为：12VD 2k2622 2h其中22009527 81103，从而k26e2639 82.29 0110922c6610 12.D212V17.8110339.8210122

10、0001062k26 222h- 8.026C / m0.01夹持边界，厚度方向上的电位移为：D222V2h由于 k2627.8110 31 ，所以夹持情况和自由情况下的D2 几乎一样大。从 ANSYS 计算出的电位移图中可知， 即使在几何边界处有类似应力集中的效应，但是整体上石英片厚度方向的电位仪是一致的，大小为 - 0.80510- 5 C/m2，与理论值十分接近。-7精选文库二、 Y 切石英板的自由振动2.1问题描述计算对象为一块石英板，其尺寸为：长度方向（ x 轴）：L=0.1m厚度方向（ y 轴）：H=2h=0.01m宽度方向（ z 轴）：W=0.02m用 ANSYS 计算石英板在厚

11、度方向振动的模态，并输出前几阶频率。该石英板的材料常数为：密度 =2650kg/m3弹性常数矩阵：压电常数矩阵：介电常数矩阵：-8精选文库2.2理论解控制方程和边界条件为：解出圆频率为：nc66 nc66,n01 232H,.h因此频率为：fn c66 ,n 0,1,2,3.22H其中，c66 c66 1 k262 ， k262e26222c662.3 ANSYS 求解(1) 建模与常数输入与 Y 切石英板静态求解部分一样，这里不再赘述。值得注意的是，由于 Y 切石英板自由振动为厚度方向的剪切振动，因此厚度方向需要划分更密集的网格。在静态求解时仅划分了两层网格，因此求解厚度方向自由振动时不能再

12、采用静态求解时的网格。本例中对石英板的几何尺寸和网格大小都做出了调整，使得求解所得模态可视性很强，结果也更为准确。-9精选文库(2) 网格划分本例中石英板的几何尺寸为：长度方向（ x 轴）：L=0.1m厚度方向（ y 轴）：H=2h=0.01m宽度方向（ z 轴）：W=0.02m将厚度方向划分为 10 层，则单元长度为 0.001m。由于长度方向和宽度方向比较次要，因此在这两个方向上单元长度可以更大，从而总节点数不至于太大使得求解困难。另外单元也不可过于狭长，因此长度方向和宽度方向的单元长度可以为厚度方向的两倍，即 0.002m。APDL 如下：lsel,s,length,H1 !选出长度为

13、H1 的线段 elesize=1e-3LESIZE,ALL,elesizelsel,inve ! 选出剩下的线段 LESIZE,ALL,2*elesizeVmesh,1 !对实体 1经行网格划分-10精选文库(3) 模态求解对于模态分析，推荐使用 Block Lanczos 法（缺省）求解。指令为：/soluantype,modal!选取模态分析modopt,lanb,5 ! 定义模态分析选项，展开 5 阶模态 mxpand,5(4) 约束与耦合考虑到建模为三维模型，而实际上关注的是 Y 方向的厚度剪切振动，因此为了使求解效率更高，不产生过多的不需要的模态，这里需要加一些约束与耦合。在厚度方向

14、上的剪切振动中， 每一个 XZ 平面运动都一致， 且只能沿着长度方向（ X 方向）运动，因此约束及耦合如下：nlayer=H1/elesize! 厚度方向的层数h=0!初始平面所在的高度*do,i,1,nlayer! 进入循环， i 从 1 到最后一层nsel,s,loc,y,h! 选出高度为 h 的所有节点d,all,uy,0! 约束住 y 方向上位移d,all,uz,0! 约束住 z 方向上位移cp,i,ux,all!耦合住 x 方向上位移h=h+elesize! 更新当前平面为下一平面*enddo!结束循环某一层加上约束及耦合后的情况如下：-11精选文库耦合列表如下：耦合共 10 组，即

15、厚度方向10 层全部分别耦合住。(5) 通用后处理后处理我们关心求出的各阶频率及对应的振型，命令流如下：/post1set,list! 列表查看所求模态set,first !查看第一阶模态 plnsol,u,sum,2 !anmode,10,0.1 !动画查看振型set,next! 查看下一阶模态-12精选文库pldisp,2!anmode,10,0.1列表查看求出的频率：第一阶为刚体位移：第二阶开始的振型分别为：-13精选文库N=1 阶频率为 166302HzN=2 阶频率为 332245Hz-14精选文库N=3 阶频率为 498518HzN=4 阶频率为 665262Hz2.4结果对比(8

16、) 各阶振型图对比ANSYS 计算出的各阶振型图如上节所示。理论计算出的振型图如下图所示：-15精选文库二者计算结果一致。(9) 各阶频率对比理论计算各阶频率由下式决定：fnc66 ,n0,1,2,3.22H其中，c66c66 1k262， k262e26222c66带入数值得：2e26239 820.095 21097 8110 3k2622c6610 1229 01.c66c66 1 k26229.011091 7.81 10- 329.24 109nc66n29 24109f.,n0,1,2,3.22H20 012650.将 ANSYS 计算的值与理论解汇入下表：阶数 n理论解（ Hz

17、）ANSYS 解（ Hz ）误差11660771663020.135%23321543322450.027%34982314985180.058%46643096652620.143%可以看出 ANSYS计算出的结果具有十分高的精度。-16精选文库三、（大作业）压电复合梁的静态分析3.1问题描述如下图所示，一上下对称的压电复合梁。梁厚度方向极化，上下各有一个电极板，板上加有电压。该复合梁长 l=0.2m，中间层的弹性梁厚度为 he=0.02m，上下两个电极板厚度都为 hp=0.005m。本例为一平面应力问题。用 ANSYS 对其经行静态分析， 并与论文 ”Analysis of Beams w

18、ith Piezoelectric Actuators”中的结果进行对比。中间层的弹性梁E=2.1E11Pa，=0.3。压电板的物理参数如下：极化方向为 z 方向极化，但是考虑到ANSYS 建二维模型时为XY 平面内，因此可将本例中的XZ 平面在 ANSYS 中看做 XY 平面，即 y 方向极化。因此各材料常数矩阵如下：弹性常数矩阵为：81.30.4320.3290cpqE0.43264.50.4320109 N/ m20.3290.43281.3000025.6压电常数矩阵为：00002eip52 1508 00 C/ m.介电常数矩阵为：Si17精选文库3.2论文

19、结果论文中计算了两个例子，第一个例子为上下极板电压都为 V 1=V2=100v，悬臂梁、简支梁对应的中线（ z=0）处的扰度曲线为：第二个例子为上极板 V 1=100v，下极板 V 2 =200v，悬臂梁、简支梁对应的中线处的扰度曲线为：厚度方向应力、应变的变化曲线：中线（ z=0）处的位移：-18精选文库3.3 ANSYS 求解(1) 建模本例为平面应力问题，可以选择 plane13 单元。通过选择 KEYOPT(1)=7来激活 plane13单元的 volt 自由度。用 TB、TBDATA 来输入弹性常数矩阵、 压电常数矩阵和相对介电常数矩阵时都需要注意 ANSYS 中的顺序规定。 详情可

20、见第一部分， Y 切石英的静态求解，这里不再赘述。用于指定单元类型和参数输入的命令流为：! 压电单元et,1,plane13,7,2 !"2" 表示平面应力c11=81.3E9c12=0.329E9c13=0.432E9c33=64.5E9c44=25.6E9e31=-5.2e33=15.08e15=0g11=1300g33=1475mp,perx,1,g11 !相对介电常数mp,pery,1,g33tb,ANEL,1 ! 刚度系数矩阵tbdata,1,c11,c13,c12tbdata,7,c33,c13tbdata,12,c11tbdata,16,c44tb,PIEZ,

21、1 !压电（应力）系数矩阵-19精选文库tbdata,2,e31tbdata,5,e33tbdata,10,e15! 弹性单元et,2,plane182mp,EX,2,2.1E11mp,NUXY ,2,0.3本例的几何模型十分简单， 可以用 RECTNG 命令直接生成。 RECTNG 命令的格式为：RECTNG, X1, X2, Y1, Y2因此对应的命令流为：! 几何尺寸L=0.2He=0.02Hp=0.005! 创建几何体rectng,0,L,-0.5*He-Hp,-0.5*Herectng,0,L,-0.5*He,0.5*Herectng,0,L,0.5*He,0.5*He+Hpaglu

22、e,all其中 aglue 是将生成的三个矩形变成一个整体，从而不会出现多余的面、线、点。(2) 网格划分首先可以用命令 ESIZE 来制定单元大小。然后要分别对三个矩形进行网格划分。其中下极板对应的面的编号为 1，材料编号为 1；中间的弹性梁对应的面的编号为 4，材料编号为 2；上极板对应的面的编号为 5，材料编号为 1。因此命令流如下：! 进行网格划分 esize,0.5*Hp !指定单元大小mat,1$type,1amesh,1 !下极板mat,2$type,2amesh,4 !中间层mat,1$type,1amesh,5 !上极板(3) 静态求解静态求解的步骤比较简单。对于本例而言，要

23、注意边界条件的施加。对-20精选文库于悬臂梁，需要约束x=0 处的所有结构自由度，命令流如下：! 悬臂梁nsel,s,loc,x,0d,all,ux,0d,all,uy,0对于简支梁，需要约束中线（ z=0）两端 Y 方向的位移和左侧端点 X 方向的位移，命令流如下：! 简支梁nsel,s,loc,x,0nsel,r,loc,y,0d,all,ux,0d,all,uy,0nsel,s,loc,x,Lnsel,r,loc,y,0d,all,uy,0对于上下极板的电压的施加类似。最后还要注意极板与弹性梁接触的两个面是接地的，因此电压为 0，命令流如下：! 定义上极板电压nsel,s,loc,y,0

24、.5*He+Hpd,all,volt,V !V=100! 中间层上下面接地nsel,s,loc,y,0.5*Hed,all,volt,0nsel,s,loc,y,-0.5*Hed,all,volt,0! 定义下极板电压nsel,s,loc,y,-0.5*He-Hp! 注意选择算例 case1和 case2 !d,all,volt,V !case1 d,all,volt,2*V !case2(4) 通用后处理简支梁和悬臂梁的位移云图 (以 case2为例 )：-21精选文库为了得出与论文结果相应的几副图，需要将某些路径上的值以图的形式输出。通过命令 PATH， PPATH， PDEF 和 PLPATH 可以实现这一点。对于算例 1，中线（ z=0）上的扰度需要输出为图的形式。其命令流如下：!取出路径 z=0path,z0,2ppath,1,0,0ppath,2,L,0!显示该路径上挠度wpdef,w,u,y-22精选文库/axlab,x, x(m) !x 轴名称/axlab,y, w(m) !y 轴名称plpath,w相应的，悬臂梁和简支梁的图如下：上图： case1，悬臂梁，极值 -9.751E-7，论文结果接近但小于 -1E-6上图： case1，简支梁，极值2.452E-7，论文结果接近2.4E-7对于