微分方程的数值解法

科学研究的基础微分方程的数值解法,计算材料学讲座,基本思路：将连续的偏微分方程及初边值条件离散为线性方程组，并加以求解。,由于离散的出发点不同、所依据的理论不同，因此产生出各种不同的数值方法。,基本方法：有限差分法、有限元法、谱方法、有限分析法。,数理方程数值解法概述,数学基础：数值逼近,不同的数值微分和数值积分方法、不同的函数插值方法，就产生了不同的有限差分法与不同的有限元法。,有限差分法数值微分,有限体积法（也可与有限差分法统称为有限差分法）数值积分,有限元法函数插值,其它数学基础：,数理方程、数值代数、最优化理论与方法等,自然界与工程技术中的许多现象，其数学表达式归结为常微分方程定解问题。,初值问题的数值方法,这里只介绍初值问题。,第一部分常微分方程初值问题的数值解法,常微分方程数值解法主要分为两大部分：,一些偏微分方程问题可以转化为常微分方程问题来（近似）求解。,常微分方程的数值解法为偏微分方程的数值解法提供了可供借鉴的思路。,边值问题的数值方法,目的：建立一阶常微分方程初值问题的数值解法。,以下设不变，记为h。,设初值问题的解析解(理论解)用表示，数值解法的精确解用表示。其中n=1,2,，。,常微分方程初值问题的数值解是求上述初值问题的解u(t)在点列上的近似值。,模型,常微分方程初值问题的数值解法分为：,单(一)步法：计算时，只用到和，即前一步的值。,显式单步法的一般形式为,多步法：计算时，除用到和以外，还用到和,即用到前k步的值。,对单步法与多步法，有显式与隐式方法之分：,显式、隐式多步法的一般形式类似。,隐式单步法的一般形式为,以显式欧拉(Euler)方法为例,上述公式称为（显式）Euler公式。,取h的线性部分，并用表示的近似值，,得,一低精度单步公式建立的基本思想,在附近把根据Taylor级数展开,由给定的初值，则可得到点列上解析解的近似值。,即（显式）Euler公式:,用差商(或数值微分)代替微分方程,略去余项，并用表示的近似值，得到,中的导数。,即,对微分方程两端积分。,用左矩形数值公式计算积分，,运用推导Euler公式的方法可以构造若干公式。,略去余项，并用表示的近似值，则得到（显式）Euler公式:,即,得到高阶方法的直接想法是用Taylor展开。,二高精度单步（Runge-Kutta）公式建立的基本思想,假定初值问题的解u(t)及函数f(t,u(t)充分光滑，则,当h充分小时略去余项，则可导出p阶公式：,其中的各阶导数可表示为,缺点：计算高阶导数十分困难。,Runge-Kutta公式建立的基本思想,用不同点的函数值作线性组合构造近似公式；,再将近似公式与解析解的Taylor展式相比较；,要求二者前面的若干项吻合，从而使近似公式达到一定的阶数。,其中为常数。,R-K法的一般形式为,选取这些常数的原则,要求第一式右端与处Taylor展式中的项重合，从而构造p阶方法。,常用的二阶R-K方法有：,Euler预估校正格式,经典四阶Runge-Kutta公式,常用的四阶R-K方法有：,其中为常数，不全为零。,三线性多步公式建立的基本思想,其中为常数，不全为零。,利用前面多步的信息计算，以获得较高精度的数值公式。,设，的近似值为,并记，k步线性多步方法一般形式为,构造线性多步公式主要有数值积分法和Taylor展开法。,基于数值积分的构造方法,要获得的近似值，只要选用不同的数值方法近似计算其中的积分项，则可导出不同的计算公式。,将两端从到积分，得：,通常对取等距节点构造插值多项式，从而对构造数值积分公式。,常用四阶Adams公式,显式,隐式,Taylor展开更具有一般性。,四阶Adams公式也可用Taylor展开的方法构造。,基于Taylor展开的构造方法,将线性公式写为,比较，要求与有p+1项重合，则可构造p阶方法。,并将其中的函数及(i=1,2,k-1)在点按h的幂展开，然后与微分方程的解u(t)在的Taylor展开式,（1）收敛性、截断误差估计与稳定性,截断误差估计：估计近似替代使数值解与真解之间产生的误差。,稳定性：最初产生的误差在以后各步的计算中是否会无限制扩大。,收敛性：当步长h取得充分小时，数值解能否足够精确地逼近初值问题的真解。,四其它问题,（2）其它方法,改造的方法：外推法、混合法等。,预测校正系统：将同阶的显式公式与隐式公式联合使用。,第二部分偏微分方程的有限差分方法,基本思想：使用离散的、只含有限个未知数的差分方程去近似代替连续变量的微分方程及边值条件，并将相应的差分方程解作为(初)边值问题的近似解。,二阶线性偏微分方程的一般形式为：,对于变量和给定的值和,若,为椭圆型偏微分方程,若,若,为双曲型偏微分方程,为抛物型偏微分方程,二阶偏微分方程的基本分类方法，可以推广到含有两个以上自变量的非线性高阶偏微分方程。,波动方程属于双曲型的方程，热传导方程、扩散方程以及与含时间的薛定谔方程属于抛物型方程，拉普拉斯方程和不含时间的薛定谔方程属于椭圆型方程。,波动方程,扩散方程,拉普拉斯方程,在三维直角坐标系中，其基本形式如下：,泊松方程,双曲型和抛物型偏微分方程通常描述的是非稳态情形，即时间相关问题。在解非稳态问题时必须定义初始条件。,初始条件，就是在起始时间时给定状态变量的数值及其导数值。,对波动方程，初始条件相当于和;,如果没有给定约束条件把解限制在某个特定的空间范围，即，则为初值问题。,对扩散方程或热传导方程，初始条件相当于。,椭圆型方程表述的模型与时间无关，因而描述的是稳态情况。在解稳态问题时必须定义边界条件。这样的问题称为边值问题。,边界条件分为,如果给定初始条件和边界条件，则称之为初边值问题。,Dirichlet边界条件（第一类边界条件）,Neumann边界条件（第二类边界条件）,Cauchy边界条件（第三类边界条件）,一椭圆型方程,二阶线性椭圆型方程,(1.1),边界条件有三种形式：,(1.3),(1.4),(1.2),如果方程在边界上不同部分满足不同类型的边界条件，则称为混合边界条件。,1差分格式建立的基本方法,（1）区域的矩形网格剖分,两组网线的交点称为节点，记为或(i,j)。,也可引入非均匀的矩形网格。,(2)矩形域上的差分形式,在上，进行矩形剖分：,沿x和y方向，分别用中心差商代替(1.1)中的导数：,上两式分别简记为,这里，表示。,则,在(i,j)点被表示为,略去局部截断误差，就得到逼近前述方程的差分格式,(1.5),其中,这里，表示的近似值。,因为上述差分方程中只出现网格函数u在(i,j)点及其上、下、左、右四个邻点上的值，故称之为五点差分格式。,当p(x,y)=1,q(x,y)=0时，方程(1.1)就是Poisson方程：,若用正方形矩形网格剖分区域，即，则,(1.6),求解Poisson方程的差分格式(1.5)简化为,若为矩形区域，且a=b，则M=N=a/h。,利用边界条件(1.2)，则可将(1.6)写成如下形式,(1.7),I为N1阶单位矩阵，而B为N1阶矩阵；,其中,H为阶矩阵。,定义,若用(i,j)点上、下、左、右四个邻点(i-1,j+1)、(i+1,j-1)、(i-1,j-1)及(i+1,j+1)构造Poisson方程的差分格式，可以得到另一种五点格式：,(1.8),将上述两种五点格式组合，可以得到求解Poisson方程的九点差分格式：,(1.9),其局部截断误差达到。,(3)矩形域上边界条件的处理,第二类或第三类边界条件，则需在这些边界点上单独列出差分方程。,第一类边界条件，只需将直接带入在内点列出的差分方程。,如果是一般的二维区域，边界是分段光滑曲线，则临近边界的网格点不一定在网线的交点上；边界的法线方向不一定是水平或垂直方向。,第一类边界条件，有以下两种方法：,（1）直接转移,（2）线性插值,第二类或第三类边界条件，主要考虑如何逼近法向导数。这时分两种情况：,(4)非矩形域上边界条件的处理,如果外法线方向与坐标轴平行，用单侧差商逼近方向导数，也可以用处理矩形区域的方法列出邻近边界的网点上的差分方程。,用单侧差商逼近x方向和y方向的导数，然后列出边界网点上的差分方程。,（1）邻近边界的网格点位于上,（2）邻近边界的网格点不在上,如果外法线方向与坐标轴不平行，其方向余弦是。,由,可以采用直接转移法近似处理，即将边界条件用于邻近边界的网格点，然后再在该点列出差分方程。,2用积分插值法构造差分格式,3差分格式的稳定性和收敛性,4差分方程求解的一些方法,交替方向迭代法,多重网格法,预处理共轭梯度法,(2.1),二抛物型方程,上式中G为x-t平面上的某区域。,三种形式定解问题：,(2.2),初值问题：,则(2.1)、(2.2)构成初值问题（或称Cauchy问题）。,初始条件,一维线性抛物型方程,(2.3),半无界域的初边值问题：,(2.4),方程(2.1)、初始条件(2.3)、边界条件(2.4)构成半无界域的初边值问题。,(2.5),(2.6),有界域的初边值问题：,方程(2.1)、初始条件(2.5)、边界条件(2.6)构成有界域上的初边值问题。,考虑一维常系数热传导方程有界域的初边值问题：,1差分格式建立的基本方法,用平行直线族对x-t平面的求解区域做矩形网格剖分，网点简记为(j,k)，其中和h为时间步长与空间步长，且h=1/N。所有属于G的网点的集合，记为。,(2.7),(2.8),（1）最简显格式,在网点(j,k)处，利用导数与差商的以下关系：,由此得到Lu在网点(j,k)处的关系式：,设表示，表示的近似值。,这儿,略去局部截断误差，得到逼近(2.7)的差分格式：,差分方程(2.9)称为方程(2.7)的最简显格式。它所用到的网点如图所示。,其中称为网比。,(2.9),写成便于计算的形式：,差分格式(2.9)的局部截断误差阶为。,（2）最简隐格式,根据,可得到逼近(2.7)的隐式差分格式：,(2.10),差分方程(2.10)称为方程(2.7)的最简隐格式。其局部截断误差阶为。,所用到的网点如图所示。,（3）CrankNicholson格式,用图内标“”号点的值构成差商，以代替标“”号点在(j,k+1/2)上的导数值：,上两式相减，则有,舍去误差项，得到：,(2.11),差分方程(2.11)称为CrankNicholson格式，它逼近微分方程(2.7)的局部截断误差为。,仍记，则(2.11)可以写为,2差分格式的稳定性和收敛性,3交替方向隐格式及相关格式,三双曲型方程,（a）一阶常系数线性双曲型方程,（b）二阶常系数线性双曲型方程（波动方程）,其中a为常数,主要对象为,这些方程的定解条件，可仅有初始条件，也可以有初始条件和边界条件。,其中a为常数,同椭圆型方程与抛物型方程相比，双曲型方程差分格式的性质与定解问题解析解的性质有更密切的关系。,1一阶线性双曲型方程,（1）初值问题,考虑,由于u(x,t)沿x-t平面上方向为dx/dt=a的直线xat=C（C为常数）的变化率为0，即,故沿x-t平面上任一条斜率为1/a的直线xat=C，u(x,t)为常数。平行直线族xat=C就是方程(3.1)的特征线。,(3.2),(3.1),利用特征线，可以求出初值问题(3.1)、(3.2)的解：,由于u(x,t)在点处的值依赖与(x)在点的值，故初始线t=0上的点称为解u(x,t)在点的依赖区域。,与抛物型方程求解类似，对x-t平面进行矩形网格剖分，x方向的步长为h，t方向的步长为，网点简记为(j,k)。,（1）偏心格式和中心差分格式,对方程(3.1)，利用差商代替导数的方法，可得,前两个格式的局部截断误差阶为，分别称为左、右偏心格式。,第三个格式的截断误差阶为，它称为中心差分格式。,其中,即,从差分格式依赖区域和微分方程依赖区域的关系，可以得到差分格式收敛的必要条件：,差分格式的依赖区域包含微分方程的依赖区域（也称为CFL条件）。,对于左偏心格式，CFL条件为：a0，且。,对于右偏心格式，CFL条件为：a0（或a0，上网点A(j1,k)，B(j,k)，C(j+1,k)处的解值已经算出，从点P(j,k+1)作特征线，它与线段AB交于点D。,（2）Lax格式,由u(p)=u(D)，有,这样，得到Lax格式：,当，Lax格式稳定，截断误差阶为。,（3）LaxWendroff格式,对方程(3.1)，利用特征线作二次插值，即可得到LaxWendroff格式：,当，LaxWendroff格式稳定，它的截断误差阶为。,应该注意：边值条件的给法与其它两类方程不同。,如果a0，方程特征线向右倾，只能在x变化区域的左边界上给出边界条件：,（2）初边值问题,如果a0，考虑下面模型问题：,前面建立的几个显格式，都适用于这个问题。,下面建立隐格式。,连同初始条件与边界条件：,（1）最简隐格式,该格式的局部截断误差阶为。,令，格式可改写为,该格式可在0x,t内所有网点上显示地计算解之近似值。,然后用中心差商逼近这些导数值，则可得到Wendroff格式：,在点处，用,（2）Wendroff格式,连同初始条件与边界条件：,该格式的局部截断误差阶为，且无条件稳定。,令，格式可改写为,该格式可在0x,t内所有网点上显示地计算解之近似值。,2二阶线性双曲型方程（波动方程）,考察,对x-t平面进行矩形网格剖分，x方向的步长为h，t方向的步长为，网点简记为(j,k)。,用二阶中心差商代替(3.3)中的二阶导数，则得到网点(j,k)处的差分方程：,(3.3),(3.4),其中。,或,该格式稳定的充分条件为。,初始条件离散：,由,消去，得,上述差分格式与初始条件的截断误差阶均为。,取为,上述方法也可用于求解初边值问题：,3交替方向隐格式,任何模拟方法，都必须在最佳计算速度和数值精度之间寻找平衡点。,要在各种可能的求解方法中找到一种统一地适用于计算材料学领域（或其它领域）的理想方法，一般是不现实的。,由于实际问题的具体特征、复杂性以及算法自身的适用范围决定了应用中必须选择、设计适合于自己特定问题的算法，因而掌握数值方法的思想至关重要。,科学计算（数值模拟）已经被公认为与理论分析、实验分析并列的科学研究三大基本手段之一，但三者之间相辅相成。,谢谢大家！,第三部分偏微分方程的有限元方法,一边值问题的变分原理,1引论,求使得泛函,达到最大的函数。,（1）等周问题,在长度一定的所有平面封闭曲线中，求所围面积为最大的曲线。,定义：当求泛函在一个函数集合K中的极小（或极大）问题，则该问题称为变分问题。,变分问题与微分方程的定解问题有一定的联系。,（2）初等变分原理,一元二次函数的变分原理,考察J(x)的极值情况。,变分原理：,设,求，使,与求解方程Lx=f等价。,对称正定,多元二次函数的变分原理,求J(x)取极小值的驻点，其中,设,设,则J(x)可表示为：,变分原理：,设矩阵A对称正定，则下列两个命题等价：,(b)是方程的解,上述两个例子表明：,其中,求二次函数的极小值问题和求线性代数方程（组）的解是等价的。,（1）弦平衡的平衡原理与极小位能原理,2两点边值问题的变分原理,考察一根长为l的弦，两端固定在点A(0,0)和B(l,0)。当没有外力作用时，它的位置沿水平方向与X轴重合。设有强度为f(x)的外荷载垂直向下作用在弦上，于是弦发生形变。假定荷载很小，因而发生的形变也很小。用u(x)表示在荷载f(x)的作用下弦的平衡位置。,求弦的平衡位置归结为求解两点边值问题：,设弦处于某一位置u=u(x)，可得到其总位能为,极小位能原理：,其中T是弦的张力。,平衡原理,弦的平衡位置(记为)将在满足边值条件u(0)=0，u(l)=0的一切可能位置中，使位能取极小值。,弦的平衡位置是下列变分问题的解,在数学上，要将某个微分方程的定解问题转化为一个变分问题求解，必须针对已给的定解问题构造一个相应的泛函，并证明定解问题的解与泛函极值问题的解等价。,有限元方法正是利用这种等价性（边值问题与变分问题的等价性），先将微分方程定解问题转化为变分问题（或变分方程）的求解问题，然后再设法近似求解变分问题（或变分方程）。,（2）两点边值问题的变分原理,构造泛函,考察二阶常微分方程边值问题：,引入泛函算子,则,变分问题,与前述二阶常微分方程边值问题相应的变分问题是,其中,求，使,变分原理（变分问题与边值问题的等价性）,设，是边值问题,的解，则使J(u)达到极小值；,反之，若使J(u)达到极小值，则是边值问题的解。,其中,是强制边界条件，是自然边界条件，区别这两类边界条件在用有限元方法求解边值问题时很重要。,（3）虚功原理,对两点边值问题：,其中,虚功原理,，且满足变分方程：,设，以v乘方程两端，沿a,b积分，并利用，得变分方程,对任意,在力学里，表示虚功,设，则是边值问题解的充要条件是：,对于复杂的边界条件，边值问题的求解一般是困难的。若将微分方程化为相应的变分问题或变分方程，则只需处理强加边界条件，无需处理自然边界条件（自然边界条件已包含于变分问题中泛函的构造或已包含于给出的变分方程之中）。这一特点对研究微分方程离散化方法及其数值解带来了极大的方便。,3二阶椭圆边值问题的变分原理,（1）极小位能原理,模型方程,其中G是平面有界区域。,构造泛函,引入泛函算子,则,变分问题,与前述二阶椭圆边值问题相应的变分问题是,求，使,其中,变分原理（变分问题与边值问题的等价性）,对第一边值问题，无论齐次或非齐次边界条件，泛函是一样的，只是边界条件要作为强加边值条件加在所取的函数类上。,设，是二阶椭圆边值问题的解，则使J(u)达到极小值；,反之，若使J(u)达到极小值，则是二阶椭圆边值问题的解。,其中,对第二、三类边值问题，无论齐次或非齐次边界条件，二次泛函形式相对于第一边值问题有所改变，但函数类的选取与边界条件无关。,（2）虚功原理,问题,其中,设，以v乘方程两端后在G上积分，并利用Green公式，得变分方程,虚功原理,在力学里，表示虚功,设是边值问题的解，则对任意，满足变分方程。,反之，若，且对任意满足变分方程，则为边值问题的解。,与极小位能原理类似，第一类边界条件为强加边界条件，第二、三类边界条件为自然边界条件。,虚功原理比极小位能原理应用更广。,目的：求解相应的变分问题或相应的变分方程。,Ritz方法是近似求解变分问题(即二次泛函极小值)的算法。Galerkin方法是近似求解变分方程的算法，这两种算法统称为Ritz-Galerkin方法。,Ritz-Galerkin方法的基本思想,以下用V表示等Sobolev空间，L表示微分算子，(u,v)为由L及边值条件决定的双线性泛函。,4Ritz-Galerkin方法,用有限维空间的函数代替变分问题(或变分方程)中无限维空间的函数，从而在有限维函数空间中求变分问题(或变分方程)的近似解，并要求当有限维空间的维数不断增加时，有限维近似解逼近原变分问题(或变分方程)的解。,由极小位能原理得出的变分问题为：,Ritz方法：求变分问题的近似解。,（1）Ritz方法,求，使,其中，,设是V的n维子空间，是的一组基底(称为基函数)。中任一元素可表示为,即选择适当的，使取极小值。,求，使,Ritz方法：,展开,令,则满足,解出代入，则得,Ritz方法步骤为：,根据最小位能原理构造相应于微分方程或物理问题的变分问题；,取作为的一组基底，即用近似代替无穷维空间V；,求解关于的线性代数方程组。,由虚功原理得出的变分方程为：,Galerkin方法：求变分方程的近似解。,（2）Galerkin方法,设是V的n维子空间，是的一组基底(称为基函数)。中任一元素可表示为,即选择适当的，使取极小值。,Galerkin方法：,求，使对，满足,由的任意性，取作为v，则得,将代入变分方程，则,解出代入，则得,Galerkin步骤为：,根据虚功原理构造相应于微分方程或物理问题的变分方程；,取作为的一组基底，即用近似代替无穷维空间V；,求解关于的线性代数方程组。,取作为v，将代入变分方程，得到满足的方程组：,有限元法广泛应用的原因,Ritz-Galerkin方法应用的困难,基函数选取必须满足强加边界条件，因此选取困难；,计算量、存储量巨大；,方程组求解病态严重。,充分发挥了变分形式和Ritz-Galerkin方法的优点；,摆脱了传统的基函数取法；,各种问题的结构程序格式统一。,有限元方法基于变分原理，又具有差分方法的一些特点，并且适于较复杂的区域和不同粗细的网格。,二椭圆型方程的有限元方法,差分法解偏微分方程，解得的结果就是准确解u在节点上的近似值；,Ritz-Galerkin方法得到近似的解析解，但对一般区域，却往往难以实现。,有限元方法与传统Ritz-Galerkin方法的差别在于有限维函数空间的构造方法。Ritz-Galerkin方法选用的基函数在整个定解区域上整体光滑，有限元则取分段或分片连续且局部非零的基函数。,考虑两点边值问题：,1一维问题的线性元,将区间a,b分割为n个子区间。,第i个单元记为，其长度。,（1）试探函数与试探函数空间,设,则称为试探函数空间，称为试探函数。,（2）用单元形状函数表示试探函数,设在节点上试探函数在节点上的一组值为,最简单的试探函数空间由分段线性函数组成。,在第i个单元上的线性插值函数为,即,当时，的（线性）插值公式称为（线性）单元形状函数。,把每个单元形状函数合并起来，就得到整个区间a,b上都有定义的函数：,为使分段插值标准化，通常用仿射变换,显然,把变到，令,则,变为,或,定义基函数系,（3）用节点基函数表示试探函数,线性无关，它们可组成试探函数空间的基，常称为节点基函数。,几何形状如图,任一试探函数可表示为,用这类插值型基函数，可以构造出适合各种边界条件的试探函数。,若借助前述放射变换,节点基函数可用变量表示为,直接形成有限元方程,(a)把表达式代入泛函；,（4）从Ritz方法出发形成有限元方程,(b)将泛函表达式中积分区间a,b变到0,1；,(c)由达到极小值的条件,得到含的有限元方程,这儿,(d)解出有限元方程的数值解，就得到使二次泛函取极小的近似函数(有限元解),有限元方程可用矩阵表示为,其中,称为总刚矩阵。,工程中形成有限元方程时，通常先在每个单元上形成单元矩阵（称为单元刚度矩阵），然后由单元刚度矩阵形成总刚度矩阵（称为总体合成）。,用单元刚度分析形成有限元方程,(a)把按单元组织，则在第i个单元上，令,其中称为单元刚度矩阵。各元素可计算得到。,再把扩展成nn矩阵，使其第i1行、第i行和第i1列、第i列交叉位置的元素就是单元刚度矩阵的四个元素，其余全为零（只是第一行，第一列元素非零）。即,记,则,其中称为总刚矩阵。,(b)由达到极小值的条件,(c)解出有限元方程的数值解，就得到使二次泛函取极小的近似函数(有限元解),得到有限元方程。,（5）从Galerkin方法出发形成有限元方程,把表达式代入变分方程,对前面的两点边值问题，变分方程变为,其中,与Ritz方法相比，Galerkin方法形成的有限元方程其系数矩阵就是总刚矩阵。,该方程即为Galerkin法形成的有限元方程。,由Galerkin方法推导有限元方程更加方便直接，且适用面广。,若希望在每个单元上提高逼近的精确度，则可通过提高插值多项式次数来实现，,在单元上可构造一、二、三及高次插值多项式，其方法有两种：,2一维问题的高次元,整个问题计算的全过程除分析单元插值外，均与前面框架类似。,Lagrange型：在单元内部增加一些插值节点。,Hermite型：在节点引进一阶、二阶乃至更高阶导数。,线性元(Lagrange型),要求：在每一个单元上是一次多项式，在单元节点处连续。插值条件：在单元的两个端点取指定值。,二次元(Lagrange型),要求：在每一个单元上是二次多项式，在单元节点处连续。插值条件：在单元的两个端点及单元中点取指定值。,三次元（Hermite型）,要求：在每一个单元上是三次多项式，在单元节点处连续。插值条件：在两个端点取指定的函数值和一阶导数值。,采用高次元，有限元方程形成的方法和线性元类似，但工作量增加。一是计算积分的复杂性增加，二是矩阵的带宽增加。,高次元的主要优点是收敛阶高，且提高了函数逼近的光滑性。,假定区域G可以分割成有限个矩形的和，且每个小矩形（单元）的边和坐标轴平行。,3二维问题的矩形元,通过仿射变换,采用矩形剖分后，任一个矩形,总可变成单位正方形,如果在上造出单元形状函数，就可得到试探函数。而上的形状函数可通过先在上造出形状函数，再通过仿射变化而得到。,在上构造形状函数，也采用Lagrange型和Hermite型插值。,Lagrange型：根据若干插值节点处的函数值决定插值函数。,Hermite型：根据若干插值节点处的函数值、一阶偏导数乃至更高阶偏导数决定插值函数。,（1）Lagrange型公式,双一次插值,插值条件：给定顶点上的函数值,求：双线性函数,满足,设,令,由为双线性函数，可求得,令,则,通过仿射变换消去、，就得到上的形状函数。把这些函数按单元叠加，即对所有单元求和，就得到G上的试探函数。,实际计算时，并不消去中间变量、，因为计算刚度矩阵元素（定积分）用、作自变量更为方便。,插值条件：给定II上九个插值节点(0,0)、(1/2,0)、(1,0)、(0,1/2)、(1/2,1/2)、(1,1/2)、(0,1)、(1/2,1)、(1,1)的函数值。,求：双二次函数,满足,双二次插值,故,通过仿射变换消去、，就得到上的形状函数。,令,由为二次函数，可求得,设,插值条件：给定II上十六个插值节点(见图)。,求：双三次函数,满足,设,双三次插值,故,令,由为三次函数，可求得,可以在四个顶点分别给定函数值、两个一阶偏导数的值和二阶混合偏导数的值(共十六条件)，确定一个双三次多项式的十六个系数。,（2）Hermite型公式,Lagrange型公式中不出现导数，这样的试探函数只属于。为了得到属于的试探函数，需要Hermite型插值公式。,双三次多项式含有十六项：,简单且常用的是不完全的双三次多项式插值。它去掉双三次多项式中的项。,插值条件：给定II上四个插值节点。,求：不完全双三次函数,满足,四个顶点处的函数值等于在该点的函数值；四个顶点处的值等于在该点的值；四个顶点处的值等于在该点的值。,根据仿射变换,则可将原插值问题转化为II上的插值问题。,满足,四个顶点处的函数值等于在该点的函数值；四个顶点处的值等于在该点的值乘以x；四个顶点处的值等于在该点的值乘以y。,插值条件：给定II上四个插值节点(0,0)、(1,0)、(0,1)、(1,1)。,求：不完全双三次函数,类似于Lagrange型公式的构造，可以求得上的形状函数。,在三角形元的有限元方法中，先将定解区域G化分为若干个小三角形（称作单元）。然后在每个单元上构造插值型函数，并用分片函数（但整体连续的函数）代替变分问题或变分方程中所需求解的函数。,4二维问题的三角形元,用有限元求解二维椭圆边值问题时，应用最广的是三角形元。,（1）三角剖分,将定解区域化分成若干个小三角形单元时应注意：,为了保证有限元解的精确度和收敛性，并避免其离散后代数方程组系数矩阵的病态性，网格剖分中疏密的过渡不要太陡。,错误,为了保证有限元解有较好的精度，每个单元中应尽量避免出现大的钝角。,应避免,单元顶点的编号顺序可以任意，但节点编号顺序将影响有限元方程组系数矩阵的结构（带宽）。,为了方便构造插值型函数，要求每个单元的顶点是相邻单元的顶点。,（2）面积坐标及有关公式,在三角