广义有限元方法研究进展

第26卷　第1期2009年3月 应用力学学报Vol.26　No.1Mar.2009 CHINESEJOURNALOFAPPLIEDMECHANICS 文章编号:1000-4939(2009)01-0096-013 广义有限元方法研究进展 ＊ 李录贤　刘书静　张慧华　陈方方　王铁军 (西安交通大学　710049　西安) 摘要:广义有限元方法是常规有限元方法在思想上的延伸,它基于单位分解方法,通过在结点处引 入广义自由度,对结点自由度进行再次插值,从而提高有限元方法的逼近精度,或满足对特定问题的特殊逼近 要求 对教师党员的评价套管和固井爆破片与爆破装置仓库管理基本要求三甲医院都需要复审吗 。基于广义有限元方法对单元形状函数构造理论的深入研究,具有任意内部特征(空洞、夹杂、裂纹等)及外部特征(凹角、角点、棱边等)的复杂问题,都将在简单、且与区域无关的有限元网格上加以求解。本文主要介绍广义有限元方法的基本思想、主要特征及对重要细节的处理策略,包括线性相关性的处理、局部逼近函数的获取、区域上的数值积分技术以及边界条件的处理。与扩展有限元方法和有限覆盖方法比较,分析它们各自的特点。综述广义有限元方法的研究现状、应用,展望广义有限元方法的未来发展。关键词:常规有限元方法(CFEM);单位分解法(PUM);广义有限元方法(GFEM);扩展有限元方法(XFEM);有限覆盖方法(FCM) 中图分类号:O34;O242.21　　　文献标识码:　A 二维的三角形单元或三维的四面体单元)或局部坐 标(例如二维的四边形单元或三维的六面体单元)的多项式插值,插值精度通过进一步细化网格(h型), 有限元方法的理论基础是变分原理,以单元和 结点为最基本要素的常规有限元方法,其核心是以所求场变量在结点上的值作为待定参数(自由度),在单元上进行插值逼近,并因此具有诸多卓越的性能,如良好的系统方程性态(即对称性、稀疏性和带状性)、稳定的收敛性、适用于任意复杂的区域和边界条件以及非线性问题,易于处理边界条件等等。相对于其它数值方法而言,有限元方法在理论和实现步骤上都已日臻完善,并诞生了许许多多有限元商用软件,为许多领域大型复杂结构分析提供了重要手段,因此,有限元方法的每一点点发展,都将具有深远的意义。 有限元法在单元内部一般采用整体坐标(例如 或在单元的边、内部增加结点(p型),或二者的结合 h-p型来得到进一步提高。另外,多项式函数在描述不连续特性上的不足,要求有限元方法的网格必须能够描述区域的几何特征、材料变化等,只有足够细的网格才能给出所期望的精度。因而,有限元方法在应用中还存在一些实际困难,如含有成百上千个微小夹杂、空洞和/或裂纹的复杂内部结构问题,又如含有凹角等不光滑边界的区域问题。如何利用有限元法高效高精度地解决上述问题,不仅仅是一个削减网格剖分工作量问题,实质上需要思维方式的变革。近年来有限元方法的重大发展,出发点之一就是克服常规有限元在这方面的不足。本文主要介绍广义有限元方法,并分析它与扩展有限元法和 1　引　言 ＊基金项目:国家自然科学基金(10472090,10572109);国家973项目(2007CB707705);教育部新世纪优秀人才计划(NCET-04-0930)来稿日期:2007-12-04　修回日期:2008-09-22 ,,,,,——。 第1期　　　　　　　　　　　　　　　　　　　李录贤,等:广义有限元方法研究进展97 有限覆盖方法的异同。 广义有限元方法的思想曾零星地出现过[1-8],但只有单位分解方法研究 [12-15] [9-11] 上的形状函数可由下列方式构造 Sτ=φβ×{1,u1,u2}={φβ,φβu1,φβu2},α (β=α,α+1) 其显式 表 关于同志近三年现实表现材料材料类招标技术评分表图表与交易pdf视力表打印pdf用图表说话 pdf 达式为 Sτ={φα,φα+1,φαu1,φα+1u1,φαu2,φα+1u2}(6)α 也就是说,通过这种具有单位分解特性的Lagrange型有限元形状函数与局部逼近u1和u2的积就可构 造广义有限单元的形状函数(加入单位1的原因将在2.2节解释),这样,单元τα共有2×3=6个形状函数,每个结点3个。当然,可以增加u1、u2的数目,进一步拓展空间Sτα的维数。 (5) 出现以后,才开始了其系统 ,并有必要重新审视以前这方面的工作。 在理论上,广义有限元方法是单位分解方法和常 规有限元方法的混合产物,它包括两个主要步骤:一是采用与区域无关的网格,二是引入特定的局部逼近函数构造单元的形状函数。下面予以详细介绍。 2　广义有限单元形状函数的 [12] 构造 2.1　广义有限元单元形状函数构造的基本思路以一维线性有限单元为例,介绍从单位分解法到广义有限元方法的基本思想。 研究定义在区域Ψ R上的场函数u。构造开集 TN=∪ωα α=1 盖,即 TN (2) 设xα为定义在支集ωα上的局部逼近函数空间,u的局部逼近uα为它的元素。假定uα能很好地逼近u,则 局部逼近经某种组合可给出u的整体逼近uhp,并将使uhp与u之间的差在给定范数下以局部误差|u-uα|为界。在单位分解法中,这种整体逼近采用定义 α上的函数φα加以实现。α具有下列特性在支集ωφ α∈C0(αφω),(s≥0;1≤α≤N) s N (1) 图1　一维广义有限元单元分解 使得Ψ的闭集被N个中心位于点xα的支集ωα所覆利用有限元形状函数所具有的单位分解特性,可以很容易地证明:(6)式所定义的形状函数的组合能“再生”局部逼近u1和u2,即 1×(φαu1)+1×(φα+1u1)=u1(φα+φα+1)=u1 (7a) 1×(φαu2)+1×(φα+1u2)=u2(φα+φα+1)=u2 (7b) 换句话说,有 u1,u2∈span{Sτ} 　　借助于支集上所定义函数的单位分解特性将局部逼近空间粘连在一起构造形状函数,再通过线性组合再生支集ωα上的局部逼近,是广义有限元方法形状函数构造的基本思想。这些特征将在2.3节详细讨论。 2.2　广义有限元方法的形状函数及误差估计 考虑n维框架,仍然采用式(1)和式(2)中的记号和关系式。xα(ωα)=span{Liα}i∈I(α)为定义在ωα上的局部空间,其中I(α)为指标集,Liα表示像u1和u2(见2.1节)一样的局部逼近函数,这种函数可按如下选择:假定具有单位分解特性的有限元形状函 ＊数φα仅具有线性插值特性,并且 ααPp-1(ω) Xα(ω),(α=1,…,N)(8)其中Pp-1表示次数低于或等于p-1的多项式空 (3)(4) ∑φ=1,( x αα ∈Ψ) 函数φα称为开集TN上的单位分解。 在有限元类单位分解中,支集ωα就是共有一个顶点结点xα的有限单元的并集 [9,16] 。此时,广义有 限元方法的实现在思想上与标准的有限元程序相 同,主要区别在于下面将要阐述的形状函数的形式。广义有限元方法仍然采用有限单元结点支集的单位分解,因而,它可使用常规有限元法中已发展成熟的程序结构等。 图1表示一维的广义有限元离散。单位分解函数φα就是通常的有限元形状函数,即以结点xα为中心的经典“帽子函数”,因而,支集ω就是单元τ与τα的并集,或者说,单元τα是支集ωα和ωα+1的交集, α α-1 并以两个支集的中心xα和xα+1为结点。现在考察拥有结点xα和xα+1的单元τα,假定与该单元关联的支集{u,＊文[12]中要求φ必须为线性函数似不准确,因为局部坐标系α中的线性函数在整体坐标系中并非就是线性函数,因而,这里改成 98应用力学学报 N(h) 第26卷 间,那么,p次有限元的形状函数定义为F p N ={φ=φαLiα,　(α=1,…,N);　i∈I(α))} (9) 注意,局部空间xα的选取有相当大的自由。xα α i ‖u-uhp‖E(Ψ)≤C((α,h,p,u)∑ε (13) α=1 其中C是一个与h、p、u等无关的常数。 2.3　广义有限单元形状函数实例分析及性质 本节以二维三角形单元为例,对广义有限单元的形状函数进行实例分析,并说明存在于广义有限单元形状函数中的线性相关性。 对于三角形单元,利用线性三角形形状函数与线性单项式的乘积,二次广义三角形有限单元的形状函数可构造为 φα×{1,ξ,η},　(α=1,…,N) (14) 其中φα为常规三角形有限元的线性形状函数(即面积坐标)。ξ= ααα ,η=,x=(xα,yα)hαhα 的基可直接选取能很好逼近光滑函数的多项式函数,此时,广义有限元方法在本质上与常规有限元方 i∈I(α)法相同。对每一个α,由于基函数{Liα}(α=1, …,N)的数目(即I(α)的元素个数)可以不同,因 而,广义有限元网格的每个顶点可选取不同的多项式阶数,并能实现广义有限元与常规有限元间的无缝连接,并不需要像文[5]那样,必须约束广义有限元结点的某些自由度。逼近还可以是各向异性的,例如,不同的方向具有不同的多项式阶数。在广义有限元方法框架内,不再需要常规p型有限元方法中的边结点、内部结点等概念。 在很多情况下,边值问题的解不是一个光滑函数。对于这些情形,像常规有限元方法一样利用多项式构造逼近空间并不是最优,除非使用精心设计的网格,否则可能导致对解u的低劣逼近。在广义有限元方法中,可以使用任何对解的已有知识,以选取较好的局部空间xα。例如,可以求解在部分区域具有点或线奇异性的边值问题,然后使用这样的局部空间xα构造广义有限元的形状函数以反映这些奇异性,比多项式函数更有效。 下面简单介绍广义有限元方法的误差估计。 假定Ψ划分成有限单元、形成单位分解,并满足二维或三维问题对单元及单位分解函数的要求。记 hα=diam(ωα)h=α=max(h)1,…,N α 是结点α的坐标。如此构造可使舍入误差最小[10]。 根据上节中的记号,有 xα(ωα)=span{Liα}=span{1,ξ,η},　　　　　　　　　(α=1,…,N) (15) 现考察具有x1、x2、x3三个结点的三角形单元 τ,由式(14),该单元的广义有限元二次形状函数为 Sτ=φα×{1,ξ,η},　(α=1,2,3) (16) 因而,每一个二次的三角形单元具有3×3=9个形状函数,而不再是常规有限单元中的6个(分别对应于6个基函数:1、x、y、xy、x、y)。应该注意的是,所有的形状函数都定义在单元的顶点上,犹如一个线性单元,没有引入边、面或内部结点自由度等概念,并不像常规高阶二维有限单元的做法。而且,高阶形状函数的支集与线性形状函数{φ1,φ2,φ3}的支集相同,这对刚度矩阵保持其优异的结构特性具有重要意义。 具有“再生”三次多项式的广义三角形形状函数可通过下列方式构造: 如局部空间xα使用基函数 22 xα=span{1,ξη,ξ,η},　(α=1,2,3) 2 2 (10a)(10b) 及X hp =span{φ},　(α=1,…,N;i∈I(α)) (10c) α i αi 其中广义有限单元的形状函数φ已在式(9)中定义。此外,为了保证逼近函数空间的完备性(2.3节将会看到)及广义有限元向常规有限元的退化(3.1节将会看到),在局部逼近函数空间中还需要加入单位常函数1,即 1∈xα(ωα),　(α=1,…,N) 的ε,使得 ‖u-uα‖E(Ψ∩ω)≤ε(α,h,p,u),　(α=1,…,N)α (12) ,hphp (17) 那么,具有x1、x2、x3结点的单元τ的形状函数为 22 Sτ=φ α×{1,ξη,ξ,η},　(α=1,2,3)(18) 　　如式(17),xα的基函数选取时去除了ξ,η的成 (11) 分,一方面由于它们本身可由单位分解函数φα(α=1,2,3)再生得到,另一方面,它们在再生二次多项式时的作用也由于局部逼近函数的二次特性不再需要。然而,这些做法并不足以避免形状函数的线性相关性,因为三次形状函数的基只有10个,即1、x、y、x22、2233,而　　根据局部逼近特征,存在一个与α,h,p,u有关 第1期　　　　　　　　　　　　　　　　　　　李录贤,等:广义有限元方法研究进展99 为3×4=12。具体分析如下。 定理:令Sτ和 Sτ如式(16)和式(17)所定义,那么i)　span{Sτ}=span{1,x,y,xy,x2,y2}=P2ii)　span{Sτ }=span{1,x,y,xy,x2,y2,x2y,xy2, x3,y3}=P3 证明:为了符号简单,将xα中的Liα写成其等价形式{1,x,y},因而 Sτ=φα×{1,x,y},　(α=1,2,3)(19) 　　由于φα为线性三角形单元的形状函数,那么, x存在常数aα,ayα,使得 x∈τ,满足 3α=13α=1 3α=1 x φαa∑α=x 应该小于或等于7-1=6。这样,定理的第i)部分得证。第ii)部分可同法得证。 上述定理说明,二次和三次广义有限元形状函数中总存在线性相关性。文献[5]对广义有限元的 这种特性没有提及,相应的数值结果也就值得商榷。 3　广义有限元方法的基本思想 3.1　广义有限元方法的基本思想 以二维问题为例,常规有限元方法一般将所研 (20)(21)(22) 究的区域划分成三角形或四边形单元组成的网格,并且在划分网格时就考虑了区域的边界及其内部的所有细节(如图2a),对于具有复杂内部结构(如含有成百上千个夹杂或空洞)的问题,网格的剖分就变得非常困难,即使勉强能够实现,也显得非常笨拙。广义有限元方法则绕过所有的内部细节进行单元剖分,生成被内部结构切割的网格。图2给出了两 {1,x,y}∈Sτ 种方法在网格剖分上的比较元方法具有更强的灵活性。 (23)(24)(25) [14] ∑φ=1 α ∑φa y αα=y 因此有 　　根据式(20)～式(22),很容易得到 3 x α α α=13α=1 3 x αα α=13 y α ,很明显,广义有限 φx)=x∑a∑a( φ=x2 与灵活的网格剖分对应的是对区域内部具体细节的处理,这是通过广义有限单元形状函数强大的 可构造性、特别是在以顶点为中心的支集上引入针对特定问题的局部逼近加以实现的。实际上,对含有凹角或裂纹等奇异性的问题,常规有限元插值函数的逼近精度和效率都很低,也需要进行改进。 根据单位分解方法,广义有限单元上的场变量可插值为[13-14] N N k (k)(k)jj FEM ∑a (φαx)=x α=1 ∑a y ααφ=xy 相似地,对于y2项,可得到相应的形式。因此 P2 span{Sα} 的维数为6。 利用式(20)表示的函数φα的单位分解特性,有 3α=1 αφx)=x∑( 3α=1 α=x∑φ 　　下面再证明Sτ是一个线性相关的集合,而且它 (26) uGFEM= k=1 ∑φ(∑a j∈I(α) ψ)+ k=1 φ∑b k k (30) 　　考虑到(21)式与上式,得到 3 α α=1 3 x αα α=1 其中:第一项是仅与单元顶点(或称为外自由度)有关的插值;第二项包括由于边结点和内部自由度的存在 (27) 而产生的高阶插值部分,是将高阶常规有限元拓展至 广义有限元时残留的;NFEM是单元边界自由度和内部自由度的总数;φk和φ k都是常规的有限元插值函数,通常以母单元上的局部坐标(ξ,η)表示。 (k) 式(30)中的ψj就是广义有限元方法在顶点k处具有特色的局部逼近函数,通常以实际单元上的整体坐标(x,y)表示,函数形式则根据具体问题而 (28) j定,将在3.2.1节讨论。可以看到,ψ中必须包含单j位常函数1,而且,当ψ仅包含单位常函数1时,式(30)即退化为常规有限元法的插值形式。 (k) (k) φx)-∑φa∑( =0 　　对y亦有相似的关系式,说明x、y分别有两种构造方法,即参与构造它们的形状函数是线性相关 的。因此,span{Sτ}的基的维数应该小于或等于9-2=7。 此外,还有 3 x α α α=1 3 x αα α=1 φy)=y∑aφ=yx∑a( 3 x α α 3 yα α α=1 考虑到式(24)及上式,得到 α=1 φy)-∑a(φx)=0∑a( (29) 3.2　广义有限元方法中特殊逼近函数的获取及应用3.2.1　特定逼近函数的获取 说明xy也有两种构造方法,即参与构造xy的形状, 100应用力学学报第26卷 研究问题的微分方程类型、几何形状以及边界条件等而不同。例如,对于含M个椭圆形空洞的二维Laplace问题,其控制方程为 Δu=0　　on　Ψ (31 ) 方法予以展示,但对刚度矩阵的零模态数目计算表明[13],它们之间确实存在线性相关性,因而必须采取恰当途径给予解决,具体参阅4.1节。3.2.2　特定局部逼近函 数的应用范围 广义有限元方法中的局部逼近函数直接取自于规定问题(例如3.2.1节提到的椭圆形空洞问题)的封闭解,但同时必须注意选择恰当的使用区域。例如,在含多个空洞情形下,首先将空洞附近的网格进行足够细化,以使每 图4　椭圆形空洞及特殊函数应用的3层顶点 个单元的结点的支集不会同时包含两种或以上的特殊逼近函数。每种特殊逼近函数的使用区域按照下述分层方式确定。 第0层:指支集包含空洞边界的单元顶点。 　　映射将用于获得式(34)。如图3所示,第m个空 洞的区域为Ψm。边界为Γm= Ψm,其上作用Neumann边界条件 Γ=g u·n|m或齐次Dirichlet边界条件 第1层:指支集包含第0层顶点的单元顶点。第2层:指支集包含第1层顶点的单元顶点。特殊逼近函数应用的上述3层顶点示意性地表示在图4中。针对具体问题的理论分析表明,同时对第0层和第1层顶点引入特殊逼近函数,可以显著提高求解精度。 [15] (32)(33) u|Γ=0m 　　因而,第m个空洞对应的局部逼近函数为 amlnrm,　　　　　　　(j=1)ψj,mRe(zm+bmzm),　　　(j=2l) Im(z+cz),　　　(j=2l+1) m)其中zm=rmem,i=-1为虚数单位,(rm,θ为 与第m个空洞相关联的极坐标。式(34)是含椭圆形 iθl m -lmm (k) l -l 4　广义有限元方法的实施策略 4.1　广义有限元方法线性相关性的处理 正如2.3节和3.2.1节所提及的,由于局部逼 近函数大多存在的线性相关性,广义有限元方法的刚度矩阵一般都具有奇异性。例如对于式(34)表示的简谐函数,通过计算刚度矩阵零特征值的数目,Strouboulis等分析了不同(n×n)网格和多项式阶数为p时刚度矩阵的秩的特性,并列于表1。 对于该Laplace问题,常规有限元方法仅包含一个零特征值,而且通过固定一个或更多的点可以消除。然而,在广义有限元方法中,零特征值的个数一般并不知道,因此需要求解以半正定矩阵A为系数的线性方程组,即 Ac=b(35) 　　虽然上述方程的解是不唯一的,但却是存在的。有两种基于直接Gauss消去的方法可供选择,它们都非常有效,而且,额外的计算机耗时可以忽略。具体如下。[13] (34) 空洞问题的解析解,常数am、bm、cm根据式(32)或 式(33)的边界条件确定。 同法可获取其它类型的特殊逼近函数,一些具体实例可参考文献[17]。 从以上获取特殊逼近函数的步骤可以看出,广义有限元方法的特殊函数取自于规定问题的封闭解,并包含根据边界条件确定的系数。 若特殊函数的解析形式无法获取,广义有限元方法还允许引入数值式函数。 (i) 应该注意到,由于ψj的形式及其个数并不是 按照基函数的标准选取,式(30)中的各个形状函数(i)φiψj之间并不一定相互独立。对于多项式情形,已在2.3节中进行了示例分析;对于由复杂局部逼近, 第1期　　　　　　　　　　　　　　　　　　　李录贤,等:广义有限元方法研究进展101 令Aε=A+εI,其中ε是一个大于零的小摄动参数,例如取为10 -10 ,I为单位矩阵。这样,Aε就是正定 ③　Step2,判断 ii ciAci 是否满足设定容差:若满 的,因而也是非奇异的。c可按下列步骤予以计算。①　Step0,赋初值 i=0,　c0=Aεb,　r0=b-Ac0,　i=i+1 ②　Step1,计算 i -1 足,停止;若不满足,继续。④　Step3,更新初值 i 　vi=Azi,　ri+1=r0-∑vi,　i=i+1 j=1 ⑤　Step4,进行迭代,转至Step1。 最后求得的ci就是待求的c。 zi=Ari-1,　ci=c0+∑zj -1 ε j=1 表1　广义有限元方法刚度矩阵零特性分析实例[13](p阶简谐多项式插值,n×n正方形网格) Elementsn21416642561024409616,384 Vertices 2 (n+1)4925812891089422516641 Nullity/TotaldegreesoffreedomP=1P=2P=35/127/2711/7519/24335/86767/3267131/12675259/49923 9/2013/4521/12537/40569/1445133/5445261/21125517/83205 12/2816/6324/17540/56772/2023136/7623264/29575520/116487 P=414/3618/8126/22542/72974/2601138/9801266/38025522/149769 P=516/4420/9928/27544/89176/3179＊140/11979268/46475524/183051 ＊笔者总结了后三列数据之间的规律,认为原来的28611有误,应为3179。 2)　Gauss消去法 利用Duff等在文献[18]提出的稀疏对称不定系统的多波前直接Gauss消去法,可以计算出c。该法已形成Harwell子程序库。 4.2　广义有限元方法中的数值积分技术 对于大多数单位分解方法,数值积分的准确性是困扰研究人员的最主要问题之一,它直接决定着整个数值方法的计算精度。常规有限元方法的插值逼近是逐单元地以映射多项式形式进行构造,相应的单元积分在单元上通过一定阶数的Gauss求积加以实现。但是,基于两方面原因,在广义有限元法中不宜直接使用Gauss求积法则:一是刚度矩阵和载荷向量的计算需要在复杂几何形状上的积分,例如被空洞分割的单元;二是积分所涉及的特殊逼近函数(或导数)常常不光滑或具有奇异性。 为了保证数值积分的精度,Strouboulis等通过自适应地细化各个(子)单元,经历了从嵌入式子划分方法[14]到快速细化方法[15]的发展过程,很好地解决了数值积分问题。笔者认为,这是广义有限元方法研究取得的特色性进展之一。 以空洞问题为例,令τ为至少被一个空洞横穿的单元,若求积分 I[f]= fdτ ∫ τ 令n=1为单元当前的零级子单元数,τk(k=1～n)为该级子单元。②　Step1,零级子单元划分 对被空洞横穿的零级子单元进行划分,使生成的一级子单元,或者全部位于求解域内部、或者全部位于求解域外部、或者满足图5所示的单元划分中止准则。 ③　Step2,一级子单元划分 对满足中止条件的一级子单元,按图6所示的划分策略生成二级子单元。④　Step3,初始估计 对于被空洞横穿的零级子单元τk,使用母单元上的7阶内嵌法则,按二级子单元积分、相加得到一级子单元的积分,或直接在一级子单元进行积分,然后将一级子单元的积分相加,获得整个零级子单元 τ上的积分Iτk和误差Ek。 对于其它零级子单元τk,使用母单元上的7阶内嵌法则以及外推法,得到积分Iτ和误差Eτ。kk 计算总积分I= I ∑I τ k 和总误差E= ∑E τ k 。 ⑤　Step4,判断与控制 dowhile 2 rel>ε (36) 在所有单元中寻找最大误差Emax=maxk(Eτ)。k 如果取得Emax的零级子单元含有一级、二级子,;; 将采用如下快速细化求积方法[15]。 , 102应用力学学报第26卷 置零级子单元的数目n-1+4※n。⑥　重复Step1～Step3,重新计算I和E。⑦　enddo 上述算法主要针对单元内空洞边界为多边形的情形。若单元内的空洞边界为曲线,则借助于混合函数法[15,17] 。 DECUHR算法[19]是计算这类积分的有效途径之一,可在非均匀的单元子划分上进行,与其它自适应方法相比,计算次数更少 [13] 。 图6　快速细化方法中单元τk的再剖分策略 图5　快速细化方法中初始单元的划分中止准则。对于空洞或外边界,求解域(即阴影部分)总在图中箭头所指边界方向的左侧;对于裂纹,两侧都为求解域 4.3　广义有限元方法中边界条件的处理[15] 由于广义有限元方法网格划分的灵活性,插值逼近网格有可能覆盖求解域的边界,因而,必须研究对 Dirichlet边界条件或Neumann边界条件的处理。 如果Dirichlet边界横穿某结点的支集,则对该结点施加相应的Dirichlet边界条件。这样施加边界条件,解的收敛性是可以保证的,但若在逼近函数中引入满足 对于裂纹和凹角周围使用奇异函数进行逼近的情形,则需采用特殊的求积技术[13]。例如,对于如图7所示的L-形区域,ΓD是齐次Dirichlet边界条件,ΓN上及其它剩余边界是Neumann边界条件,网格剖分成均匀的正方形单元。其精确解为u0=rsin) 3 因而,广义有限元的插值函数可表示为 3 (37) ψ (A) A =[r3sin)]φ 3 (38) Dirichlet边界条件的特殊函数,如角点或棱边函数(参考图8),可以大大 提高求解精度。 图7　L-形区域及广义有限单元网格剖分(48个正方形单元);(r,θ)为原点位于A点的极坐标 其中φA是顶点A处的帽子函数。 考虑到有限元刚度矩阵与插值函数的关系,必 3 第1期　　　　　　　　　　　　　　　　　　　李录贤,等:广义有限元方法研究进展103 Newmann边界条件式(32),将以边界积分形式出现在载荷向量中,即 Fk|τ= 5.2　广义有限元方法与扩展有限元方法的异同 广义有限元法与扩展有限元法[30]都是常规有限元法的延伸,它们之间具有相似之处。 1)　相同或相近的网格剖分思想。两种方法在划分 用于插值逼近的有限元网格时,都不考虑区域的内部结构细节,允许空洞、夹杂或裂纹等内边界横穿单 图8　Dirichlet边界条件的处理:在○处可以引入角点特殊函数;在□处可以引入棱边特殊函数 ∫ gφkdΓ τ (39) 　　上述边界积分分成n τ=段来处理就成为Fk|nj=1 nj=1 元。并且,广义有限元法的网格还允许区域的外部轮廓横穿单元,这一思想大大增强了有限元方法网格剖分的灵活性。 2)　相同的理论基础。广义有限元法和扩展有限元法均根据基于单元的单位分解方法发展而来,在单元上对场变量进行局部插值逼近。由于二者均源于与常规有限元方法相同的插值逼近理论,因而它们保持了常规有限元方法的优良性能。3)　处理结构细节的思想相同。这两种方法都寻求在插值函数中引入特殊函数(从而增加结点自由度)的思想来考虑真实结构中各种细节的存在,并发展了与之对应的算法和策略(如数值积分技术等)。4)　相同或相近的单元分解思想。由于形状函数中引入了特定的局部逼近函数,广义有限元法和扩展有限元法在进行数值积分时都采用了单元分解思想。以二维问题为例,广义有限元法常常将单元子划分成四边形单元,扩展有限元法则划分成三角形单元,但这种为了数值积分所进行的单元分解,完全不同于常规有限元法中的网格细化,不会增加问题的总自由度数、扩大问题的求解规模。 广义有限元法与扩展有限元法之间主要有两个区别:①两种单元结点处自由度的物理涵义不同。在引入局部逼近函数的同时,广义有限元方法在结点处引入了广义自由度,其物理意义可解释为以该结点为中心的支集上的局部逼近函数(封闭解的各项)的权重;扩展有限元法附加自由度可解释为结点处增强基函数的权重,是对结点常规结点自由度的补充和扩展。例如对于空洞问题,当采用常规有限元网格时,扩展有限元的附加自由度会自动消失(或通过去除附加自由度或置附加自由度为零实现),很自然地退化成常规有限元方法。但广义有限元方法必须改变局部逼近函数,首先在广义自由度上反映这种变化,才能退化为常规有限元方法;②形状函数的特性不同。虽然这两种方法在形状函数中都引入了针对具体问题的局部逼近函数(对于广义有限元法)或增强函数(对于扩展有限元法),但其特性有着本质差异,是这两种方法的最根本区别所在。以g(s)φ(s)ds=∑∫ s k j g(t)φ(tdt∑∫dt s j k (40 ) 图9　横穿单元的Neumann边界条件及处理方式 其中t为定义在母线段 s上,值位于[-1,1]间的参数(参考图9)。数值积分将利用恰当阶的分片一维Gauss法则来实现。当然,如果被积函数是高度震荡函数或具有奇异性,还需采用自适应求积。 5　广义有限元法的特点及与扩展有 限元法和有限覆盖方法的异同 5.1　广义有限元方法的特点 由于广义有限元方法中的局部逼近函数取自于满足边界条件的规定问题的封闭解,因而它能像常规有限元方法一样,精确施加本质(Dirichlet)边界条件,这正是困扰其它类单位分解方法(例如无网格法[20]、扩展有限元方法[21-22]以及数值流形方法[23])的一个主要问题。 通过对线性求解器的改进,简单而有效地解决了广义有限元方法中系统方程组所固有的线性相关问题(详细参阅4.1节)。类似问题也存在于扩展有限元方法[24]和有限覆盖(数值流形)方法[23,25-26]中,但没有加以认真对待。 广义有限元方法可以对复杂局部逼近函数的数值积分精度实现自适应的控制(详见4.2节),以使它们不影响形状函数的插值精度。其它类单位分解方法,例如无网格法和扩展有限元方法[27] [28-29] ,对 104应用力学学报第26卷 扩展为特征,着眼于内部几何或材料发生变化时表现出的最基本性质,利用整体坐标表示的符号距离函数表征界面两侧位移梯度(应变)所出现的跳跃(详见文献[30]),因而,扩展有限元方法的增强函数具有插值特征,也不随空洞的形状和边界条件而发生变化。对同样问题,广义有限元法中的局部逼近函数是形如式(34)的封闭解的各项,并不具有插值基函数的特性;而且,形式因结构细节、载荷和边界条件而不同。 另外,这两种方法所采用的数值积分技术也不同。由于形状函数通常不再是母单元坐标的多项式函数,被积函数在单元上常常是不光滑的、不连续的、甚至是奇异的,因而必须发展相应的数值积分技术。扩展有限元法在扩展单元内部通过单元分解(例如二维问题,子划分成三角形和四边形),目前都采用高阶Gauss求积进行数值积分,更深入的研究尚未见报道。已有研究 [12-15] 等[26]又将满足协调性和完备性(又称再生性条件)的物理覆盖称为有限覆盖。 在数学覆盖MI上定义权函数wI(x),并满足 wI(x)≥0,　x∈MI wI(x)=0,　x MI 以及在整个数学域上 N M (41) I=1 x)=1∑w( I (42) 其中NM为数学覆盖的总数目。式(41)表明权函数wI(x)是MI上的支函数,式(42)就是逼近连续性要求的单位分解。 在物理覆盖PI上定义覆盖函数uI,这样,整 体位移函数由下式逼近 N M [α][α] n I u(x) u(x)= h I=1i=1 x)u∑∑w( I [α] I (43) 其中:nI是与数学覆盖MI相关联的物理覆盖的数目;α=α(I,i)是I与i的函数,描述数学覆盖MI与物理域Ψi的共有情况。 在有限覆盖方法中也可仿照有限元方法定义诸 如“单元”和“结点”等几何要素。从网格剖分角度,若数学覆盖与物理覆盖重合,有限覆盖方法就退化为常规有限元方法;再若数学覆盖选为有限元结点的支集,则二者完全相同,即就是,“单元”是物理覆盖(或结点支集)的公共区域,“结点”是数学覆盖(或结点支集)的“星”。 广义有限元方法的网格与实际的物理域也可以不相同,犹如有限覆盖方法的数学域与物理域之不同。但是二者仍有较大区别:①出发点不同:广义有限元方法,一方面通过采用规则的网格,克服具有复杂内外结构问题在网格剖分上的困难;另一方面,加入规定问题的封闭解,既易于边界条件的处理,又提高了求解精度。数值流形方法则针对岩土工程中块体结构之间存在的不连续性和大变形问题,通过在物理覆盖上设置待定函数,予以有效解决。两种方法各有特色;②核心问题不同:由于研究对象不同,这两种方法面临的核心问题也不相同。广义有限元方法必须获得对应问题的解析解;数值流形方法则必须很好解决物理覆盖间的接触问题。应该说,两种方法都对各自的核心问题提出了较好的处理策略。另外,这两种方法都面临数值积分和线性相关性两个主要问题。 由于数学域与物理域可以不同,即使采用规则形状的数学覆盖,也可能遇到不规则单元上的数值表明,广义有限元的数值 积分精度与求解精度关联,并且,为了获得更高的积 分精度,发展了新的快速细化自适应数值积分技术,或者建议采用DECUHR等自适应方法等,进行数值积分。 应该提及的是,广义有限元方法基于严格的数学基础,其误差可根据局部逼近函数的误差事先估计;而对扩展有限元方法,这方面的工作难度较大,研究进展也相对较少。5.3　有限覆盖方法与广义有限元方法的异同 1991年和1992年,石根华[31-32]先后撰文,提出了数值流形方法(MM)。1997年,裴觉民翻译、整理出版了有关数值流形方法的专著[23]。1998年Chen等 [25] 建立了高阶数值流形方法的 公式 小学单位换算公式大全免费下载公式下载行测公式大全下载excel公式下载逻辑回归公式下载 。2003年, [26] Terada等将数值流形方法发展成为概念上更加 清晰的有限覆盖方法。这类方法在处理不连续问题及大变形问题时具有独特优势,目前在岩土工程中得到了深入研究和广泛应用。 有限覆盖方法以数学域和物理域开始对方法的描述。所谓数学域就是指与物理无关的数学函数的定义区域;物理域则是物理量的定义区域Ψ,可根据内部的物理或边界特征划分成不同的几个区域Ψi。数学域由部分或全部重叠的有限分片相并而成,这些分片称为数学覆盖;利用有限覆盖方法在描述问题时,数学域(即数学覆盖的并集)不需与物理域一致,但必须完全覆盖物理域。将数学覆盖MI与物理域Ψi的公共区域称为物理覆盖,记做PI,这与数在[α] [26] 第1期　　　　　　　　　　　　　　　　　　　李录贤,等:广义有限元方法研究进展105 理[23],广义有限元方法则采用自适应方法保证积分的精度[14-15]。但对于常规有限元方法,数值积分只需采用标准的高斯积分就可解决,是个相当简单的问题。 为了提高逼近精度,有限覆盖方法常采用高阶覆盖函数;在广义有限元方法中,则对结点自由度通过引入广义自由度进行高阶插值。实际上,这两种做法中的插值基函数均具有线性相关性[26]。正如4.1节所描述,广义有限元方法对基函数相关性导致的系统矩阵的秩的特性进行了详细研究;但在数值流形类方法中,只有文献[26]提及了该问题,并建议采用广义有限元方法中的处理策略,初期的工作[23]、甚至对高阶数值流形方法的专门研究中[26],也没提及线性相关性这一重要问题。 当结构从实体部件到柔性部件过渡时,需要非常多的单元;②为了捕捉角点和区域棱边处的应力集中和奇异性,常规有限元方法的多项式逼近需要使用大量单元,但考虑到总的自由度数目以及网格剖分技术等原因,在这些区域的网格并没有足够细化,以致得到的结果不能令人满意;③生成的三角形或四面体单元,只能是各向同性的,即不同方向不能具有不同的插值多项式阶数。为了从理论上加以改进,与常规有限元方法对照,该文详细介绍了广义有限元方法插值函数的构造思想,并以三维问题为例说明了所构造插值函数本身固有的线性相关性。 文献[15]根据文献[13-14]中提出的思想,以曲多边形域上的Laplace边值问题为例,详细阐述了广义有限元方法的设计思想及实现步骤。研究表明,广义有限元方法具有两个重要优点:①能够在不依赖于区域几何形状的网格上建立整体的插值逼近,即可采用与区域无关的网格;②能够在逼近构造中融入任何特定的局部逼近函数,即可采用特定函数进行局部逼近改进。 事后误差估计是用来评判计算结果的主要工具,因而也是数值方法理论研究的一个重要 内容 财务内部控制制度的内容财务内部控制制度的内容人员招聘与配置的内容项目成本控制的内容消防安全演练内容 。文献[34]对广义有限元方法的事后误差估计进行了研究,并以含597个空洞、外边界为Neumann边界条件的Laplace问题为例,展示了在广义有限元具体计算中的预估算子及其有效性。 考虑到广义有限元方法能很好适用于区域内含大量内部特征(如空洞、夹杂、裂纹等)的问题,Strouboulis等利用从离散求解域得到的“数值式”函数作为局部顶点逼近函数,在相对于内部特征分布而言较粗糙的网格上,能够获得较高的数值精度。该方法的难点就是要解决从求解网格上的“数值式”函数向整体网格“数值式”函数的变换。 广义有限元方法建立在对局部逼近解(函数)信息掌握的基础之上,然而,局部解的信息是不知道的、或者说是模糊的,一般能得到的也只有它所包含的函数空间,在给定情形下到底该使用什么样的逼近函数?在什么情形下常用的多项式逼近最佳?这些都需要解决,以期选取性能优异的局部逼近。Babuska等引入了逼近函数性能判别的定量性指标,并进行了实例分析。结果表明,在一维情形中,如果我们能得到的信息是局部函数位于与x无关的Sobolev类空间,那么多项式局部逼近的性能最佳。该方法可望用来研究高维问题中逼近函数的选取。 [36] [35] 6　广义有限元方法的研究现状、 应用及发展展望 6.1　广义有限元方法的研究现状 广义有限元概念于1983年提出[1],1995年开始得到重视,2000年Babuska研究组结合常规有限元法和单位分解法,首次对广义有限元方法进行系统研究,将广义有限元方法的主要特征界定为:①可以像常规有限元方法一样精确处理本质边界条件,而这恰恰是其它类单位分解方法所面临的一个主要问题;②刚度矩阵和载荷向量中的数值积分精度可以得到自适应地控制,以使积分的误差不影响插值函数的逼近精度,该问题在其它类单位分解方法中没有凸现出来;③通过对线性求解器的简单修改,很好解决了线性方程组中存在的线性相关性问题。以具有复杂角点和裂纹的Laplace问题为例,展示了广义有限元方法在处理复杂几何形状问题、计算误差及计算机耗时诸方面的优势。 同年,他们在文献[14]中给出了用广义有限元方法求解含多个椭圆形空洞Laplace问题的详细实施过程,并认为广义有限元方法的突出优点在于:①能很容易包含反映局部特征的特定函数;②能使用覆盖区域边界的网格;③能克服其它方法(如无网格法)在数值积分和施加Dirichlet边界条件上的困难;④可很容易并入已有的有限元程序中,当单元被不同材料或区域边界横穿时,只需增加特定函数的有关程序语句即可。 文献[12]认为常规有限元所使用的网格剖分, [2,33] [13] 106应用力学学报第26卷 中的多项式定义,并不随所使用的背景网格和单元类型而变,Duarte等[37]在广义有限元方法中提出了一种简单有效、且与网格无关的p-型正交局部逼近函数(增强函数),并用粗糙的四面体网格对具有边界层和内部界面层的问题进行了分析。该方法不需要修改就可用于任意维问题的任何类型单元,而且并不像常规有限元方法一样只能在所使用的网格和单元方向上进行p型增强。 在多年研究基础上,Babuska等[38]以求解线性椭圆方程为例,利用变分原理,从数学角度对无网格法和广义有限元方法这两种最接近的数值方法进行了概括。第二年,他们[39]又对广义有限元方法进行了全面总结,给出了广义有限元方法实施的四部分主要内容。认为广义有限元方法是常规有限元方法在思想上的“推广”,可以涵盖许多种有限元方法。当形状函数取多项式及其它特殊函数时,就可得到扩展有限元方法(文献[30])。值得一提的是,广义有限元方法在解决复杂问题时非常有效,但对简单的标准测试问题,并不能体现出它的优势。 就在最近,Duarte等[40]采用基于Ck(k为任意大数)有限单元的单位分解,建立了具有以下特征的广义有限元形状函数:①任意光滑性;②适用于非结构性网格并能推广到高维;③与Lagrange型形状函数相同的支集;④可使用与常规有限元方法相同的数据结构。这种CGFEM形状函数可用于对光滑性要求较高的问题(如受集中力和集中矩等分布数据的问题),还去除了C∞GFEM形状函数对支集的凸性要求。由于具有相同的支集,CGFEM与C0GFEM(原来的Lagrange型)一样都可使用Gauss求积法则,但需要较多的求积点。鉴于此,他们又提出了一种技术,仅在需要高光滑性的区域引入CGFEM,从而结合使用CGFEM与CGFEM两种单元,既简单、又通用。6.2　广义有限元方法的应用 广义有限元方法的研究目标就是以有限元方法为基础,高效、高精度地解决复杂实际问题,因而在方法本身研究和发展的同时,广义有限元方法已在不同领域得到了应用。 Plaks等利用广义有限元法模拟了磁纳米颗粒的自组装问题,使复杂问题数值模拟的困难,从网格剖分和防止扭曲上转移到了数值积分上。 Simone等 [43][42] k k [41] k k 利用插值函数的单位分解特性嵌入到有限元方法中,因而,有限单元的网格不需要与多晶的拓扑一致。该法与文献[44]中处理同类问题的扩展有限元方法比较,只在构造局部逼近函数的原理上有区别,却更容易处理文献[29]中同样的分支和相交裂纹问题。 Strouboulis等[45]借助于广义有限元方法,通过在顶点处的插值逼近中引入平面波函数,分析了高波数(k≤64)的Helmholtz方程。 Pereira等[46]根据与求解方法及有限元离散无关的序列映射,利用广义有限元方法,提出了从环路积分法(CIM)、截断函数法(CFM)和J积分法中提取应力强度因子的实施 方案 气瓶 现场处置方案 .pdf气瓶 现场处置方案 .doc见习基地管理方案.doc关于群访事件的化解方案建筑工地扬尘治理专项方案下载 。结果表明,这三种方法的收敛速度快于广义有限元解的能量范数误差,而且CIM和CFM比J积分法更好。 Barros等[47]在广义有限元方法基础上,采用平衡单元残差法,定义了整体和局部的误差预估算子,并以平面弹性问题为例,展示了预估算子的有效性。在文献[48]中将广义有限元方法延伸至处理由于累计损伤产生的材料非线性问题,并因此提出了一个特定的p-型自适应策略。 Duarte等利用广义有限元方法进行了三维动态裂纹扩展模拟,不需对网格做任何重构。广义有限元形状函数中的单位分解函数,综合了不连续的Shepard单位分解和有限元的单位分解(称之为FE-ShepardPU)特性而构造,以发挥各自在裂纹几何建模和数值积分方面的优点;局部逼近函数可以是多项式函数或用户定义函数,均能有效近似裂纹尖端奇异场。当裂纹扩展时,只需修改沿裂纹表面的单位分解函数,并不需要对网格的连续重构或对不同网格上解之间的映射。数值算例展示了该法在动态裂纹扩展模拟时的主要特征和计算效率。6.3　广义有限元方法的发展展望 准确地说,广义有限元方法是一种数值求解思想。广义有限元方法的形状函数由具有单位分解特性的函数和对具体问题的局部逼近函数相乘而构造,只与问题的物理本质和几何特征有关,均独立于所采用的网格,因而可以采用非常简单直观、相对粗糙的规则单元和网格,以利于进行数值积分。但是,由于在简单的标准问题上难以很好展示其优势,广义有限元方法的有效性一般首先通过严格的数学证明给出,因而涉及的都是给定的具体问题、并针对明确的边界条件,较难被广大工程技术人员接受。展望广义有限元方法的未来发展,尚有诸,。 [49] 利用广义有限元方法,对多晶体分 析中的晶体颗粒边界进行了明晰的处理。认为晶粒, 第1期　　　　　　　　　　　　　　　　　　　李录贤,等:广义有限元方法研究进展107 1)　广义有限元方法和扩展有限元方法分别为不同学派所提出,应各取二者之长。广义有限元方法将局部逼近网格(背景网格)与求解域分离,实现了有限元方法与单位分解方法间的巧妙结合,可以解决复杂内部结构问题和不规则外部边界问题。这些特点,极有可能促进扩展有限元方法的发展或二者之间的结合。2)　广义有限元方法严格的数学理论,是对有限元理论的重大发展,基于此,有望创造出新类型的有限单元,并对过去有限元方法的体系进行重新梳理。3)　复杂形状函数插值精度的判断。对于多项式插值,广义有限元方法通过再生性研究,可以预估形状函数的插值精度。但是,广义有限元方法和扩展有限元方法在形状函数中都引入了复杂的函数形式,对新的形状函数的插值精度,有待进一步研究。4)　广义有限元方法形状函数的线性相关性,说明了形状函数构造的不唯一性。这种不唯一性对数值计算结果的影响需要研究,以更加明确。 5)　通用的广义有限元方法的框架有待搭建。为此,需要开展相关研究,以建立完善的局部逼近函数库,提出数值积分、线性相关性处理以及边界条件处理策略选取准则等。 6)　借鉴数值流形方法处理覆盖接触的思想,发展广义有限元方法,以应用于诸如块体结构等具有较大变形的复杂结构问题,发挥两种方法的各自优势 [50] [43] mentmethod:basictheoryandapplications[J].ComputerMethodsinAppliedMechanicsandEngineering,1996,139:289-314. [10]DuarteCA,OdenJT.AnH-Padaptivemethodusingclouds [J].ComputermethodsinAppliedMechanicsandEngineer-ing,1996,139:237-262. [11]BabuskaI,MelenkJM.Thepartitionofunitymethod[J].In-ternationalJournalforNumericalMethodinEngineering,1997,40:727-758. [12]DuarteCA,BabuskaI,OdenJT.Generalizedfiniteelement methodsforthree-dimensionalstructuralmechanicsproblems[J].Computer&Structures,2000,77:215-232.[13]StrouboulisT,BabuskaI,CoppsK.Thedesignandanalysis ofthegeneralizedfiniteelementmethod[J].ComputerMeth-odsinAppliedMechanicsandEngineering,2000,181:43-69. [14]StrouboulisT,CoppsK,BabuskaI.Thegeneralizedfiniteel-ementmethod:anexampleofitsimplementationandillustra-tionofitsperformance[J].InternationalJournalforNumeri-calMethodsinEngineering,2000,47:1401-1417. [15]StrouboulisT,CoppsK,BabuskaI.Thegeneralizedfiniteel-ementmethod[J].ComputerMethodsinAppliedMechanicsandEngineering,2001,190:4081-4193. [16]OdenJT,DuarteCA,ZienkiewiczOC.Anewcloud-basedhp finiteelementmethod[J].ComputerMethodsinAppliedMe-chanicsandEngineering,1998,153:117-126. [17]SzaboBA,BabuskaI.FiniteElementAnalysis[M].Wiley, NewYork,1991:169-173. [18]DuffIS,ReidJK.Themultifrontalsolutionofindefinite sparsesymmetriclinearsystems[J].ACMMathematicalSoftware,1983,9:302-325. [19]EspelidTO,GenzA.DECUHR:Analgorithmforautomatic integrationofsingularfunctionsoverahyperrectangularre-gion[J].NumericalAlgorithms,1994,8:201-220. [20]KrongauzY,BelytschkoT.Enforcementofessentialbounda-ryconditionsinmeshlessapproximationsusingfiniteelements[J].ComputerMethodsinAppliedMechanicsandEngineer-ing,1996,131:133-142. [21]JiH,DolbowJ.Onstrategiesforenforcinginterfacialcon-straintsandevaluatingjumpconditionswiththeextendedfi-niteelementmethod[J].InternationalJournalforNumericalMethodsinEngineering,2004,61:2508-2535. [22]MoesN,BechetE,TourbierM.ImposingDirichletboundary conditionsintheextendedfiniteelementmethod[J].Interna-tionalJournalforNumericalMethodsinEngineering,2006,67:1641-1669. [23]石根华.数值流形方法与非连续变形分析[M].裴觉民,译.北 京:清华大学出版社,1997. [24]DolbowJ,MoesN,BelytschkoT.Discontinuousenrichment infiniteelementswithapartitionofunitymethod[J].FiniteElementsinAnalysisandDesign,2000,36:235-260. ],Y,Thigher-man- 。 参　考　文　献 [1]　BabushkaI,Osborn.Generalizedfiniteelementmethods:their performanceandtheirrelationtomixedmethods[J].SIAMJour-nalofNumericalAnalysis,1983,20(3):510-535. [2]　梁国平,何江衡.广义有限元方法-一类新的逼近空间.[J]力学 进展,1995,25(4):562-565. [3]　栾茂田,田荣,杨庆.广义结点有限元法[J].计算力学学报, 2000,17(2):192-200. [4]　田荣,栾茂田,杨庆.高阶形式广义结点有限元法及其应用[J]. 大连理工大学学报,2000,40(4):492-495. [5]　邵国建,刘体锋.广义有限元及其应用[J].河海大学学报:自然 科学版,2002,30(4):28-31. [6]　栾茂田,田荣,杨庆.平面广义四结点等参元GQ4及其性能探 讨[J].力学学报,2002,34(4):578-584. [7]　张振宇,董兴平.弹塑性力学问题的广义有限元法[J].盐城工 学院学报:自然科学版,2004,17(4):19-21. [8]　彭自强,李小凯,葛修润.广义有限元法对动态裂纹扩展的数值 模拟[J].岩土力学与工程学报,2004,23(18):3132-3137.[9], 108应用力学学报 ActaNumerica,2003,12:1-125. 第26卷 ifoldmethod[J].InternationalJournalforNumericalMethodsinEngineering,1998,43:685-712. [26]TeradaK,AsaiM,YamagishiM.Finitecovermethodforlin-earandnon-linearanalysesofheterogeneoussolids[J].Inter-nationalJournalforNumericalMethodsinEngineering,2003,58:1321-1346. [27]DolbowJ,BelytschkoT,NumericalintegrationoftheGaler-kinweakforminmeshfreemethods[J].ComputationalMe-chanics,1999,23:219-230. [28]XiaoQZ,KarihalooBL.ImprovingtheaccuracyofXFEM cracktipfieldsusinghigherorderquadratureandstaticallyadmissiblestressrecovery[J].InternationalJournalforNu-mericalMethodsinEngineering,2006,66:1378-1410.[29]DauxC,MoesN,DolbowJ,etal.Arbitrarybranchedandin-tersectingcrackswiththeextendedfiniteelementmethod[J].InternationalJournalforNumericalMethodsinEngineering,2000,48:1741-1760. [30]李录贤,王铁军.扩展有限元及其应用[J].力学进展,2005,35 (1):5-20. [31]ShiG.Manifoldmethodofmaterialanalysis[C]//Trans9th ArmyConfonAppliedMathematicsandComputing,Minne-apolis,Minnesota,1991:57-76. [32]ShiG.Modelingrockjointsandblocksbymanifoldmethod. RockMechanics[C]//Proc33rdUSSymposium,SantaFe,NewMexico,1992:639-648. [33]MelenkJM.Ongeneralizedfiniteelementmethods[D].Uni-versityofMaryland,1995. [34]StrouboulisT,ZhangL,WangD,etal.Aposteriorierrores-timationforgeneralizedfiniteelementmethods[J].ComputerMethodsinAppliedMechanicsandEngineering,2006,195(9-12):852-879. [35]StrouboulisT,ZhangL,BabuskaI.Generalizedfiniteelement methodusingmesh-basedhandbooks:applicationtoproblemsindomainswithmanyvoids[J].ComputerMethodsinAp-pliedMechanicsandEngineering,2003,192(28-30):3109-3161. [36]BabuskaI,BanerjeeU,OsbornJE.Onprinciplesforthese-lectionofshapefunctionsfortheGeneralizedFiniteElementMethod[J].ComputerMethodsinAppliedMechanicsandEngineering,2002,191(49-50):5595-5629. [37]DuarteCA,BabuskaI.Mesh-independentp-orthotropicen-richmentusingthegeneralizedfiniteelementmethod[J].In-ternationalJournalforNumericalMethodsinEngineering,2001,55(12):1477-1492. [38]BabuskaI,BanerjeeU,OsbornJE.Surveyofmeshlessand generalizedfiniteelementmethods:aunifiedapproach[J]. [39]BabuskaI,BanerjeeU,OsbornJE.Generalizedfiniteelement methods-mainideas,resultsandperspective[J].InternationalJournalforNumericalMethodsinEngineering,2004,1(1):67-103. [40]DuarteCA,KimDJ,QuaresmaDM.Arbitrarilysmoothgen-eralizedfiniteelementapproximations[J].ComputerMethodsinAppliedMechanicsandEngineering,2006,196:33-56. [41]EdwardsHC.Aparallelinfrastructureforscalableadaptivefi-niteelementmethodsanditsapplicationtoleastsquaresC1collocation[D].TheUniversityofTexasatAustin,Austin,TX,USA,1997. [42]PlaksA,TsukermanI,FriedmanG,etal.Generalizedfinite-elementmethodformagnetizednanoparticles[J].IEEETransactionsonMagnetics,Part1,2003,39(3):1436-1439. [43]SimoneA,DuarteCA,VanderGiessenE.AGeneralizedFi-niteElementMethodforpolycrystalswithdiscontinuousgrainboundaries[J].InternationalJournalforNumericalMethodsinEngineering,2006,67(8):1122-1145. [44]SukumarN,SrolovitzDJ,BakerTJ,etal.Brittlefracturein polycrystallinemicrostructureswiththeextendedfiniteele-mentmethod[J].InternationalJournalforNumericalMeth-odsinEngineering,2003,56:2015-2037. [45]StrouboulisT,BabuskaI,HidajatR.Thegeneralizedfinite elementmethodforHelmholtzequation:Theory,computa-tion,andopenproblems[J].ComputerMethodsinAppliedMechanicsandEngineering,2006,195(37-40):4711-4731. [46]PereiraJP,DuarteCA.Extractionofstressintensityfactors fromgeneralizedfiniteelementsolutions[J].EngineeringA-nalysiswithBoundaryElements,2005,29(4):397-413. [47]BarrosFB,ProencaSPB,deBarcellosCS.Onerrorestimator andp-adaptivityinthegeneralizedfiniteelementmethod[J].InternationalJournalforNumericalMethodsinEngineering,2004,60(14):2373-2398. [48]BarrosFB,ProencaSPB,deBarcellosCS.Generalizedfinite elementmethodinstructuralnonlinearanalysis-ap-adaptivestrategy[J].ComputationalMechanics,2004,33(2):95-107. [49]DuarteCA,HamzehON,LiszkaTJ,etal.Ageneralizedfi-niteelementmethodforthesimulationofthree-dimensionaldynamiccrackpropagation[J].ComputerMethodsinAppliedMechanicsandEngineering,2001,190(15-17):2227-2262. [50]裴觉民.数值流形方法与非连续变形分析[J].岩石力学与工程 学报,1997,16(3):279-292.