弹性动力学的相似边界元法

第22卷第6期重庆建筑大学学报Vol 22 No 62000年12月Journal of Chongqing Jianzhu University文章编号∶1006-732X20006-0001-04弹性动力学的相似边界元法程玉民1,彭妙娟2(1.上海大学上海市应用数学和力学研究所,上海200072;2.上海大学土木工程系,上海200072)摘要时论了弹性动力学边界元法中边界单元相似时单元之间的一些矩阵关系建立了相似边界元法的公式。在一组相似单元中只要求得一个单元的相应矩阵通过比例关系即可求得其它单元的相应矩阵然后通过迭加建立代数方程组系数矩阵。与通常的每个单元都各自进行积分计算相比本文方法可大幅度减少计算量。关键词弹性动力学;边界积分方程;相似单元;相似边界元法中图分类号:O343文献标识码:A在边界元法中要建立最终可求解的代数方程组需要在所有单元上进行大量的积分运算。当单元数目较多时单元上积分的计算量将大大增加计算时间。在有限元法的研究中文1时论了单元相似时相似单元之间单元刚度矩阵的关系。在无限元法的研究中虹2地讨论了单元相似时的情况。在边界元法中若能建立相似单元之间相应矩阵的关系则可以减少大量的积分运算为此本文对弹性动力学问题的边界元法讨论了边界单元相似时单元之间的一些矩阵关系建立了弹性动力学的相似边界元法。在一组相似单元中只要求得一个单元的相应矩阵通过比例关系即可求得其它单元的相应矩阵然后通过迭加建立线性代数方程组的系数矩阵。与通常的每个单元都各自进行积分计算相比本文方法可大幅度减少计算量。弹性动力学边界元法对于均匀、各向同性的线弹性材料在小变形条件下运动微分方程为Ccf-c2 ) u; i ( x,t)+c2u; (x,)+f(x,)=u(x,t)xEn其中c1和c2分别为弹性体内膨胀波和畸变波的传播速度;为体力分量为弹性体所在的域,其边界为rx=(x1x2,x3对三维问题咸x=(x1,x2〔对二维问题应力和位移满足如下边界条件x,t)=rin;=p(x,)x∈r(2)u x t)= g(x t)∈和初始条件u(x0+)=;(xx∈g(x0+)=t(x)本构关系为Ii (x ,t)=ol(ci-2c2 )um m( x , t )um m( x,THa中国煤化工(x,t))x∈CNMHG(4)其中r和r分别为已知面力和位移的边界r∪T=Fδ为 Kroneker delta收稿日期2000-05-01基金项目国家自然科学基金作者简介程玉民1965-)男山西人教授博士导师博士主要从事计算力学、大跨空间结构研究。重庆建筑大学学报第22卷对式1)-(4)作用 Laplace或 Fourier变换即得相应的变换域中的方程再应用加权残数法,得到弹性动力学问题在变换域中的边界积分方程Ck Q )u(o)=ukit dr-t hil dr+ufida(5)其中Q为边界点;k和t为基本解:C(Q是自由项与Q点处边界的几何特征有关;;和;分别表示变换域中的位移和面力分量在无体力的情况下边界积分方程5河离散为C(Q)(Q)=221,Nm从5d7t后N,),)d其中N为边界单元数每个单元均为有M个节点的等参单元边界节点总数为K;"n,"分别表示第n个单元第m个节点的i方向的位移和面力洲N()=12灬,M是形函数沃Jacobi行列式K安,)对嵌6进行数值求解即得变换域中边界节点的位移和面力分量然后由数值反变换求得时间域中的解。将边界积分方程6泻成矩阵形式〔CIU=Σ[ GIT)]-∑ HIUn=1其中U〕={u1},i=12…,,K为节点位移向量fUn]={n"},i=12灬…,M为单元n的节点位移向量fTn〕=〔t"〕,i=12灬…,,M为单元n的节点面力向量C〕为与〔U对应的自由项矩阵;〔Gn〔H玢分别为与TnUn]对应的系数矩阵进一步方程(8)写为〔C〕+〔H〕U)=〔GT〕其中T]={t1}i=12灬…,,K为节点面力向量fH〔G是分别由Hnl〔Gn组装而成的对应于〔U〔T)的系数矩阵。2相似边界元法由方程9冋知只要求得矩阼〔H〕和〔G〕,求解方程组即可求解。而〔H〕和〔G〕是分别由〔HnGn咀组装而成的。要求得HnGn〕需要做大量的如下形式的积分thiN(s n) s n sdn=0(10)uNm(,)安dn(11)在非奇异单元上采用Gaus积分在奇异单元上需对奇异积分做消除奇异性或降低奇异性阶数处理或利用刚体位移特解来求奇异积分相似边界元法的基本思想是将弹性体所在的TH为减少形如01的积分的计算量本文V中国煤化工之了相似边界元法。CNMH(条件或几何形状划分为若干局部区域然后将每个局部区域剖分为若干相似单元建立相似单元上矩阵〔Hn〕及〔Gn〕的关系。在一组相似单元中只要求得某个单元上的Hn)和〔Gn)其它单元上的〔H和Gn即可由比例关系得到而不需再做积分运算。下面来推导相似边界元法的公式。第6期程玉民等弹性动力学的相似边界元法将一局部区域剖分为一组相似单元即有i= ax(12)z=az"-(如果需要的话)其中:x"〃和z"为此局部区域中第n个单元的第i个节点的坐标a为比例常数。那么J(13)可得〔Hn]=a2Hn-1(14)〔Gn)=a2Gn-1(15)这样就建立了相似单元之间矩阵〔Hn〕及〔Gn的关系。则在一局部区域中不论有多少单元只要通过积分求得〔H1)〔G1〕由式(14厢和式15)求得此局部区域中所有单元上的〔Hn]和〔Gⅷ而不必再做大量的奇异或非奇异积分运算从而大幅度減少计算量。对二维问题有Jn aJ(16)〔Hn〕=a[Hn-1(17)3算例以下对突加载荷σ作用下的二维矩形中心裂纹板的应力强度因子作了计算。矩形中心裂纹板的模型如图1所示板长4cm宽2cm裂纹长度为2aa=0.24cm。材料的弹性模量E=2.0×105MN/m2泊松比ν=0.3质量密度p=0.005MNs2/m4对I型对称裂纹其动态应力强度因子的计算公式为k(1)=212x(4k1)-)(19)其中对平面应变问题k=34y对平面应力问题k=(3v)(1+v);Kt厢和v(t)分别为裂纹尖端单元上B点和C点在t时刻y方向的位移。计算所得的正则应力强度因子(即K1(t)K1,K1为静态应力强度因子)其与前人所- BEM-PLT做工作的比较见图2。可以看出本文方法的计算结果与其它方法基本吻合。参考文献HE中国煤化工力强度因子CNMHG[1] Leung AYT and Su RKL. Mode I Crack Problems by Fractal Two Level Finite Element Methods J). EngineeringFracture Mechanics 1994 A86)847-856[2]应隆安.无限元方法(M〕北京北京大学出版社1992[3]嵇醒臧跃龙程玉民.边界元法进展及通用程序〔M〕.上海同济大学出版社997重庆建筑大学学报第22卷Similar Boundary Element Method in ElastodynamicsCHENG Yu-min, PENG Miao-juan2(1. Shanghai Institute of Applied Mathematics and Mechanics, Shanghai University Shanghai 200072, Chian 2. Department of Civil Engineering Shanghai University Shanghai 200072 China)Abstract: For the boundary element method of elastodynamics some properties of matrices are discussed in case of similar boundary elements and the similar boundary element method is presented. In aseries of similar boundary elements when the corresponding matrices of a boundary element are obtained the ones of other boundary elements in the series can be obtained by proportion. Then the coefficient matrix of the last system of linear algebraic equations can be obtained by the method of superposition. Compared with the general boundary element method the computing speed can be raised bythe similar boundary element method given in this paperKeywords elastody namics i boundary integral equation similar boundary elements i similar boundaryelement methodH中国煤化工CNMHG