电磁弹性材料Ⅲ型界面断裂分析
2019-01-03周震寰徐旺邓子辰徐新生徐成辉
周震寰, 徐旺, 邓子辰, 徐新生, 徐成辉
(1.大连理工大学 工程力学系 工业装备结构分析国家重点实验室, 辽宁 大连 116024;2.西北工业大学 力学与土木建筑学院, 陕西 西安 710072)
随着材料科学的发展,具有多物理场耦合特性的电磁弹性复合材料已经广泛用于智能器件的开发与制造中[1-2]。该类材料往往由两相或多相压电和压磁材料按照一定的方式复合而成。由于不同材料相间的失配性,电磁弹性复合材料在制造和使用过程中会不可避免地出现界面裂纹。这些裂纹会在外荷载作用下产生应力集中,进而引发结构失效破坏,造成安全事故。因此,研究电磁弹性复合材料的界面断裂问题具有重要的理论和实际意义。
在现有研究中,解析方法往往仅限于求解无限大或半无限大电磁弹性材料中的裂纹问题,如积分方程方法[3]、Stroh变换方法[4]。对于有限几何尺寸的含裂纹电磁弹性材料还主要依赖于数值方法,如有限元方法[5]、边界元法[6]和无网格方法[7]。然而,大部分数值方法仅适用于单相电磁弹性材料,无法直接应用于双相电磁弹性复合材料的界面断裂分析。为解决上述问题,本文提出一种适用于电磁弹性复合材料反平面界面断裂分析的辛离散有限元方法。该方法将一类基于哈密顿体系的解析方法[8]与传统有限元方法相结合,能够简单、高效地获得相关断裂参数,并同时获得裂纹尖端附近奇异物理场的显式表达式。目前,该方法已经成功应用于弹性材料与压电材料的断裂分析中[9]。
本文提出的辛离散有限元方法可以分为2步:第一,根据裂纹尖端的位置将整体结构分为2类区域,即包含奇异性的近场和无奇异性的远场;第二,通过在裂纹尖端附近引入解析的辛本征解函数,实现该区域内的未知量变换,将数量庞大的节点未知量转化为少量的辛本征解待定系数。该方法可以避免传统有限元方法在计算断裂参数时的网格敏感性和路径依赖性,直接提高计算效率和计算精度。
1 电磁弹性材料反平面有限元列式
考虑如图1所示的含界面裂纹的电磁弹性复合材料。上下层材料分别记为材料1(M1)和材料2(M2),z轴选取为电磁弹性介质的极化方向,坐标原点位于裂纹尖端。结构承受反平面剪应力τ0,平面内电位移D0和磁感应强度B0。曲线Γ0将整体结构划分为2类区域,即近场区域和远场区域。
在直角坐标系下,电磁弹性复合材料在反平面荷载作用下的基本方程可以表示为:
本构方程
(1)
几何方程
(2)
平衡方程
(3)
图1 含裂纹电磁弹性复合结构
本文使用的电磁弹性材料单元为八节点四边形单元。单元内任意点未知量由节点未知量表示
(4)
式中,Nj是形函数。单元势能可以表示为
(5)
式中,L(i)为应变能密度函数
(6)
将(1)式、(2)式和(4)式带入(5)式可得
(7)
根据最小势能原理δΠ(i)=0,可得
(8)
式中,Ke是单元刚度矩阵。根据节点编号,组集整体刚度矩阵、节点位移向量和整体载荷向量,可得电磁反平面有限元列式为
(9)
式中,K为刚度矩阵,f为荷载向量,下标N和F分别代表近场和远场区域。
2 电磁弹性材料Ⅲ型界面断裂辛方法
在近场区域,定义全状态向量为
Ψ(i)={q(i),p(i)}T=
(10)
则电磁弹性复合材料Ⅲ型界面断裂问题的哈密顿控制方程可以表示为
(11)
式中,H(i)为哈密顿矩阵,见附录(A-2)。界面连续条件为
q(1)|θ=0=q(2)|θ=0,R(1)gq(1)|θ=0=R(2)gq(2)|θ=0
(12)
裂纹面边界条件为
R(1)gq(1)|θ=π=0,R(2)gq(2)|θ=-π=0
(13)
利用分离变量法求解哈密顿方程(11),并结合(12)式和(13)式得到问题的辛本征值和本征解,见附录(A-3~A-7)。因此,双材料电磁反平面的位移,电势和磁势可以表示为
(14)
式中,M是所取的本征解项数,cj, k为待定系数。
3 辛离散有限元方法
由(14)式,近场区域内节点未知量可以表示为
uN=Φc
(15)
式中,c={c0, 1,c0, 2,c0, 3,c1, 1, …,c1, 6,c2, 1, …,cM, 6}T为辛展开函数的待定系数向量
(m=1, 2, …,NΩN;j=1, 2,…, 3M+6),NΩN为近场区域内总节点数,(rm,θm)是第m个节点的极坐标。将(15)式代入(9)式,则电磁反平面断裂问题的辛离散有限元表达式为
(16)
4 强度因子和能量释放率
由(1)式和辛本征解可知,在裂纹尖端(r=0)[10],广义应力强度因子可以表示为
(17)
式中,K3,KD和KB分别为应力强度因子、电位移强度因子和磁感应强度因子,是由反平面剪应力τ0,平面内电位移D0和磁感应强度B0耦合作用的结果。能量释放率为
(18)
5 数值算例
为验证本文提出方法的正确性与有效性,数值算例计算了3种典型裂纹对应的断裂参数。本节中所有计算数据均采用无量纲形式[10],涉及的材料参数如表1所示。
表1 材料参数
图2 含边裂纹电磁弹性矩阵板
5.1 算例1
表2 不同高度和宽度对应的强度因子和能量释放率
5.2 算例2
图3为1个含有平行裂纹的电磁弹性材料矩形板,其裂纹长度分别为2a1和2a2。该板两端受到均布荷载τ0=1,D0=1和B0=1作用。令W=1,H/W=0.5,H1/W=0.5,表3给出了不同裂纹长度对应的裂纹尖端A处的强度因子与能量释放率。从表中数据可以看出,对于给定的a2,各强度因子与能量释放率总是随着的a1增加而增大;而对于给定的a1,断裂参数表现出相反的变化趋势。表4给出了裂纹长度相同时(a1=a2)断裂参数随H/W的变化规律。从表中可以看到,各断裂参数随着H/W单调变化,并逐渐趋近于某一固定值。该现象可以通过圣
维南原理解释,当H/W增加到一定值时,边缘效应不再对裂纹尖端物理场分布产生影响。
图3 含2条平行裂纹电磁矩阵板
a2/Wa1/W0.20.40.50.60.70.80.2K*30.977 41.068 01.131 51.216 71.344 61.570 4K*D0.979 41.070 11.133 61.218 61.346 11.571 3K*B1.055 41.148 61.210 91.288 81.401 61.606 3G0.342 90.818 91.148 91.594 12.271 23.540 60.5K*30.810 20.946 11.040 71.155 91.308 61.552 6K*D0.817 10.951 41.044 91.159 11.310 81.553 9K*B1.068 71.145 21.201 61.277 71.392 31.600 6G0.235 70.642 60.971 91.438 72.151 23.460 90.8K*30.694 90.818 90.916 71.047 01.226 61.504 1K*D0.705 10.827 70.924 41.053 11.231 01.506 6K*B1.076 91.151 31.204 91.277 51.389 61.597 5G0.173 40.481 60.754 31.180 61.890 43.248 0
表4 裂纹尖端A的强度因子和能量释放率(a1=a2)
5.3 算例3
考虑如图4所示的1个含有折线裂纹的电磁弹性材料矩形板。裂纹与水平方向的偏离角度分别为θA和θB。令W=H,a=0.4H,τ0=1,D0=1和B0=1,表5和表6分别给出了的裂纹尖端A和B处的强度因子和能量释放率随θA和θB的变化规律。从表5中可以看出,所有断裂参数随θA的增大均表现出先增大后减小的趋势。当θA=0时,应力、电位移强度因子和能量释放率达到最大值。表6表现出与表5相似的变化趋势。从上述现象可以看出,合理控制外部电磁场可以有效改变电磁弹性材料裂纹尖端的能量释放率,从而阻滞裂纹扩展。
图4 含有折线裂纹电磁矩形板
θBθA-π/4-π/100π/10π/4-π/6K*30.757 30.968 51.017 10.982 10.789 0K*D0.771 90.979 51.025 10.987 10.788 9K*B1.301 41.383 71.323 81.171 80.781 7G0.412 00.673 60.742 80.692 50.446 90K*30.798 21.019 51.073 61.043 00.855 3K*D0.805 41.022 61.073 61.039 90.847 4K*B1.061 61.135 91.073 70.926 10.554 5G0.457 50.746 10.827 40.780 80.524 9π/6K*30.788 71.017 91.076 81.050 80.869 8K*D0.788 71.013 41.069 11.040 10.854 8K*B0.784 00.853 30.791 30.650 50.300 9G0.446 60.743 60.832 10.792 30.542 7
表6 裂纹尖端B强度因子和能量释放率随θA的变化
6 结 论
本文将辛离散有限元方法成功应用于电磁弹性复合材料在反平面荷载作用下的界面断裂分析中。该方法将传统有限元方法与哈密顿体系辛方法有机结合,可以直接计算出高精度的应力、电场、磁场强度因子和能量释放率,并同时获得裂纹尖端附近区域各物理场的显式表达式。相比其他数值方法,该方法有三方面优势:①在裂纹尖端奇异区域内引入具有辛展开形式的解析解,将该区域内传统有限元对应的大量节点未知量转化为少量辛本征解待定系数,解决了传统有限元方法中断裂分析的网格敏感性问题,并同时大幅提高了计算效率;②引入的辛本征解函数能够直接表征裂纹尖端的奇异性,避免了传统有限元方法中断裂参数计算的路径依赖性问题,直接提高了应力强度因子计算精度;③无需额外引进新的单元以及后处理程序,可以直接获得裂纹尖端附近奇异物理场的显式表达式,有利于该类结构的前期优化设计与后期安全评估。
附录A
(A-1)
(A-2)
式中
(A-3)
(A-4)
(A-5)
(A-6)
(A-7)
式中,φj=rμjsin(μjθ),φj=rμjcos(μjθ);μj=j/2是辛本征值;ρ=(R(2))-1R(1),n=1, 2, 3, …。