（工程力学专业论文）弹塑性力学问题中的无单元伽辽金法.pdf

摘要 s i 4 4 2 2 无单元伽辽金法是最近几年发展起来的与有限元相似的一种数值方法。它 采用移动最小二乘法构造形函数，从能量泛函的弱变分形式中得到控制方程， 并用罚函数法施加本质边界条件，从而得到偏微分方程的数值解。该法只需节 点信息，不需将节点连成单元；此外，还有精度高、后处理方便等优点。 无单元伽辽金法的数学基础是移动最小二乘法。用移动最小二乘法构造形 函数时，只需在求解的区域内布置一系列的节点。因此，无单元伽辽金法可以 不需单元。但是，移动最小二乘法的近似函数不一定精确地通过计算点，除非 使用奇异的权函数。因此，本质边界条件的施加和集中载荷的处理变的复杂。 但与这种方法带来的优势相比，是微不足道的。 无单元伽辽金法现已成功地应用于弹性力学问题。本文将无单元伽辽金湮 推广应用于求解弹塑性力学问题。首先，利用变分原理推导了微分形式求解群 塑性问题的离散方程。另外，进一步讨论了如何用无单元伽辽金法处理集中 的问题。求解实际中，可用增量形式代替微分形式，并且采用n e w t o n r a p h s o i 增量迭代法进行计算。然后，编制了相应的计算程序，并给出了算例。所有簧 例的结果都与a n s y s 的结果进行了比较，两种计算的结果非常吻合，这验词 了本文理论的可靠性。同时，算例还表明：无单元伽辽金法求解弹塑性问题帝 仍具有稳定性好，收敛快的优点。最后，对无单元伽辽金法的优点和一些关锤 问题也进行了讨论。 关键词：移动最小二乘法，无单元伽辽金法，罚函数法，弹塑性问题， 牛顿一拉普森 a b s t r a c t e l e m e n t f r e eg a l e r k i nm e t h o d ( e f g m 、，s i m i l a rt of i n i t ee l e m e n tm e t h o d ，i s n e wn u m e r i c a lm e t h o dd e v e l o p e dr e c e n t l y i ne f g m ，i no r d e rt og e tan u m e r i c a s o l u t i o nf o rap a r t i a ld i f i e r e n t i a le q u a t i o n ，t h es h a p ef u n c t i o ni sc o n s t r u c t e db ! m o v i n gl e a s ts q u a r e s ( m l s l ，t h ec o n t r o le q u a t i o ni sd e r i v e df r o mt h ew e a kf o r mo v a r i a t i o n a le q u a t i o na n de s s e n t i a lb o u n d a r yc o n d i t i o n s a r e i m p o s e db yp e n a l t , f u n c t i o nm e t h o d t h ea d v a n t a g e so fe f g m a r e ：f 11 o n l yn o d a ld a t a sa r en e c e s s a r y ， iet h e r ei sn on e e dt oi o i nn o d e si n t oe l e m e n t s ：( 2 1h i g ha c c u r a c yc a nb ea c h i e v e d ： ( 3 ) p o s t p r o c e s si se a s y ，e t c t h em a t h e m a t i c a lb a s i so fe f g mi sm o v i n g1 e a s ts q u a r e sm e t h o dt ou s em l ! i ti so n l yn e c e s s a r yt oc o n s t r u c ta na r r a yo fn o d e si nt h ed o m a i nu n d e rc o n s i d e r a t i o n j u s tb e c a u s eo ft h i s e f g mi sc o m p l e t e l yf r e e m o v i n gl e a s ts q u a r e si n t e r p o l a n t sd ( n o tp a s st h r o u g ht h ed a t ab e c a u s et h ei n t e r p o l a t i o nf u n c t i o n sa r en o te q u a lt ou n i t va t h en o d e su n l e s st h ew e i g h tf u n c t i o n sa r es i n g u l a r t h i si so fd i s a d v a n t a g ei ne f g f o ri t c o m p l i c a t e s t h e i m p o s i t i o n o fe s s e n t i a l b o u n d a r yc o n d i t i o n s a n dt hc a p p l i c a t i o no fp o i n tl o a d s h o w e v e r ，t h e s ed i s a d v a n t a g e sa r eh e a v i l yo u t w e i g h t e db ! i t sa d v a n t a g e s e f g mh a sb e e ns u c c e s s f u l l yu s e dt os o l v ee l a s t i cp r o b l e m si nt h i s p a p e r ， e f g mi s a p p l i e dt os o l v ee l a s t o p l a s t i cp r o b l e m s f i r s t l y ，t h ed i s c r e t ee q u a t i o n si r d i 髓r e n t i a lf o r m ，w h i c hi su s e dt os o l v ee l a s t o p l a s t i cp r o b l e m s i sd e r i v e du s i m v a r i a t i o n a lp r i n c i p l ei na d d i t i o n ，i ti sd i s c u s s e dh o wt os o l v ep r o b l e m si n v o l v i m p o i n t 1 0 a d sw i t h e f g m d u r i n gt h es o l v i n gp r o c e s s t h e i n c r e m e n tf o r r nc a l a p p r o x i m a t et h e d i f f e r e n t i a lf o f i na n dn e w t o n r a p h s o ni t e r a t i o n t e c h n i q u e s a r ( i n t r o d u c e di n t ot h e c o m p u t a t i o n t h e nt h ec o r r e s p o n d i n gc o m p u t e rp r o g r a mi ! d e v e l o p e da n ds e v e r a le x a m p l e sa r eg i v e n t h er e s u l t so f a l le x a m p l e sa r ec o m p a r e c w i t ht h o s eo fa n s y sa n dt h e ys h o wg o o da g r e e m e n t t h ec l o s e n e s so ft h er e s u l t ： o b t a i n e db yt h e s et w om e t h o d sv e r i f i e st h er e l i a b i l i t yo ft h et h e o r yi nt h e p r e s e n p a p e ra n d a l s os h o w st h a ti ns o l v i n ge l a s t o p l a s t i cp r o b l e m s ，e f g ms t i l l p o s s e s s e , , s o m ea d v a n t a g e ss u c ha sg o o ds t a b i l i t ya n dh i g hr a t eo f c o n v e r g e n c e a tl a s t t h ( a d v a n t a g e so fe f g m a n ds o m e k e yi s s u e sa r ea l s od i s c u s s e di nt h i sp a p e r k e yw o r d s ：m o v i n g l e a s ts q u a r e s ，e l e m e n tf l e eg a l e r k i n m e t h o d ，p e n a l wf u n c t i o r m e t h o d ，e l a s t o p l a s t i cp r o b l e m s ，n e w t o n r a p h s o n 陈莘苹弹塑性j j 学问题的九啦兀伽辽金法( e l ? o m ) 1 1 引言 第一章绪论弟一早珀下匕 工程实际中许多力学问题都可以归结为在给定的边界条件下，求解一组偏 微分方程组。在理论上，这种边值问题有唯一确定的解，但一般难以求得解 折解。这是由于边界的几何形状或问题本身的一些特性很复杂。克服这些困 难的补救办法是对问题作较多的简化假设，使问题能够求解，但是这样做的 结果往往导致精度太差，有时甚至得出错误的解答。 以前，在得不到解析解的时候，人们或者采用差分法，按差分格式离散以 获得数值解；或者按问题的特点，选取试函数，采用里兹法或伽辽金法等近 似方法来获得近似解。这些近似法总有这样或那样的缺点而不能令人满意。 现在由于电子计算机的迅速发展和计算方法的新进展，可以在保留问题复杂 性的前提下设法去寻找它的近似解。 有限元法的思想最早出现在c o u r a n t1 9 4 3 年所发表的一篇著作中。当时 由于受到一些客观条件的限制而未能得到很快的发展。到5 0 年代，由于工程 分析的需要，计算工具和计算方法都已具备了一定的条件，有限元法在分析 复杂的航空结构中最先得到应用，并且逐步显示出了巨大的优越性，迅速被 众多的科学家和工程师所接受。有限元法把差分法的离散改造成更为灵活的 有限元离散，把里兹法全域内的试函数近似换成局部区域( 即单元) 内的插值函 数近似，以变分原理作为推导的根据，并充分利用电子计算机的计算能力， 从而开拓了现代数值方法的广阔领域。 有限单元法是对某些工程问题求得近似解的一种数值分析方法。这种方 法是将所要分析的连续物体或工程结构分割为很多较小的区域( 称为单元或 元素) ，这些单元的集合体就代表原来的物体或工程结构，然后建立每个单元 的有关特性的关系式，再组合起来就能求得相应物体或工程结构问题的解答。 湖南大学坝! j 学 ? r 造文 这是一种从部分到整体的方法，分析过程人为简化。从数学角度来说，彳限 单元法是从变分原理或加权残数法出发，通过区域剖分和分片插值，把微分 方程的边值问题化为等价的一组线性代数方程来求解。 在有限元法中，最终求解的是线性代数方程组，它的系数矩阵总是对称 的，对于正定的变分问题，有限元离散化后保持了正定性，而且有限元法的 系数矩阵是稀疏的。有限元法不仅适应复杂的几何形状和边界条件，而且很 容易通过对不同的单元规定不同的性质，成功地用于多种介质和非均匀连续 介质的问题。这是其它数值方法最难于处理的问题。人们已用它来求解各种 力学和非力学问题，线性和非线性问题，均取得很好的成效。有限元法特别 适合于求解大型复杂结构的静力学和动力学问题。有限元法还允许把求出各 种问题的程序纳入到一个程序系统以形成通用程序包。现在功能齐全的大型 通用程序包已经商品化，在科学研究和工程应用中起到了愈来愈大的作用。 尽管有限元法所取得的成就与目俱增，但有限元法还不是十全十美的， 改进有限元法的努力一真进行着。但有限元法的某些不足是有限元法固有的， 是无法克服的。冽如：有限元法不大适合求解无限边界场域边值问题，而只 能求解有界问题，因为用有限个单元离散无限域显然是不可能的。因此对无 限域问题只能人为地截取有限域。用有限元法难于处理的另一类问题是域内 具有应力奇异的问题。在固体力学问题中，这类应力奇异通常发生在不规则 的凹角或孔洞附近。由于应力奇异可能引起断裂扩展，因此在奇异点附近能 否得到一个较为精确的解答，有时就显得十分重要。但是在奇异点处，理论 上应力为无限大，用有限元法可能产生毫无意义的分析结果。类似的情况也 可能发生在有集中载荷处，如在势论问题中有点源存在处。在这些情况下， 不论单元划分得多么细小，用有限元法得到的结果通常不可能反映出它们在 奇异点附近的迅速变化。有限元法需要在整个求解域上进行离教，虽然目前 有一些网格生成器，但要对形状复杂的三维体进行网格剖分，仍然不是一件 轻松的事，而且导致问题的自由度和原始信息量大。有限元法的离散技术本 身也存在着缺陷，它把本来是连续的介质用仅在节点处连接的有限个单元的 集合来模拟，这样不仅带来了离散的误差，而且在单元间的连续性要求较高 ， 【：缸苹芋：弹塑性力学问题的无o p 元伽辽金法( e f g m ) 时，单元的构造很困难。另外，用位移型钉限元法求解日i 的应力精度低立 移的粘度以及对于不可压缩物体存在体积自j 锁现象等。 为了弥补有限元法的不足，许多科学工作者不断地致力于研究新的数值 方法。边界元法是继有限元法之后的一利r 别具特色的新的数值方法，它是将 描述弹性力学问题的偏微分方程边值问题化为边界积分方程并吸收有限元的 离散化技术而发展起来的。将弹性力学问题归结为求解一组边界积分方程， 若在边界上已知了三个位移分量和三个面力分量中的三个分量，则由边界积 分方程可以确定其余三个未知分量，而任意内点的位移和应力可由6 个边界 分量通过边界积分来确定，这就是边界积分方程方法。边界积分方程有奇异 性，解析求解非常困难。有限元法所取得的成就吸引人们对边界积分方程在 边界上划分单元进行离散，然后由全部边界节点的三个已知边界分量求出全 部边界节点的另外三个边界分量，这就是边界元法的由来。边界元法中包含 有限元法的思想，它把有限元法的按求解域划分单元离散的概念移植到边界 积分方程方法中，但边界元法不是有限元法的改进或发展，边界元法与有限 元法存在着本质的差异。 边界元法具有有限元法所没有的优点。由于边界元法的离散处理仅涉及 边界，整个域内不再出现待求参数，因此，其待求参数的数目可以比有限元 法( 需同时将全域和边界离散) 所用的少很多，使方程规模缩小，故边界元法可 以用较少数量的未知数分析有限元法同样的问题。边界元法的这个特性使得 它在三维问题中特别具有吸引力，因为在三维问题中，求解区域的外表面对 体积的比值是很小的。由于边界元法能使问题降一维，并且分析同样一个问 题比有限元法简化了输入数据的准备工作，因此，边界元法一般能节省计算 讥内存、机时和人的数据准备工作量，使解题较为经济。边界元法在得出边 界近似解之后，虽得不到解析显式，但可以逐点计算域内点的近似解，而有 限元法则必须同时对所有域内的结点联立求解。因此，当只需对个别点求解 时，边界元法较简便。 由于边界元法采用无限域的基本解自然满足场域无限远处的条件，用边 界元来求解有限元难以求解的无限域问题是非常合适的。另外，对应力奇异 问题也特别适用，易于处理，边界元洲、在理论上能够计算任意点处的解答， 不论这些点在非常遥远处或是在离奇异点为任意小的距离处，因为边界元实 际上摆脱了有限元法中存在的这样一个约束，即必须在一个给定的网格各点 上寻求问题的解答。 实际上，任何一种方法都不是十全十美的，边界元法也有其自身难以克 服的缺点。用边界元法求解边值问题需要找到控制微分方程的个基本解或 控制微分方程在无限空间上的格林函数，这对于某些问题是十分困难的。为 了保证边界元法求解的控制方程为常系数偏微分方程，一般还要求求解区域 是均匀介质的。因此，边界元法比较难于求解控制微分方程为非线性的问题 和含有非均匀介质的问题；虽然边界元法可以把问题的维数减少一维，但所 得的影响系数矩阵是一个非对称、非带状和非稀疏的满秩矩阵。特别是在求 解非线性问题时，不可避免地包含非线性项的域积分，因而就必须在域内划 分单元或网格用以计算非线性项的域积分，而这些积分的精度在非线性问题 的收敛性中起支配作用，因此域内单元或网格的划分就不能太粗糙。这样边 界元法可以把问题的维数减少一维的优越性就不显著了。 有限元、边界元等方法是求解边值问题强有力的数值方法。这类方法 都是以单元为基础的，因此存在个共同的缺点是，每次计算前都要剖分阚 格，数据准备工作量大，尤其是对三维问题。当这类方法用于自适应计算或 模拟裂纹扩展时，一般要不断更新网格( r e m e s h i n g ) 。虽然目前已有一些网格生 成器，但人们还是觉得准备数据占用机时多，不方便】。因此，人们呼唤一 种不划网格的数值方法，即无网格法( m e s h t e s sm e t h o do rg r i d l e s sm e t h o d l 。 1 2 无网格法的研究现状 无网格法就是采用对所考虑问题域内随机分布点的变量的值( 或名义节点 值) 的局部插值函数作为试函数，来满足数值求解的局部性要求。无网格方法 具有灵活，容易实施数值计算，求解精度高，在离散模型中不需要划分单元( 在 边界上) 和网格( 在域内) ，在未知变量急剧变化的地方，只需增加节点，对求 4 - 陈竿早：弹塑肚力学问题的无单元伽辽金法( e f g m l 解复杂边界问题极具灵活性，特别在工程应用中容易实现智能和自适应算法 等优点，所以近些年来，这种方法已被广为推广和应用u - 4 s l 。 追溯起来，无网格数值方法的研究已有2 0 多年的历史。1 9 7 7 年，l u c y 4 1 提出了。种新的数值方法一一光滑粒子水动力学方法( s m o o t h e d p a r t i c l e h y d r o d y n a m i c s ：s p h ) 。s p h 是一种纯拉格朗日方法，不需网格，在天体物理 领域里得到了成功应用。近几年，s w e g l e l 2 1 ，d y k a 1 3 1 等人提出了s p h 方法不 稳定的起因及稳定化方案；j o n h s o n 和b e i s s e l 1 4 1 等人提出了一些改善应变计算 的方法；l i u l i5 j 等人提出了对核函数的修正方案。 另外一条构造无网格方法的途径是采用移动最小二乘法( m o v i n 2l e a s t s q u a r e s m e t h o d ，简记为m l s ) 进行近似。m l s 最早由l a n c a s t e rp 等【1 6 1 提出， 用于构造插值函数来拟合曲线和曲面，n a y r o l e s 等【5 1 在研究有限元法的过程 中，提出使用一种被称为“弥散单元”( d i f f u s ee l e m e n t ) 的新的单元类型。事 实上，他们是将移动最小二乘法近似用于g a l e r k i n 方法中，并将之称为弥散 单元法( d i f f u s ee l e m e n tm e t h o d s ，简称d e m ) 。在这种方法中位移函数的形成和 区域积分的实现都可以脱离单元的概念。b e l y t s c h k o 等( 1 9 9 4 ) x 寸1 f 1 2 做了进一步 改进，这些改进包括 1 2 3 1 9 1 ： ( 1 ) 对形函数导数考虑得更全面；( 2 ) 采用高阶高 斯积分完成区域积分：( 3 ) 引入l a g r a n g e 乘子法施加本质边界条件。这些改 进使得d e m 求解精度更高，更具发展前途【1 , 3 1 。b e l y t s c h k o 等称改进后的d e m 为无单元伽辽金法( e f g m ) ，也有的学者将之简称为无单元法。无单元伽辽金 法与有限元法极为相似，都是基于伽辽金公式，采用局部插值函数作为试函 数来求得近似解。它们之间关键的不同之处在于插值方法、积分方式及本质 边界条件的施加方法。 无单元伽辽金法就是无网格法的一种，但它不是真正的无网格法，因为 积分时仍需背景网格。所幸的是这种网格与节点无关，可以非常粗糙地划分。 鉴于无单元伽辽金法的众多优点，近年来它吸引了大量研究人员的注意，做 了大量深入的研究： 湖南大学帧10 位论文 1 9 9 8 年，周维垣等17 , 1 8 1 对其在平面弹性连续体问题中的应用做了探讨。 他在前人： 作的基础上，做了些改进：( 1 ) 对高斯权函数、样条权函数作了研 究，提日1 了一种新的权函数：( 2 ) 采用罚函数法引入边界条件，对面约束、点 约束以及各种面力的处理都提出了解决方案；( 3 ) 对积分方法、支持域、高斯 点等一些关键问题进行了探讨，并用更丰富的工程实例说明了无单元伽辽金 z 去( e f g m ) 的应用。 周小平、周瑞忠等【20 】对无单元法的插值函数进行了专门的论述，提出使 用动态单元法( d y n a m i c - e l e m e n t m e t h o d ) 来理解无单元法，并着重讨论了移动 最小二乘函数中的矩阵a ( x ) 的相关问题，给出了无单元法比有限元法具有更 高的精度、更高次连续性的直观解释。 庞作会等讨论如何用无单元伽辽金法求解集中力问题，同时在前人工 作 2 , 1 9 的基础上，提出了无单元伽辽金法的一种点积分形式，并给出了相关算 例。但是，点积分实现e f g m 时，积分点多少个合适，没有给出理论上的证 明。从算例来看，点积分的总积分点数与高斯积分的总积分点数大致相当即 可。 由于无单元伽辽金法的近似函数不通过结点变量，即吼( x ，) 占无单元 伽辽金法的一个难点是本质边界条件的引入。在目前的发展状况下，通常使 用的方法有：( 1 ) l a g r a n g e 乘子法 3 2 ”。这种方法是通过l a g r a n g e 乘子引入本 质边界条件，因此不要求近似函数满足本质边界条件。l a g r a n g e 乘子法最大 的缺点是它引入了新的未知量，并且使离散方程的系数矩阵不再具有正定、 带状的特点。但是，l a g r a n g e 乘子法是引入本质边界条件最精确的方法，因 此对于规模较小的问题，这种方法是非常适用的；( 2 ) 修正的变分原理方法 2 6 】。 这种方法是将l a g r a n g e 乘子用相应的物理量代替，这样就可以避免由乘子未 知量产生的不良影响，还可以保持方程的带状特性，从而求解工作量比 l a g r a n g e 乘子法要减小许多：( 3 ) 罚函数法1 82 8 舶】。这种方法具有实施简单， 不引入新的未知量等优点。但是，本质边界条件只能近似地得到满足，罚数口 坠羔兰! 堂型堂塑兰i 型些璺查望! ! 塑坚竺堡! l ：”竺塑! 越大，本质边界条件的满足就越好。实际计算中罚数a 不可能取得无穷大， m 只能取为较大的有限值；( 4 ) 与有限无耦合法【2 7 】。将有限元法和无单元伽辽 金法耦合起来使用的方法具有很大的优点【2 7 2 9 1 。通过与有限元法的鹈合，我 们可以将本质边界条件区域用有限元法近似，这样就可以方便地施加本质边 界条件了。 k r y l 等p3 i 对无单元伽辽金法在薄板弯曲问题中的应用做了探讨，张伟星 等州对其在钢筋混凝土筏板中的应用作了研究，张建辉等【3 5 】又对其在筏板基 础中的应用作了研究，所有算例均表明无单元伽辽金法在解决要求插值函数 c 连续的板弯曲问题是合理可行的，其优势是明显的。 n a g a s h i m a 和o u a t o u a t i 等口”3 2 】对无单元伽辽金法做了些改进，然后将其 推广应用于结构的模态分析，得到了很好的结果。这表明无单元伽辽金法在 动力学分析中的应用也是合理可行的，值得深入研究。 陈建、吴林志等i s 6 对无单元伽辽金法在功能梯度材料( f g m ) 断裂行为方 面的应用进行了一定的尝试。算例计算结果表明，该方法具有较高的计算精 度，其分析f g m 材料断裂行为的有效性、灵活性且不受材料特性参数随坐标 连续变化的影响，可以方便、准确地得到弹性模量梯度变化对应力强度因子 的影响规律。而通常有限元法在分析该类问题时，往往需要大量节点和单元， 处理困难，难以得到满意的结果。可以预计，无单元伽辽金法在f g m 之类的 非均匀材料的力学行为分析方面是非常合适的，还可以得到更为广泛的应用。 刘欣、朱德懋等 3 7 1 针对平面弹性问题发展和推导一种显式后验误差指示 公式，对平面裂纹实例进行了h 型，p 型，h p 型三种不同类型的的无网格自 适应分析。 庞作会、朱岳明等 2 5 1 探求采用无单元伽辽金法求解接触问题。该文采用 的方法是将k a t o n a 界面单元引入无单元伽辽金法。在无单元伽辽金法中引入 界面单元，总的原则与有限元( f e m ) 中引) , # n 单n - - r 4 ，具体过程有些不同。 即t 根据点对状态修改总刚及荷载矩阵时，点对支持域内其它节点上的数据 湖南大学坝士学位论文 也耍修改。造成这些不同的原因z l ：了。无单元伽辽金法是通过移动最小二乘法 构造出位移函数。这是无单元伽辽金法在求解非线性问题应用中的尝试，但 是无单元伽辽金法在求解材料非线性、几何非线性问题的应用都尚未见报道。 a t l u r i 等人 4 4 , 4 6 】提出了另外- f e e 无网格法无网格局部p e t r o v g a l e r k i n 方法。这种无网格法在积分时也不再需要背景网格，是- f e e 真正的无网格方 法。a t l u r i 等人首先将这种方法应用于求解调和算子的拉普拉斯方程和泊松方 程。龙述尧【4 5 l 又将这种无网格法推广应用于求解弹性力学平面问题。这种方 法采用移动最小二乘近似函数作为试函数，并且采用移动最小二乘近似函数 的权函数作为加权残值法的加权函数：同时这种方法只包含中心在所考虑点 处的规则局部区域及其边界上的积分，所得系统矩阵是一个带状稀疏矩阵。 算例结果证明：该方法具有收敛快，精度高等优点。 但是，无网格局部p e t r o v g a l e r k i n 方法的积分是在包含中心在所考虑点 处半径为矗的圆域q 。及其边界上进行的，这种圆域q 。就相当于无单元伽辽 金法中的背景网格，而且影响求解结果的局部域半径n 的选取又没有具体的 公式可遵循。从理论上讲，只有所有子域并集覆盖了整体域q ，即u q 3 d ， 求解的结果才能满足整体域及其边界上的平衡方程和边界条件。这样，计算 量显然比无单元伽辽金法要大，计算效率比无单元伽辽金法要低。 我们可以看出，无单元伽辽金法与无网格局部p e t r o v g a l e r k i n 方法都是 与有限元相似的数值方法。无网格局部边界积分方程方法日8 。4 1 则是从边界元 的基础上发展起来的。z h u 和z h a n g 等 3 8 1 首先将无网格局部边界积分方程方 法应用于求解调和算子的拉普拉斯方程和泊松方程，z h u 和z h a n g 等【4 0 l 用其 求解位势非线性问题，龙述尧等 4 3 1 又将该方法推广应用于求解弹性力学平面 问题，s l a d e k 和a t l u r i 等1 4 2 】进一步用该方法求解各向异性材料的弹性力学平 面问题。这种方法采用局部边界积分方程来表示所考虑点的未知函数值，并 包含该点的支持域内其它点的值。在局部边界积分方程中，用于近似函数的 试函数的连续性要求可以大大降低，在构造系统刚度矩阵时，通常不需要形 陈苹莘：弹塑性力学0 u 题的无单元伽辽金法( e f g m l 函数的导数；当采用非插值型的移动最小二乘函数用作试函数h j ，本质边界 条件也可以直接施加进去，这比无单元伽辽金法及无网格局部p e l t o v - g a l c l k i n 方法都有优势。而且，这种方法积分时也不再需要背景网格，同样是一种1 真 萨的无网格法。 但是，无网格局部边界积分方程方法也有其固有的缺点。在用无网格局 部边界积分方程方法求解由微分算子所控制的边值问题时，除了需要所求问 题微分算子在无限空间中的基本解外，还需要所求问题微分算子在球域上( 对 三维问题) 或圆域上( 对二维问题) 的“友解”( c o m p a n i o ns o l u t i o n ) 。但是 有些问题的友解不容易求得，这就给无网格局部边界积分方程方法的应用带 来了局限性。无网格局部边界积分方程方法的积分是在包含中心在所考虑点 处半径为的圆周崛上( 有些问题包含域积分) ，这种圆域q 。就相当于无单 元伽辽金法的背景网格，而且影响求解结果的局部域半径r o 的选取同样没有 具体的公式可遵循。另一个令人感到棘手的是奇异性问题。虽然a t l u r i 和 s l a d e k 等口“4 1 1 对奇异性问题进行了讨论，但所得的计算公式不够简便。 此外，参照无单元伽辽金法的理论形成了一种被称为“流形元法”的新 的数值方法4 7 删1 ，这也是一种无网格方法，有些研究人员正在致力于研究。 1 3 课题的研究意义及论文的研究内容 随着计算机技术的迅速发展，求解力学问题的数值解法得到了不断的发 展。从有限元法到边界元法，再到无网格法，各有各的优点，但各有各自不 可克服的缺点。所以对各种数值解法的研究一直不间断，各种数值解法也得 到不断的完善。 虽然有限元，边界元都已经发展的十分完善，而且一般形式的无网格方 法在目前的发展状况下还不如有限元计算效率高( 当然无单元伽辽金法也不例 外) ，但是无单元伽辽金法在场函数近似和对局部特性的描述、不需划分单元、 后处理方便、可消除体积闭锁现象等诸多方面都具有有限元法不可比拟的优 q 湖南人学硕l 学f 十论义 点。无单元伽辽金法的这些优点使得它在处理需不断增加h 点、不断更新嘲 格的问题时独具魅力。因此，无单元伽l 、辽金法必将对计算力学和结构分午i 生具大的影响，无单元伽辽金法的应用研究具有重大的实用价值。 另外，我们还可以考虑将有限元法、边界元法以及无单元伽辽金法等4 ； 同数值方法耦合起来，从而形成功能更加强大通用程序包。这样就能充分发 挥各种数值方法优势互补的作用，能够根据具体问题的需要自动选择最适合 的数值方法来计算。因此，对无单元伽辽金法的研究也将对完善数值计算通 用软件( 如a n s y s ) 的功能起到重大作用。 目前，国内外学者对无单元伽辽金法做了大量的研究工作，但是用其来 求解弹塑性问题的研究还未见报道。本论文就是将无单元伽辽金法( e f g m ) 推 广应用于求解弹塑性问题，最后用算例验证了本文方法的可行性。 本文工作主要有： ( 1 ) 在阅读大量文献的基础上，系统地综述了无网格方法的发展现状及各 种无网格方法的优缺点。 ( 2 ) 采用移动最小二乘法来构造位移函数，是目前各种无网格方法最流行 采用的方式。本文对移动最小二乘法的物理意义、爿( 工) 的可逆性、权函数的 选取、支持域半径的确定以及插值函数的性质等关键性问题都进行了讨论。 ( 3 ) 在对有限元求解弹塑性问题的理论进行分析的基础上，结合无单元伽 辽金法的特点，提出了用无单元伽辽金法求解弹塑性问题的原理，并对一些 关键性问题提出了解决的方案。 ( 4 ) 编制了无单元伽辽金法求解弹塑性问题的程序。最后，利用已编的程 序做了大量的算例分析，并将所得结果与a n s y s 的计算结果进行了l l 较5 7 析，从而验证了程序的正确性和本文方法的可行性。 ( 5 ) 在大量试算的基础上，对无单元伽辽金法求解弹塑性问题的特点进行 了讨论，得出了一些结论。 2 1 引言 陈莘莘弹塑性力学问题的无t 一九伽辽金7 ，( e f g m ) 第二章移动最小二乘法 无网格方法中位移函数的形成和区域积分的实现都可以脱离单元，是因 为采用了与有限元、边暴元等数值方法不同的插值函数。无网格方法的插值 函数一般有三种：( 1 ) 核函数近似方法 5 1 1 ：( 2 ) 移动最小二乘法( m l s ) ；( 3 ) 单位 分解法p “。目前，各种无网格法最流行采用的是移动最小二乘法f m l s l 。 移动最小二乘法最早用于构造插值函数来拟合曲线和曲面。移动最小 二乘法是用随机分布的节点上未知变量的值来表示试函数，且试函数与有限 元法的插值函数样，具有局部性质。 本章首先以标量函数“为例说明了移动最小二乘近似法的原理，然后对移 动最小二乘法的物理意义、爿( x ) 的可逆性、权函数的选取、支持域半径的确 定以及插值函数的性质等关键性问题都进行了讨论。 2 2 插值公式的推导 考虑点x 的邻域( 子域) 力。，它位于全域力内，为了近似函数“在子域力一 的分布，在有限个随机分布的节点x ，i = 1 。2 一， 上函数“的近似式 “( x ) ，v x 力：，可以定义为 “1 x ) = p r ( x ) 口( x ) ，v x 力。，x = i xy z 7( 2 1 ) 式中p 7 ( x ) = p ( x ) p ：( x )p 。( x ) 】是m 次完备单项式的基，( 工) 是包含 系数日，( x ) ，= l ，2 ，脚的向量，这些系数是空间坐标x = i xy ：】7 的函数。 例如，对于二维问题 湖南火学础i 七彳位论文 线性基p 7 ( x ) = 1 x y 】，= 3 二次基p v ( x ) = 1 x y x 2 x yy2 ，= 6 三次基p 1 ( x ) = 【1xyx 2x yy 2 x 3 x 2 yx y 2y 3 ，m = 1 0 ( 22 a ) ( 22 b ) f 22 c ) 定义加权离散l 模为 ，( x ) = w ，( x ) e p 7 ( x ，) 口( ) 一二，】2 = p 口( x ) 一二 7 w pn ( x ) 一二( 2 3 ) f = l 式中w ，( x ) 是节点f 的权函数，对于在w ，( j ) 的支持域内的所有x ，有w ，( 工) 0 ， x ，表示节点i 的x 值，h 是域 2 ，内权函数w ，( x ) 0 的节点数，而矩阵p 和w 分别为 且 p = p 1 ( x ，) p 1 ( x 2 ) p 7 ( _ ) j = 降 p l ( x 1 ) p l ( x 2 ) p ! ( x 1 ) p 2 ( x 2 ) p l ( x 。) p 2 ( x 。) 0 矗7 = 眩如二。1 p 。( x i ) p 。( x 2 ) f 24 、 ( 2 5 1 ( 2 6 1 在这里应该特别注意：在式( 23 ) 和( 2 6 ) 中的，i = 1 , 2 ，n 是名义节点值，通常 它不是未知试函数“( x ) 的节点值( 见图2 一1 ) 当o ) 取极小值时，由孚望：o 得到 d ( 工) a ( x ) a ( x ) = b ( x ) f i 式中矩阵爿( x ) 和曰( x ) 分别为 ( 27 ) 陈莘卑弹塑性力学问越的尢单元伽l 辽金法( e i ? ( ，m ) 4 ( x ) = p t w p = b ( x ) p = 1 0 ( x ) ，( 卫，) p 1 ( x ，) ( 28 ) 嚣( x ) = p w = w i ( x ) p ( x 】) w 2 ( x ) p ( x2 ) b o u n d a r yn o d e 图2 1 试函数u “( x ) 中的节点“，和名义节点值五的区别 只有在方程( 2 7 ) 中的矩阵4 是非奇异时，才可定义移动最小二乘近似函 数。而要4 非奇异，必须使得p 的秩等于m 。因此定义移动最小二乘近似函 数的必要条件是，对于每个样点x 6 2 至少有m 个权函数不为零f 即”m ) ， 而且在力，中的节点不能有特殊形式的分布，例如分布在一条直线上。这里所 说的样点既可以是所考虑的节点也可以是高斯积分点。 由式( 2 7 ) 解出口( 工) 并代入式( 21 ) 得 式中 “6 ( x ) = 西1 ( x ) f i = 妒，( 工) 西， “6 ( 一) 驯，五，工力。( 2 1 0 ) j e l 垂7 ( x ) = p r ( 工) 4 一| ( 工) b ( x ) ( 21 1 ) 型l ! 垒! 型! 土兰些堡兰 妒，( x ) = p ，( x ) 爿。1 ( 戈) 曰( 工) 】，( 2 1 2 ) = l 蛾( x ) 是类似于有限元法的插值形函数，我们称之为节点x 的移动最小二乘近 似的形函数。由式( 2 9 ) 及( 2 1 2 ) 可以看出，当w ，( x ) = 0 时，纪( x ) ：0 。在实际 应用中，通常选择w ，( x ) 在节点x 的支持域上不为零，而节点x 的支持域通常 取中心在x 点，半径为l 的圆。因此，对于不在节点x 的支持域内的x ，有 够( x ) = 0 ，这使得移动最小二乘近似法具有如有限元法一样的局部性质。 移动最小二乘近似函数“并不对节点进行插值，因为“一( x ，) ；“，且 妒，( 工，) 吒。形函数移( x ) 的光滑性由基函数和权函数的光滑性所确定。令 c k ( 力) 是阶连续可微函数空间，如果w ，( x ) c ( 力) ，p ，( x ) c l ( 力) ，其中， i = 1 ，2 ，+ ，n ：= l ，2 ，一，m ，则妒，( 戈) c ( 力) ，r = m i n ( k ，) 。 在计算时，我们需要用到形函数p ，( x ) 的偏导数，其表达式如下f 3 】 纪女= z i p ，女( 爿一1 曰) ，+ p ，( 爿b ，+ 爿j 1 b ) 。( 2 13 ) j 2 l 式中爿j 1 = ( 4 “) ，。表示爿的逆矩阵对x 的偏导数，给出为 爿j 12 4 “爿女a r 2 1 4 ) 式中( ) ，表示偏导a ( ) 玉7 。 2 3 插值函数的物理意义 根据式( 2 1 2 ) ，可以很容易地将插值函数分成相互独立的3 个部分，即 p ( x ) , a - i ( x ) 和占( x ) 。其中，p ( x ) 包含了插值点的位置信息；a t ( x ) 是由所 有插值基点的坐标进行包括加权在内的简单运算后叠加得到的。由于采用了 以距离为单变量的权函数，使得与某一插值点有关的插值基点被限制在以 这一插值点为圆心的圆形区域，即支持域内，故4 1 ) 是由插值点x 的支持域 内所有插值基点的坐标形成的；口似) 中的每列可看成是表示对应的插值苯 点与插值点x 之间的分配关系。因此，对于无单元伽辽金法的插值函数，可 一 堕羔兰! 翌望型生兰塑望塑垄望兰塑坚垒堕! 兰里! l 一以这杼理解，即对于插值点x ，先由其支持域内的所有插值基点形成代表这 一支持域整体的，4 ( x ) 并求逆，再利用召( ) 将这一整体信息分配到域内的各插 值基点，就可得到插值点与这些基点问的插值函数值。不难理解，虽然无单 元伽辽金法不需要单元信息，但在计算爿“) 时已在无形中构造了一种虚拟的 单元。这种单元随着插值点的移动而移动，并随着影响半径的变化而伸缩， 且在互相邻近的插值点上出现( 如图2 - - 2 ) 。因此，将无单元伽辽金法称为动 态单元法( d y n a m i ce l e m e n tm e t h o d ) 似乎更为贴切。 支持域 图2 - 2 支持域示意图 插值基点 插值点 分析有限元法的插值方法，单元内任一点处的场函数值是由该单元所有 的结点插值得到的。以三角形单元为例，它的插值基点只有3 个。而对于无 单元伽辽金法，即使采用与有限元法相同的结点布置，对同一点来说，只要 取足够大的影响半径，则有效的插值基点数将是任意的。因此，采用数量更 多的、分布更广的插值基点，并使用县有高阶连续性的插值函数，这就是无 单元伽辽金法比有限元法具有更高的精度、更高次的连续性的直观解释。 2 4 a ( x 1 的n - - i 逆性 由式( 2 1 1 ) 和式( 2 1 3 ) 可以看出：每次求形函数及其导数时，都涉及到求 a ( x ) 的逆，因此如果a ( x ) 不可逆，将是整个计算无法进行下去。而要a ( x ) 可 逆，必须有：对于每个样点x q 至少有m 个权函数不为零( 即n m ) ，即每 个样点工力的支持域内至少有m 个插值基点。现给出证明如下： 1 5 - 塑塑叁兰：塑土堂垡篓兰 由式f 24 ) 可得：p 是一个月m 的实矩阵；e t i a ( 25 ) 可得：w 是一个n n 的实矩阵；由式( 28 ) 可得：a ( x ) 是一个m 阶方5 1 ：。 要爿( x ) 非奇异，即爿( x