华科大压电材料

一、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.3ANSYS求解(1)单元选择根据该石英板的材料属性，可以选择solid226单元进行求解。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中需要输入的是相对介电常数，因此需要用介电常数除以真空中介电常数ε0=12-108548-.。TBOPT=0则指定恒应变介电系数矩阵[εS]；TBOPT=1则指定恒应力介电系数矩阵[εT]，后者将被转换为恒应变介电系数矩阵[εS]。值得指出的是，solid5单元只能用MP命令制定PERX，PERY和PERZ三项，所以本次计算不采用。ANSYS中输入相对介电系数矩阵的各项时也需要注意顺序，其顺序如下：从而该石英的介电系数矩阵可以按如下输入：tb,DPER,1tbdata,1,4.4285,4.4974,4.5642tbdata,5,0.0971TB和TBDATA指令用于输入矩阵中的各项。其中DPER用于指定该矩阵为各向异性相对介电系数矩阵，另外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,-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]与各项异性弹性刚度矩阵[c]相联系，而压电应变矩阵[d]与弹性柔度矩阵[s]相联系。ANSYS中压电应力矩阵[e]中常数的顺序为：因此命令流可以输入如下：tb,PIEZ,1tbdata,1,0.171!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,SPACE,KFORC,LAYER1,LAYER2,KYNDIV划分后网格如下图：(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!定义上表面电压!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)通用后处理输入“/post1”进入通用后处理部分。本例采用的命令流为：/post1pldisp,2!显示变形前的轮廓和变形后的图形plnsol,d,y!显示厚度方向的电流密度save,YshearStaticClamped,dbfiniPLDISP命令用来显示结构的位移，格式为：PLDISP,KUND。PLNSOL命令用来显示节点解。其格式为：PLNSOL,Item,Comp,KUND,Fact,FileID后处理结果如下：变形图（位移矢量和）：XZ平面（厚度方向）显然产生了剪切。厚度方向电位移图：1.4结果对比自由边界，厚度方向电位移为：hVkD21222262其中391226622226226108171001291082390950....cek，从而26123222262mC100268-010200010823910817121/....hVkD夹持边界，厚度方向上的电位移为：hVD2222由于1108173226.k，所以夹持情况和自由情况下的D2几乎一样大。从ANSYS计算出的电位移图中可知，即使在几何边界处有类似应力集中的效应，但是整体上石英片厚度方向的电位仪是一致的，大小为5-108050-.C/m2，与理论值十分接近。二、Y切石英板的自由振动2.1问题描述计算对象为一块石英板，其尺寸为：长度方向（x轴）：L=0.1m厚度方向（y轴）：H=2h=0.01m宽度方向（z轴）：W=0.02m用ANSYS计算石英板在厚度方向振动的模态，并输出前几阶频率。该石英板的材料常数为：密度ρ=2650kg/m3弹性常数矩阵：压电常数矩阵：介电常数矩阵：2.2理论解控制方程和边界条件为：解出圆频率为：...,,,n,cHnchn321026666因此频率为：...,,,n,cHnf32102266其中，22666661kcc，6622226226cek2.3ANSYS求解(1)建模与常数输入与Y切石英板静态求解部分一样，这里不再赘述。值得注意的是，由于Y切石英板自由振动为厚度方向的剪切振动，因此厚度方向需要划分更密集的网格。在静态求解时仅划分了两层网格，因此求解厚度方向自由振动时不能再采用静态求解时的网格。本例中对石英板的几何尺寸和网格大小都做出了调整，使得求解所得模态可视性很强，结果也更为准确。(2)网格划分本例中石英板的几何尺寸为：长度方向（x轴）：L=0.1m厚度方向（y轴）：H=2h=0.01m宽度方向（z轴）：W=0.02m将厚度方向划分为10层，则单元长度为0.001m。由于长度方向和宽度方向比较次要，因此在这两个方向上单元长度可以更大，从而总节点数不至于太大使得求解困难。另外单元也不可过于狭长，因此长度方向和宽度方向的单元长度可以为厚度方向的两倍，即0.002m。APDL如下：lsel,s,length,,H1!选出长度为H1的线段elesize=1e-3LESIZE,ALL,elesizelsel,inve!选出剩下的线段LESIZE,ALL,2*elesizeVmesh,1!对实体‘1’经行网格划分(3)模态求解对于模态分析，推荐使用BlockLanczos法（缺省）求解。指令为：/soluantype,modal!选取模态分析modopt,lanb,5!定义模态分析选项，展开5阶模态mxpand,5(4)约束与耦合考虑到建模为三维模型，而实际上关注的是Y方向的厚度剪切振动，因此为了使求解效率更高，不产生过多的不需要的模态，这里需要加一些约束与耦合。在厚度方向上的剪切振动中，每一个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!结束循环某一层加上约束及耦合后的情况如下：耦合列表如下：耦合共10组，即厚度方向10层全部分别耦合住。(5)通用后处理后处理我们关心求出的各阶频率及对应的振型，命令流如下：/post1set,list!列表查看所求模态set,first!查看第一阶模态plnsol,u,sum,2!anmode,10,0.1!动画查看振型set,next!查看下一阶模态pldisp,2!anmode,10,0.1列表查看求出的频率：第一阶为刚体位移：第二阶开始的振型分别为：N=1阶频率为166302HzN=2阶频率为332245HzN=3阶频率为498518HzN=4阶频率为665262Hz2.4结果对比(8)各阶振型图对比ANSYS计算出的各阶振型图如上节所示。理论计算出的振型图如下图所示：二者计算结果一致。(9)各阶频率对比理论计算各阶频率由下式决定：...,,,n,cHnf32102266其中，22666661kcc，6622226226cek带入数值得：391226622226226108171001291082390950....cek93-922666661024291081711001291...kcc...,,,n,..ncHnf32102650102429010222966将ANSYS计算的值与理论解汇入下表：阶数n理论解（Hz）ANSYS解（Hz）误差11660771663020.135%23321543322450.027%34982314985180.058%46643096652620.143%可以看出ANSYS计算出的结果具有十分高的精度。三、（大作业）压电复合梁的静态分析3.1问题描述如下图所示，一上下对称的压电复合梁。梁厚度方向极化，上下各有一个电极板，板上加有电压。该复合梁长l=0.2m，中间层的弹性梁厚度为he=0.02m，上下两个电极板厚度都为hp=0.005m。本例为一平面应力问题。用ANSYS对其经行静态分析，并与论文”AnalysisofBeamswithPiezoelectricActuators”中的结果进行对比。中间层的弹性梁E=2.1E11Pa，μ=0.3。压电板的物理参数如下：极化方向为z方向极化，但是考虑到ANSYS建二维模型时为XY平面内，因此可将本例中的XZ平面在ANSYS中看做XY平面，即y方向极化。因此各材料常数矩阵如下：弹性常数矩阵为：29mN1062500003