1、边界元法之介绍,姚振汉（北京清华大学工程力学系）,引言由偏微分方程建立边界积分方程求解边界积分方程的边界元法位势问题边界元法的应用例题边界元法的各种应用快速多极边界元法的大规模计算边界型无网格法简介进一步研究方向,内 容 提 要,边界元法又称边界积分方程边界元法是在有限元法 (1955，56) 之后发展起来的一种精确高效的工程与科学问题的数值分析方法。它以边界积分方程为数学基础，借鉴了有限元法的分元离散技术，形成了精确高效的数值分析方法。经过近50年的发展，它在一些领域已经成为有限元法最重要的补充。,引 言,边界元法间接法位势问题 (Smith & Pierce, 1958)弹性力学 (Mas
2、sonet，1965)边界元法直接法位势问题 (Jaswon，1963)弹性力学 (Rizzo，1967)本研究组1979年开始有关边界元法的研究，至今没有间断。,引 言,边界元法的主要优点：降低维数 (便于模拟复杂几何构型)；高精度；适于处理高梯度，甚至有奇异性问题；适于处理无限域、半无限域问题。常规边界元法的主要缺点：求解方程组具有非对称满阵，解题规模受限制；对一般非线性问题缺少高效计算方案。,引 言,近年的新进展：一方面对于原来有限元能算、边界元不能算的问题，发展了边界元法的计算方案，并验证其计算结果和有限元法吻合很好；针对不能计算大规模问题的弱点，采用快速多极算法加速，发展了快速多极边
3、界元法，已能处理甚至比有限元法处理的规模更大的问题；结合近年无网格法的发展，提出了多种边界型的无网格法，例如边界点法，杂交边界点法等。,引 言,以位势问题为例：在物理问题中有不少问题可以归结为位势问题，即控制方程为调和方程，又称Laplace方程。比如最常见的有温度场问题，电场问题，理想流体流动的流场问题，还有用翘曲函数求解的柱体扭转问题等。 在数学上，以二维问题为例，调和方程可写成它的具有二阶连续偏导数的解称为调和函数。,由偏微分方程建立边界积分方程,由偏微分方程建立边界积分方程,内部问题和外部问题,由偏微分方程建立边界积分方程,Dirichlet问题和Neumann问题,此外还有Robin
4、问题，和更一般的混合边值问题。,由偏微分方程建立边界积分方程,Green等式,设函数 u, v 及其所有一阶偏导数在闭域 上连续，并在 W 内有所有的二阶偏导数，则成立,Green第一等式,Green第二等式,由偏微分方程建立边界积分方程,Laplace方程的基本解和解的积分表达式,基本解即无限域上的单位奇异解，满足方程：,对于二维问题，该基本解为,代入Green第二等式，可得,由偏微分方程建立边界积分方程,边界积分方程的建立,在解的积分表达式中令奇异点 P 从域内趋于边界点 p,其中,这就是位势问题的边界积分方程,由偏微分方程建立边界积分方程,边界积分方程,对最一般情况，将边界分成给定函数值
5、部分和给定法向导数值部分，边界积分方程可写成,求解边界积分方程的边界元法,边界的离散和边界变量的离散,边界划分为单元,边界变量在每个单元上用插值函数近似,按插值函数不同可分：常值元，线性元，二次元等,求解边界积分方程的边界元法,边界积分方程离散化为线性代数方程组,将插值近似的边界变量代入边界积分方程，可得,近似解不可能要求对任意 p 点均满足此方程，而只能要求在加权余量法的弱形式下使方程得到满足。,求解边界积分方程的边界元法,边界积分方程离散化为线性代数方程组,通常采用加权余量法的配点格式，得到,这实际上就是对于各节点未知量的线性代数方程组。,求解边界积分方程的边界元法,线性代数方程组的矩阵形
6、式,上述线性代数方程组可以写出如下矩阵形式,U 和 Q 分别为节点函数值和法向导数值列阵，而 H和 G 两个矩阵的元素都是核函数（即基本解）和插值形函数乘积在单元上的积分,求解边界积分方程的边界元法,线性代数方程组的矩阵形式,将节点函数值和法向导数值列阵按照给定量和未知量重新排列，得到,将给定量和相应系数阵直接相乘得到常数列阵，最终得到边界元法的求解代数方程组,求解边界积分方程的边界元法,计算核函数和形函数乘积的积分是关键之一,计算量大，而且积分精度对结果的精度有直接影响。,当奇异点 p 落在积分单元上时，积分有奇异性。对于各类积分采用不同的处理办法，例如：非奇异积分可用等精度高斯积分公式；弱
7、奇异积分通过变量置换消除奇异性，然后积分；强奇异的Cauchy主值积分，可利用简单特解代入方程，然后由其它系数间接确定。,求解边界积分方程的边界元法,线性代数方程组的求解,小规模问题用高斯消去法：存储量: O(N 2)计算量: O(N 3),规模较大问题可用迭代解法：存储量: O (N 2)计算量: O (k N 2),计算量大，解题规模受到限制。,位势问题边界元法的应用例题,渗流问题 (成功大学土木系吴教授，1995),位势问题边界元法的应用例题,渗流问题 (成功大学土木系吴教授，1995),位势问题边界元法的应用例题,漫流问题 (弗吉尼亚州立大学Wu P-H，1996),位势问题边界元法的
8、应用例题,漫流问题 (弗吉尼亚州立大学Wu P-H，1996),位势问题边界元法的应用例题,铸模传热优化问题 (Davey K，2003),位势问题边界元法的应用例题,铸模传热优化问题 (Davey K，2003),位势问题边界元法的应用例题,埋管阴极保护 (Purcar M，2003),位势问题边界元法的应用例题,埋管阴极保护 (Purcar M，2003),位势问题边界元法的应用例题,多孔板传热 (Liu YJ等，2005),位势问题边界元法的应用例题,多孔板传热 (Liu YJ等，2005),边界元法的各种应用,二维弹性力学问题模拟光弹性实验 (本组, 1998),边界元法的各种应用,三
9、维弹性应力分析 (本组, 2000),滚筒轴边界元模型：5个子域，8184结点，1721单元,边界元法的各种应用,三维弹性应力分析 (本组, 2000),滚筒轴有限元模型：13000结点，64000单元,边界元法的各种应用,三维弹性应力分析 (本组, 2000),边界元解：左端局部Mises应力(单位：MPa),有限元解：左端局部Mises应力 (MPa),2000年得到的启示：,快速多极边界元法的大规模计算,快速多极算法发展概况 传统的求解 N particle 问题算法的计算量： O(N 2)静电势静电力 Barnes, Hut (1986) ： 树结构算法 (tree code) 计算量
10、：O (NlogN)J. Barnes, P. Hut. A hierarchical O (NlogN) force calculation algorithm. Nature, 324, (1986), 446-449.,快速多极边界元法的大规模计算,快速多极算法发展概况 Rokhlin, Greengard (1987)： 快速多极算法 (Fast multipole method, FMM)： O(N) 采用了： 多极展开 (multipole expansion) 局部展开 (local expansion) Carrier等 (1988) 自适应快速多极算法 (Adaptive F
11、MM)： O(N),快速多极边界元法的大规模计算,快速多极算法发展概况 Greengard, Rokhlin (1997): 3D新型快速多极算法 (New Version FMM) O(N) 引入 指数展开 (Exponential expansion) Hrycak, Rokhlin (1998)： 2D新型快速多极算法 O(N),快速多极边界元法的大规模计算,二维含大量随机分布圆形夹杂的复合材料大规模计算,快速多极边界元法的大规模计算,三维含大量随机分布球形颗粒的复合材料大规模计算沿X轴单向拉伸时颗粒/基体界面法向应力分布图,快速多极边界元法的大规模计算,三维含大量随机分布和随机方向短纤
12、维的复合材料大规模计算,快速多极边界元法的大规模计算,并行计算的计算规模,最大三维算例 (4000夹杂)：5,064,000自由度，44.775小时。,最大二维算例 (4000夹杂)：8,080,000自由度，8.29小时。,快速多极边界元法的大规模计算,界面法向面力分布（隐去基体显示）,快速多极边界元法的大规模计算,快速多极边界元法的大规模计算,用二次单元的高精度大规模计算,二维无限弹性体中3000裂纹 (900,000 DOF, p = 20, 6.15 h),快速多极边界元法的大规模计算,含大量裂纹弹性体的等效材料特性微裂纹非均匀随机分布的二维弹性体模型,快速多极边界元法的大规模计算,含
13、400裂纹交变受拉长方板疲劳裂纹扩展模拟,快速多极边界元法的大规模计算,含400裂纹交变受拉长方板疲劳裂纹扩展模拟,快速多极边界元法的大规模计算,含400裂纹交变受拉长方板疲劳裂纹扩展模拟,快速多极边界元法的大规模计算,低雷诺数流体微喷的数值模拟,快速多极边界元法的大规模计算,低雷诺数流体微喷的数值模拟 (边界元网格),快速多极边界元法的大规模计算,低雷诺数流体微喷的数值模拟 (面力分布计算结果),快速多极边界元法的大规模计算,低雷诺数流体微喷的数值模拟 (速度分布计算结果),快速多极边界元法的大规模计算,碳纳米管复合材料大规模模拟,快速多极边界元法的大规模计算,碳纳米管复合材料大规模模拟,最
14、大模型：16,000 CNT28,800,000 DOF计算时间 34 小时,快速多极边界元法的大规模计算,碳纳米管复合材料大规模模拟,快速多极边界元法的大规模计算,MEMS (AutoMEMS, Coyote System Inc. ),快速多极边界元法的大规模计算,超高频电磁波散射大规模分析,快速多极边界元法的大规模计算,超高频电磁波散射大规模分析,快速多极边界元法的大规模计算,无网格局部边界积分方程法 (MLBIE) (Zhu and Atluri et al., 1998) 无网格局部彼得洛夫-伽辽金法 (MLPG) (Atluri et al., 1998; Kim and Atlu
15、ri, 2000; Lin and Atluri, 2000) 两种方法均利用局部子域上的弱形式，采用移动最小二乘插值，得到了“真正的无网格法”。,边界型无网格法简介,边界点法 (BNM) (Mukherjee et al., 1997) 将移动最小二乘插值(MLS)和边界积分方程结合 (BIE) ，无网格，且降维。仅需要物体边界上的结点数据结构， 但在积分时还需要背景网格。,边界型无网格法简介,杂交边界点法将移动最小二乘插值和修正变分原理结合，包含三类互相独立的变量，以位势问题为例，即: 边界上的势边界上的法向流量域内的势,边界型无网格法简介,域内变量用基本解来插值，从而可以将弱形式中的域内
16、积分化为边界积分。边界变量采用移动最小二乘插值。 象边界元法那样降低维数， 又是真正的无网格法，只需要建立在边界上分布的点的信息。,边界型无网格法简介,用FM-HBNM求解大规模三维位势问题,620,708 节点，53,840 秒,含复杂形状空腔,边界型无网格法简介,碳纳米管复合材料的一种代表性胞元,边界型无网格法简介,得到的等效导热系数和碳纳米管(CNT)长度的关系,边界型无网格法简介,碳纳米管复合材料的另一种代表性胞元,边界型无网格法简介,得到的等效导热系数和碳纳米管波纹度的关系,边界型无网格法简介,边界元法和各种无网格法等在许多工程领域都属于有限元法的一个补充。它的发展既受到有限元法的启发，又受到有限元法的制约。只有在其优势领域，作为补充才有实际意义。如果只是有限元能计算的它也能算，而且结果吻合，那就证明这种补充没有必要。,进一步研究方向,快速多极边界元法已经证明可以有效地进行有限元法都难以完成的大规模计算。当务之急是为它寻找工程上有重要意义的实际应用。由于线性问题的有限元方法已经非常成熟，而边界型数值方法在非线性分析方面的研究还很少，因此快速多极边界元法在非线性分析中的应用值得进行研究。我们在弹塑性分析快速多极边界元法方面已经初步取得成功。,进一步研究方向,谢谢诸位！,
