1、任何模拟方法，都必须在最佳计算速度和数值精度之间寻找平衡点。,要在各种可能的求解方法中找到一种统一地适用于计算材料学领域（或其它领域）的理想方法，一般是不现实的。,由于实际问题的具体特征、复杂性以及算法自身的适用范围决定了应用中必须选择、设计适合于自己特定问题的算法，因而掌握数值方法的思想至关重要。,科学计算（数值模拟）已经被公认为与理论分析、实验分析并列的科学研究三大基本手段之一，但三者之间相辅相成。,偏微分方程的有限元方法,一 边值问题的变分原理,1 引论,求使得泛函,达到最大的函数 。,（1）等周问题,在长度一定的所有平面封闭曲线中，求所围面积为最大的曲线。,定义：当求泛函在一个函数集合
2、K中的极小（或极大）问题，则该问题称为变分问题。,变分问题与微分方程的定解问题有一定的联系。,（2）初等变分原理, 一元二次函数的变分原理,考察J(x)的极值情况。,变分原理：,设,求 ，使,与求解方程 Lx = f 等价。,对称正定, 多元二次函数的变分原理,求J(x)取极小值的驻点， 其中,设,设,则J(x)可表示为：,变分原理：,设矩阵A对称正定，则下列两个命题等价：,(b) 是方程 的解,上述两个例子表明：,其中,求二次函数的极小值问题和求线性代数方程（组）的解是等价的。,求弦的平衡位置归结为求解两点边值问题：,设弦处于某一位置u=u(x)，可得到其总位能为,极小位能原理：,其中T是弦
3、的张力。,平衡原理,弦的平衡位置 (记为 )将在满足边值条件 u(0)=0，u(l)=0 的一切可能位置中，使位能取极小值。,弦的平衡位置 是下列变分问题的解,在数学上，要将某个微分方程的定解问题转化为一个变分问题求解，必须针对已给的定解问题构造一个相应的泛函，并证明定解问题的解与泛函极值问题的解等价。,有限元方法正是利用这种等价性（边值问题与变分问题的等价性），先将微分方程定解问题转化为变分问题（或变分方程）的求解问题，然后再设法近似求解变分问题（或变分方程）。,（2）两点边值问题的变分原理, 构造泛函,考察二阶常微分方程边值问题：,引入泛函算子,则, 变分问题,与前述二阶常微分方程边值问题
4、相应的变分问题是,其中,求 ，使, 变分原理（变分问题与边值问题的等价性）,设 ， 是边值问题,的解，则 使 J(u) 达到极小值；,反之，若 使 J(u) 达到极小值，则 是边值问题的解。,其中,是强制边界条件， 是自然边界条件，区别这两类边界条件在用有限元方法求解边值问题时很重要。,（3）虚功原理,对两点边值问题：,其中,虚功原理,，且满足变分方程：,设 ，以v乘方程两端，沿a,b积分，并利用 ，得变分方程,对任意,在力学里， 表示虚功,设 ，则 是边值问题解的充要条件是：,对于复杂的边界条件，边值问题的求解一般是困难的。若将微分方程化为相应的变分问题或变分方程，则只需处理强加边界条件，无
5、需处理自然边界条件（自然边界条件已包含于变分问题中泛函的构造或已包含于给出的变分方程之中）。这一特点对研究微分方程离散化方法及其数值解带来了极大的方便。,3 二阶椭圆边值问题的变分原理,（1）极小位能原理,模型方程,其中G是平面有界区域。, 构造泛函,引入泛函算子,则, 变分问题,与前述二阶椭圆边值问题相应的变分问题是,求 ，使,其中, 变分原理（变分问题与边值问题的等价性）,对第一边值问题，无论齐次或非齐次边界条件，泛函是一样的，只是边界条件要作为强加边值条件加在所取的函数类上。,设 ， 是二阶椭圆边值问题的解，则 使 J(u) 达到极小值；,反之，若 使 J(u) 达到极小值，则 是二阶椭
6、圆边值问题的解。,其中,对第二、三类边值问题，无论齐次或非齐次边界条件，二次泛函形式相对于第一边值问题有所改变，但函数类的选取与边界条件无关。,（2）虚功原理,问题,其中,设 ，以v乘方程两端后在G上积分，并利用Green公式，得变分方程,虚功原理,在力学里， 表示虚功,设 是边值问题的解，则对任意， 满足变分方程。,反之，若 ，且对任意满足变分方程，则 为边值问题的解。,与极小位能原理类似，第一类边界条件为强加边界条件，第二、三类边界条件为自然边界条件。,虚功原理比极小位能原理应用更广。,目的：求解相应的变分问题或相应的变分方程。,Ritz方法是近似求解变分问题(即二次泛函极小值)的算法。G
7、alerkin方法是近似求解变分方程的算法，这两种算法统称为Ritz-Galerkin方法。,Ritz-Galerkin方法的基本思想,以下用V表示 等Sobolev空间，L表示微分算子，(u,v)为由L及边值条件决定的双线性泛函。,4 Ritz-Galerkin方法,用有限维空间的函数代替变分问题(或变分方程)中无限维空间的函数，从而在有限维函数空间中求变分问题(或变分方程)的近似解，并要求当有限维空间的维数不断增加时，有限维近似解逼近原变分问题(或变分方程)的解。,由极小位能原理得出的变分问题为：,Ritz方法：求变分问题的近似解。,（1）Ritz方法,求 ，使,其中 ，,设 是V 的n维
8、子空间， 是 的一组基底(称为基函数) 。 中任一元素 可表示为,即选择适当的 ，使 取极小值。,求 ，使,Ritz方法：,展开,令,则 满足,解出 代入 ，则得,Ritz方法步骤为：,根据最小位能原理构造相应于微分方程或物 理问题的变分问题；,取 作为 的一组基底，即用 近似代替无穷维空间V；,求解关于 的线性代数方程组。,由虚功原理得出的变分方程为：,Galerkin方法：求变分方程的近似解。,（2） Galerkin方法,设 是V 的n维子空间， 是 的一组基底(称为基函数) 。 中任一元素 可表示为,即选择适当的 ，使 取极小值。,Galerkin方法：,求 ，使对 ， 满足,由 的任
9、意性，取 作为v ，则得,将 代入变分方程，则,解出 代入 ，则得,Galerkin步骤为：,根据虚功原理构造相应于微分方程或物理问 题的变分方程；,取 作为 的一组基底，即用 近似代替无穷维空间V；,求解关于 的线性代数方程组。,取 作为v ，将 代入变分方程，得到 满足的方程组：,有限元法广泛应用的原因,Ritz-Galerkin方法应用的困难, 基函数选取必须满足强加边界条件，因此 选取困难；, 计算量、存储量巨大；, 方程组求解病态严重。,充分发挥了变分形式和Ritz-Galerkin方法的优点；, 摆脱了传统的基函数取法；, 各种问题的结构程序格式统一。,有限元方法基于变分原理， 又
10、具有差分方法的一些特点，并且适于较复杂的区域和不同粗细的网格。,二 椭圆型方程的有限元方法,差分法解偏微分方程，解得的结果就是准确解u在节点上的近似值；,Ritz-Galerkin方法得到近似的解析解，但对一般区域，却往往难以实现。,有限元方法与传统Ritz-Galerkin方法的差别在于有限维函数空间的构造方法。Ritz-Galerkin方法选用的基函数在整个定解区域上整体光滑，有限元则取分段或分片连续且局部非零的基函数。,考虑两点边值问题：,1 一维问题的线性元,将区间a,b分割为n个子区间 。,第i个单元记为 ，其长度 。,（1）试探函数与试探函数空间,设,则 称为试探函数空间， 称为试
11、探函数。,（2） 用单元形状函数表示试探函数,设在节点上试探函数 在节点上的一组值为,最简单的试探函数空间 由分段线性函数组成。,在第i个单元 上的线性插值函数为,即,当 时， 的（线性）插值公式称为（线性）单元形状函数。,把每个单元形状函数合并起来，就得到整个区间a,b 上都有定义的函数 ：,为使分段插值标准化，通常用仿射变换,显然,把 变到 ，令,则,变为,或,定义基函数系,（3） 用节点基函数表示试探函数,线性无关，它们可组成试探函数空间的基，常称为节点基函数。,几何形状如图,任一试探函数 可表示为,用这类插值型基函数，可以构造出适合各种边界条件的试探函数。,若借助前述放射变换,节点基函
12、数可用变量表示为, 直接形成有限元方程,(a) 把表达式 代入泛函；,（4）从Ritz方法出发形成有限元方程,(b) 将泛函表达式中积分区间a,b变到0,1；,(c) 由 达到极小值的条件,得到含 的有限元方程,这儿,(d) 解出有限元方程的数值解 ，就得到使二次泛函取极小的近似函数(有限元解),有限元方程可用矩阵表示为,其中,称为总刚矩阵。,工程中形成有限元方程时，通常先在每个单元上形成单元矩阵（称为单元刚度矩阵），然后由单元刚度矩阵形成总刚度矩阵（称为总体合成）。, 用单元刚度分析形成有限元方程,(a) 把 按单元组织，则在第i个单元上，令,其中 称为单元刚度矩阵。各元素可计算得到。,再把
13、 扩展成nn矩阵，使其第i1行、第i行和第i1列、第i列交叉位置的元素就是单元刚度矩阵的四个元素，其余全为零（ 只是第一行，第一列元素非零）。即,记,则,其中 称为总刚矩阵。,(b) 由 达到极小值的条件,(c) 解出有限元方程的数值解 ，就得到使二次泛函取极小的近似函数(有限元解),得到有限元方程 。,（5）从Galerkin方法出发形成有限元方程,把表达式 代入变分方程,对前面的两点边值问题，变分方程变为,其中,与Ritz方法相比， Galerkin方法形成的有限元方程其系数矩阵就是总刚矩阵。,该方程即为Galerkin法形成的有限元方程。,由Galerkin方法推导有限元方程更加方便直接
14、，且适用面广。,若希望在每个单元上提高逼近的精确度，则可通过提高插值多项式次数来实现，,在单元 上可构造一、二、三及高次插值多项式，其方法有两种：,2 一维问题的高次元,整个问题计算的全过程除分析单元插值外，均与前面框架类似。, Lagrange型：在单元内部增加一些插值节点。, Hermite型：在节点引进一阶、二阶乃至更高阶导数。, 线性元(Lagrange型),要求：在每一个单元上是一次多项式，在单元节点处连续。 插值条件：在单元的两个端点取指定值。, 二次元 (Lagrange型),要求：在每一个单元上是二次多项式，在单元节点处连续。 插值条件：在单元的两个端点及单元中点取指定值。,
15、三次元（Hermite型）,要求：在每一个单元上是三次多项式，在单元节点处连续。 插值条件：在两个端点取指定的函数值和一阶导数值。,采用高次元，有限元方程形成的方法和线性元类似，但工作量增加。一是计算积分的复杂性增加，二是矩阵的带宽增加。,高次元的主要优点是收敛阶高，且提高了函数逼近的光滑性。,假定区域G可以分割成有限个矩形的和，且每个小矩形（单元）的边和坐标轴平行。,3 二维问题的矩形元,通过仿射变换,采用矩形剖分后，任一个矩形,总可变成单位正方形,如果在 上造出单元形状函数，就可得到试探函数 。而 上的形状函数可通过先在 上造出形状函数，再通过仿射变化而得到。,在 上构造形状函数，也采用L
16、agrange型和Hermite型插值。,Lagrange型：根据若干插值节点处的函数值决定插值函数。,Hermite型：根据若干插值节点处的函数值、一阶偏导数乃至更高阶偏导数决定插值函数。,（1）Lagrange型公式, 双一次插值,插值条件：给定 顶点上的函数值,求：双线性函数,满足,设,令,由 为双线性函数，可求得,令,则,通过仿射变换消去、 ，就得到 上的形状函数。把这些函数按单元叠加，即对所有单元求和，就得到G上的试探函数。,实际计算时，并不消去中间变量、 ，因为计算刚度矩阵元素（定积分）用、作自变量更为方便。,插值条件：给定II上九个插值节点(0,0)、(1/2,0)、(1,0)、
17、 (0,1/2)、(1/2,1/2)、(1,1/2)、(0,1)、(1/2,1)、(1,1)的函数值。,求：双二次函数,满足, 双二次插值,故,通过仿射变换消去、 ，就得到 上的形状函数。,令,由 为二次函数，可求得,设,插值条件：给定II上十六个插值节点(见图) 。,求：双三次函数,满足,设, 双三次插值,故,令,由 为三次函数，可求得,可以在四个顶点分别给定函数值、两个一阶偏导数的值和二阶混合偏导数的值(共十六条件)，确定一个双三次多项式的十六个系数。,（2） Hermite型公式,Lagrange型公式中不出现导数，这样的试探函数只属于 。为了得到属于 的试探函数，需要Hermite型插
18、值公式。,双三次多项式含有十六项：,简单且常用的是不完全的双三次多项式插值。它去掉双三次多项式中的 项。,插值条件：给定II上四个插值节点。,求：不完全双三次函数,满足,四个顶点处 的函数值等于 在该点的函数值； 四个顶点处 的值等于 在该点的值； 四个顶点处 的值等于 在该点的值。,根据仿射变换,则可将原插值问题转化为II上的插值问题。,满足,四个顶点处 的函数值等于 在该点的函数值； 四个顶点处 的值等于 在该点的值乘以x； 四个顶点处 的值等于 在该点的值乘以y。,插值条件：给定II上四个插值节点 (0,0)、(1,0)、(0,1)、 (1,1) 。,求：不完全双三次函数,类似于Lagr
19、ange型公式的构造，可以求得 上的形状函数。,在三角形元的有限元方法中，先将定解区域G化分为若干个小三角形（称作单元）。然后在每个单元上构造插值型函数，并用分片函数（但整体连续的函数）代替变分问题或变分方程中所需求解的函数。,4 二维问题的三角形元,用有限元求解二维椭圆边值问题时，应用最广的是三角形元。,（1）三角剖分,将定解区域化分成若干个小三角形单元 时应注意：, 为了保证有限元解的精确度和收敛性，并避免其离散后代数方程组系数矩阵的病态性，网格剖分中疏密的过渡不要太陡。,错误,为了保证有限元解有较好的 精度，每个单元中应尽量避免出现大的钝角。,应避免, 单元顶点的编号顺序可以任意，但节点
20、编号顺序将影响有限元方程组系数矩阵的结构（带宽）。, 为了方便构造插值型函数，要求每个单元的顶点是相邻单元的顶点。,（2）面积坐标及有关公式,在三角形单元上构造插值型函数，并不简单类同于矩形单元。, 面积坐标,考虑一个面积为S的三角形单元，其顶点按反时针顺序记为i, j, k。在此单元内部任取一点p(x,y)，连接p和三个顶点 ，此单元则被分成三个小三角形它们的面积记为 和 。,记,单元内任一点p(x,y)的位置与三维数 一一对应，称 为面积坐标。,面积坐标和直角坐标之间的关系：, 面积坐标与直角坐标的关系,面积坐标与坐标系无关，这是采用面积坐标的优点,面积坐标有性质：,由于面积坐标满足 ，将
21、其代入,得：, 任意三角形到标准等腰直角三角形的变换,将 看作是某一平面的坐标，则上式表明了一种变换。可以证明它把X-Y坐标系的任意三角形单元映射为 坐标系中的标准等腰直角三角形单元。且,即该变换是仿射变换。,这个变换的Jacobi行列式,该变换除了能将三角单元仿射变换为标准三角单元，还能将三角单元上的插值型函数变换为标准三角单元上的同类型函数。因此，采用面积坐标可使计算工作简单化、标准化。,另外，利用面积坐标表示的齐次多项式在(i,j,k)上的积分也非常容易计算。即,其中p, q, r是任意非负整数。,一般在三角形元 上，构造一个m次完全多项式,（3）Lagrange型公式,两个变量x, y
22、的高次多项式可用Pascal三角形表示：,故Lagrange型插值需要 个节点的函数插值条件。,逼近u(x,y)时，该多项式 具有的项数为, 一次多项式,是一次多项式，插值节点数是3。取(1,2,3)的三个顶点为插值节点，运用待定系数法，易得, 二次多项式,是二次多项式，插值节点数是6。取(1,2,3)的三个顶点及三边中点为插值节点。设,运用待定系数法，可得,（4） Hermite 型公式, 二次多项式,在三边中点处的法向导数取指定的,根据插值条件可得,在 的顶点上取指定的 ，,设,其中 是 边长度， 是顶点i, j, k的偏心率。,由三个顶点上的函数值和两个一阶偏导数，加上在三角形形心上的函
23、数值（共十个条件）确定三次插值函数。,设,可求得三次多项式为：, 三次多项式,只利用插值节点处函数值得信息构造出来的插值函数称作是Lagrange型插值函数。,除利用插值节点处的函数值信息，还利用插值节点处的导数信息作为插值条件构造出来的插值函数称作是Hermite型插值函数。,就整体而言，Hermite型插值函数的光滑程度会有所提高。,一般情况下，高次元有限元解的精度比线性有限元解的精度高。但在形成和求解有限元方程组的过程中所花费的工作量也随之猛增。通常除特殊需要，一般不必采用高次元，而是采用线性元、网格细分、外推等措施。,除了三角形元和矩形元，还可考虑任意(甚至是曲边)的四边形单元及曲边三角形单元。,(b) 从虚功原理及节点基函数出发形成有限元方程。,5 用有限元方法解题主要步骤, 对求解域作网格。, 构造基函数（或单元形状函数）。, 求解线性代数方程组。, 形成有限元方程。,(a) 从变分原理及单元形状函数出发，先形成单元刚度矩阵，再由单元刚度矩阵叠加成总刚矩阵。, 把边值问题化成等价的变分问题（或建立等价的变分方程）。,
