
本章介绍关于有限元方法的一些数学概念和结论,目的在于使读者对于有限元解的收敛性以及单元精度问题能有确切的了解。以后各章的内容在本章提供的基础之上进行。对于有限元方法的数学研究,目前已进行得相当充分,对这方面有兴趣的读者可进一步查阅有关的专著。
πP
本章介绍的主要对象是函数:真实解是一个函数;基函数是一组函数;试探函数是某一类函数,有限元解是这类函数中使 取驻值(最小值)的那一个函数。下面讨论中的 “元素”实际指的就是函数,“空间”实际指的就是某种函数的集合,即函数空间。
§4-1线性空间(向量空间)
1、线性空间的定义
满足下列条件的空间E为线性空间
(1) x, y , zE有如下“加法”运算
(i)
(ii)
(iii)存在“零元素”θE x E有
(iv) x E存在逆元素-xE使
(2)设E中的元素与实数域的元素有“数乘”运算,即 x, y E,α,β K(实数域)
(i)
(ii)
(iii)
(iv)
若K为实数域则E称为实线性空间,K为复数域则E称为复线性空间。
例1 C[a,b]
若 、 是[a, b]上的连续函数,则 也是[a, b]上的连续函数。故定义在[a, b]上的所有连续函数组成一个线性空间。记作C[a, b]。
例2 L2(a,b)
若 、 是(a, b)上平方可积的函数,即
,
存在,则
所以 也是(a, b)上平方可积的函数。所有(a, b)上平方可积的函数组成一个线性空间,记作L2 (a, b) 。
例3 C1[a,b]
若 、 、 、 在[a, b]上连续,则
也在(a, b)上连续。所有函数本身及一阶导数都在(a, b)上连续的函数组成一种线性空间,记作C1[a, b]。
例4 Rn n维欧氏空间是线性空间,R2(二维平面), R3(三维空间)是n维欧氏空间的特例。
例5 Pn(x) [a,b]
定义在[a,b]上的n次多项式Pn(x) [a,b]C[a,b] 构成线性空间。
2、线性空间的维数
(1)线性相关与线性无关
设 为线性空间E的n个元素
(i)若存在不全为零的常数1、2…n 使得
则称 线性相关;
(ii) 若
仅当
才成立,则称 线性无关。
(2)线性空间的维数
若线性空间E满足
(i)任意n+1个元素一定线性相关。
(ii)存在着n个线性无关的元素。
则称线性空间E的维数为n。
例1 若 线性无关,则所有形式为
的试探函数组成n维线性空间。而所有形式为
的位移场则组成2n维线性空间。
例2 由于可以找出任意多个线性无关的连续函数(例如:1 ), 所以C空间为无限维线性空间。L2空间也是无限维线性空间。
3、线性空间的模(范数)
(1)模的定义
当线性空间E中的任意一个元素x可用一个非负实数与之对应,记作‖x‖(表示“大小”或“长度”)称为E空间为模线性空间或赋范线性空间,实数‖x‖称为模或范数。模的性质如下:
(i)‖x‖≥0 ,仅当x≡0时‖x‖=0
(ii)对任一常数α有‖αx‖=∣α∣•‖x‖
(iii) 对任意 x、y∈E 有
‖x+y‖≤‖x‖+‖y‖ (此式又称三角不等式)。
只要满足这些要求,均可作为一种模的定义。‖x‖描述了元素x的“大小”,‖x-y‖ 则描述了两个元素x、y之间的“距离”。设真实解为u,有限元解为uh,如果当网格无限细分时有‖u-uh‖→0,则说明有限元解收敛于真实解。模的定义不同,收敛的意义也不同。
例1 在平面R2内,向量x(x1, x2)可以有下列三种模的定义
例2 设x, y E则,x-y的模可以表示这两个元素的“接近程度”,若在R2空间中的两个元素,x (1,1), y (2, 4) 可以有如下模的定义
可见模的定义不同,其意义不同, 在实数域内,
模与绝对值是等价的。
(2)两种常用的模
下面是最常见、也是具有代表性的两种模
① 一致模 若u∈C[a, b],则u必在[a, b]上取到最大值和最小值,故可按如下方式定义一致模‖u‖∞。
② L2模 若u∈L2(a, b) 则
存在,可以按如下方式定义L2模
显然按一致模收敛是一致收敛,按L2模收敛则是平均收敛。
对于有限维空间,不同类型的模是等价的。即对某种模成立的结论,对另一种模也成立。而对于无限维空间则未必如此。
§4-2内积空间
内积空间扩展了向量点积和正交的概念。
1、内积 对于线性空间E的每一对元素u、v定义一个确定的实数与之对应,称为u、v的内积,记作(u、v),且满足
(i) ( u、v)=(v、u) (对称性)
(ii) 对任一常数α有(αu, v)=α(u, v) (齐次性)
(iii) 对于u1、u2、v∈E 有
(u1+u2, v)=( u1, v)+( u2, v) (可加性)
(iv) (u, u) ≥0 仅当u≡0时 (u, u)=0。
定义了内积的线性空间称为内积空间。
例 u(x), v(x) C2[a, b] 至少存在以下四种形式的内积
其中是[a, b]上的给定函数。
2、内积模
对于内积空间,可以直接利用内积来定义元素的模
模的定义并不是唯一的,如果确有必要,我们仍然可以再定义其他形式的模(例如一致模)。但对于内积空间用得最多的当然是由内积定义的模。
在内积空间E中,u与v之间的距离可用内积模表示为
3、正交性
内积空间与一般线性空间的不同之处是可以用内积来定义两个元素之间的正交关系,函数之间的“正交”。
若( u、v)=0 则称u、v 正交。显然,模及正交性涵义取决于内积的定义。
例1 若积分
存在,则可定义
则模
而 u、v 正交则意味着
例如在[0, π]上当m≠n时 sinmx与sinnx正交。
例2 若ρ>0且积分
存在,也可定义
则模
而 u、v 正交则意味着
即通常理解的u、v 以ρ为权正交。
4、Schwarz 不等式
设u、v 是内积空间的两个元素,t 为任一实数,则tu-v也是内积空间的一个元素,显然,它自身的内积
展开
可以证明,上式对任何实数t都成立的充分必要条件是它的判别式
即
或
这就是著名的Schwarz不等式
对于a、b两个向量的内积显然有
上式为Euclid空间的三角不等式,此式仅是Schwarz不等式的一个特例。
5、收敛性与完备性
(1)收敛性
点列E(赋范线性空间),若存在
则,称为点列的强极限,读作:强收敛于,注意模的定义不同收敛的涵义不同。
(2)完备性
若E空间中的每一个元素列收敛于E中的一个元素,则称空间E是完备的。例如
Hilbert空间---完备的内积空间。
Banach空间---完备的赋范线性空间。
是Hilbert空间的子空间。
是的子集。
许多数学物理问题的解存在于某一类函数空间中。换句话来说为了得到有意义的解,必需明确解的存在空间。即对组成解的函数的类型作一限定。值得注意的是,如果的过于严格会将一些解排除在外,过宽可能导致解无意义。有限元解的存在空间为索伯列夫空间,是Hilbert 空间的子空间。
§4-3索伯列夫空间HK
1、HK 空间的定义及实例 广义导数
当k为非负整数时,HK(Ω)表示在定义域Ω内,函数本身以及直到k阶导数都平方可积的函数全体。(当k为负数或分数时将赋于另外的意义。本章不讨论这种情况。)
下面结合几个具体实例说明有关概念。
(1)若u(x)∈H1(0, L),则积分
、
一定存在。对于图4-2所示的分段线性插值函数u(x)而言,显然
和
都存在。在点 x2、x3、x4处,按通常的意义u;不存在,左导数不等于右导数。我们不妨补充定义:在这些点上,取两侧一阶导数的平均值作为该点的导数值(广义导数)。定义了广义导数的空间就是完备的线性空间,由于在通过结点时不连续,有限跳跃量,在结点处为δ—函数。δ—函数本身可积,但平方不可积。我们可以得到这样的结论:对于分段线性插值的函数u(x)而言:u∈H1(0, L),但uH2(0, L)。
(2)设Ω为一平面二维区域,将Ω分成若干三角形的片(单元)(图4-3)。取各三角形的项点为结点,以结点处的函数值对单元内的位移场进行分片线性插值。根据第3-4节的分析可知,对于这样定义的函数u(x,y)在Ω上连续,且积分
、 、 存在。
在单元边界和结点处,通常意义下 、 不存在。我们不妨补充定义: 、 在单元边界上取边界两侧的平均值,在结点处取包围这个结点各单元的加权平均值(按各单元所占有的角度加权)。在单元边界和结点处 、 一般不连续, 、 、 为一δ函数,故对于二维分片线性插值的函数u,有u∈H1(Ω),但uH2(Ω)。如图4-3所示。
(3)在研究梁的弯曲时,曾以结点处的函数值和一阶导数值作为结点参数,采用分段三次Hermite插值来构造试探函数v(x)(图4-4)。根据第3-3节中的分析可知,v(x)、v’(x)在 [0, L]上连续,
, 和
存在。
在通常意义下,结点处v’’不存在。我们可以补充定义:结点上的v’’值取为两侧v’’的平均值。在结点处v’’’为一δ函数,平方以后不可积。因此,对于分段三次Hermite插值的函数v(x) 有vH2(a, b), 但uH3(a, b)。
综上所述,为了处理以分片插值定义的函数(这种形式为处理强制边界条件提供了很大方便)引入了广义导数的概念。我们没有涉及广义导数的数学定义,但对广义导数所做的直观上的解释(理解为某种意义下的平均值)已能够帮助我们理解一些重要的结论。HK空间中所提到的导数,都是指广义导数。不难看出,如果通常意义下的导数存在,它将与广义导数完全一致。
任何一个HK空间都是无限维线性空间。显然,H2(Ω)是H1(Ω)的子空间;H3(Ω)是H2(Ω)的子空间又是H1Ω)的子空间。H0(Ω)只要求函数本身平方可积,与L2(Ω)是一回事。
结论:一般说来,用分片插值多项式定义的函数,在单元内部的可微性都比较好。在整个区域上的可微性主要取决于插值函数在穿过单元边界时的性质。
2、索伯列夫空间的模
设Ω为一平面二维区域,当u∈H1(Ω)时,定义
当u∈H2(Ω)时, 定义
这样的定义意味着,若u为位移函数,则在讨论收敛性时,若收敛一定意味着函数本身(位移)的收敛性,也包括了函数的若干阶导数(自然包括了应力)的收敛性。
3、索伯列夫空间的半模 (积分中不含函数自身)
当Ω为一平面二维区域时,半模的定义下:
当时
当时
显然,上两式定义的半模描述了某一阶广义导数的“大小”。
4、能量模和能量内积
众所周知,弹性体的变形能总是非负的。因此,可以用一非负实数与其对应。如果所加的位移约束条件使整个系统不能做刚体运动,那么只有在系统的位移恒为零时总变形能才等于零。在这种情况下就允许我们把变形能作为一种模的定义(能量模),并且由此受到启发去定义能量内积。变形能的具体表达式依赖于问题的类型和边界条件:对于二阶问题(例如轴向受拉的杆和平面应力问题)变形能依赖于位移的一阶导数(应变);对于四阶问题(例如梁的弯曲)变形能依赖于位移的二阶导数(曲率);如果系统具有弹性支承,那么变形能还将依赖于支承刚度和支承点的位移。
为了与通常的内积定义区别,将u, v 的能量内积写成D(u, v)而u的能量模则写成 记作‖u‖H。对于图4-5所示
的轴向受拉杆,有
对于图4-6所示弹性基础上的简支梁,若支承刚度为k(x), 则能量内积和能量模的平方分别为
当k≡0时
即为通常的弯曲变形能。
对于图4-7的平面应力问题,要略为繁琐一些。若以
为表位移,以
代表另外一种位移,则能量内积
能量模的平方为 等于
总之,对于具体问题。我们总可以写出D(u, v)、D(u, u) 的具体表达式。尽管有些情况下较为繁琐,但却可以有“一劳永逸“的效果。外力f在位移场 u上作的功也可表示为内积形式(f, u)那么系统的总势能则可以十分简洁地写成
而势能驻值条件则可统一地写成
这里u泛指一般的位移场,而f 则泛指一般外载荷。这样借助能量内积和能量模的概念,我们就可以暂时摆脱问题的具体物理背景,而理论分析得到的结论将适用于任何一个具体的问题(椭圆型方程边值问题)。而当需要作某些直观解释时,只要选择一些较为简单的问题即可达到目的。鉴于上述理由,本章下面的讨论将仅限于二阶问题,对一维情况仍然以受拉杆为例,对二维情况将以泊松方程(薄膜问题)作为代表,这样在书写上要简单很多,而又不失去结论的普遍意义。
§4-4插值逼近的误差分析
1、一维情况
单元e为长h的区间,设真解u∈H3(0, h)。A1、A2、A3为曲线u(x)上的三个点,过这三个点作一条二次曲线p2(x.) (图4-8)。现在研究以p2(x) 代替u(x) 造成的误差。
设误差函数为
所以有
且 、
在[0, h]上连续。若A1、A2、A3的坐标为x1、x2、x3 , 由于A1、A2、A3在曲线u(x)及p2(x)上,必有
因而必存在两个点 x1*, x2* , x1 因此又可以找到一点x**, x1* 若将二次导数用三次导数来表示,则 对于误差函数的二次导数下式成立 利用Schwarz不等式,则 上式中为函数u对x的三次导数的模的平方,进一步可得到误差函数的二次导数的半模与函数u的半模之间的关系 类似有 以及 即 上述讨论中认为插值点上取到的函数值是精确值。 而有限元解求得的结点上的函数值(位移值)只是近似值。后面可以证明,这种差别并不会带来不利影响。 综上所述可得到如下重要结论: (1)适当增加单元内插值点个数,增加插值多项式次数,有可能提高精度,但要受到函数本身(真解)可微性的。当u∈HK(e) 、(k≥2) 时采用k―1次多项式可以达到最高的精度(函数本身的误差为hk,e)。再增加插值多项式次数,精度当然不会降低,但也未必能够提高。当u∈H1(e) 时一般采用线性插值。 (2)函数本身的精度最好,导数的精度低于函数,导数的阶数越高,精度越低。 2、一般情况 设e为一个一般性的单元(在平面情况下即为三角元或矩形元)。若函数(真解) u∈HK(e), (k≥2)。在单元内用一个多项式逼近函数u;多项式中k―1次以及低于k―1次的项是完全的,k次以及高于k次的各项或者没有,或者有缺项;多项式的系数由结点处函数u(必要时还有u的导数)的精确值定出。这样的插值多项式记作πk―1u。其中下标k―1意味着插值多项式完全到k―1次,即当u(真解)为不超过k―1次的任何多项式都有 对于三结点三角元和四结点矩形元有 k―1=1 (线性位移场) 对于六结点三角元和八结点矩形元有 k―1=2 (二次函数描述的位移场) 对于这样的插值多项式,可得到如下的误差估计: (4-4-1) 当k≥2时有 (4-4-2) (4-4-3) 当 k≥3时还有 其中h为单元直径,即单元内任意两点之间的最大距离,对三角形为最长边,对矩形为对角线。C代表与u和、h无关、只与单元形状(对三角元指内角,对矩形元指长宽比)有关的常数。从形式上看这些估计式与一维情况类似,但要证明它所涉及的数学知识要多得多。 至此我们可以得到几条量的界限: (1)当u∈H2(e) 时,采用完全到一次项的多项式,当h→0时函数误差为Ο(h2)阶,一阶导数误差为Ο(h)阶。 (2)当u∈H3(e) 时,可以采用二次插值多项式,提高收敛速度。当u∈H1(e) 时前述误差估计式不能直接用,但仍可采有线性插值,函数和一阶导数的收敛性可以加以证明,但收敛速度可能较慢。 可见有限元解的误差不仅取决于插值多项式次数,而且有赖于真实解的可微性(实际的位移场的平缓性)。 §4-5广的可微性 角点 1、有限元解的容许空间 满足强制性边界条件和协调条件、且使总势能 取驻值(最小值)的u称为相应的椭圆型方程边值问题在Ritz意义下的广。对于二阶问题D(u, u)可能包括u和u的一阶导数,为使πP存在必须要求 u∈H1(Ω)。 H1(Ω)称为二阶问题的容许空间。而四阶问题的容许空间则是H2(Ω)。这是对广可微性最起码的要求。而广的可微性(特别是在Ω的内部子域上)可能比这最低要求要好,因而采用高次插值有可能改善有限元解的精度。一般说来,导致广可微性差有下列三方面的因素:材料性质不连续、载荷不连续以及角点。下面分别讨论它们的影响和处理方法。 2、材料性质不连续的影响 假设求解区域Ω由两种材料组成,材料分界线Г1是一光滑曲线。只要在划分单元时使Г1作为单元边界,那么这种材料性质的不连续虽然降低了广在整个区域Ω上的可微性,却不会降低在单元内的可微性。(图4-11) 3、载荷不连续的影响 当分布载荷不连续时,只要取载荷不连续处作为单元边界,这种不连续性就不致降低单元内广的可微性。 例如,对于一维二阶问题(杆拉伸),允许内部载荷为集中力,当取力作用点为结点时,单元内广的可微性不会因此而降低。(图4-12) 对于二维二阶问题(膜、平面应力)不允许在Ω内施加集中力(过分奇异)。但对二维四阶问题(板弯曲)允许施加集中力。 4、角点 我们看到,对于材料性质不连续和分布载荷不连续这类因素,只要在划分单元时采取一定措施,这些因素并不会给我们带来更多的困难。而角点的情况则完全不同。所谓角点,是指边界上不光滑的点,例如图4-13(a)中的A、B、C、D、E、F各点。当不同材料的分界面Г1不是一条光滑曲线时,在不光滑的点(例如图4-13(b)中的P点)上会 出现与角点类似的情况。 在图4-13(c)中尽管边界Г本身是光滑的,但由于一部分边界自由,一部分边界施加了位移约束条件,分界点P1、P2仍然属于角点。在角点附近往往伴随着应力集中现象,应力变化急剧。尽管我们仍然可在划分单元时使分界线为单元边界,取角点为结点,但是单元内假设的多项式形式的插值函数仍然不能很好地逼近角点附近急剧变化的应力场。下面先以薄膜问题为例。在角点P附近挖出一个小扇形区,取P为极点,则在极坐标下泊松方程的形式为 设角点附近的边界为直线,沿此边界施加u=0的边界条件,扇形夹角 φ=απ(0<α≤2)设 代入方程可定出 uk(r) 则u可能包含如下形式的奇异项 , 在第二种情况下,奇异性最严重的当然是k=1的项,即 由 , 中可能出现最严重的奇异项估计u的可微性,如下表所示 的可微性估计 我们看到,对任何φ值u都没有超出二阶问题的容许空间H1 。但在φ>π即出现凹角时,u不在属于H2。特别是在φ=π时, 具有 的奇异性,对应的位移项为 平面应力场的裂纹尖端存在着类似的情况(图4-14(a))。对于固定 自由边界的分界点上,通过对称延拓,不难看出,与内角为2φ的自由边界相当(图4-14(b)) 对于凹角(特别是对于裂纹),采用大致均匀的单元网格即使是采用高次插值多项式也不会在奇点附近取得较好的逼近效果(因为在奇点附近,真实解与多项式相甚远),而且这种误差会传播到不包含奇导性的其他单元。要改进逼近效果,可以有两种方法。一种是采用逐步加密的单元网格,这种网格在奇点附近的单元取的很小,但只采用低阶插值多项式。另一种方法是在奇点附近的单元上假定具有某种奇异项的位移场(即,采用与远离奇异点的单元不同的单元),例如在裂纹附近的单元中采用除完全的一次项外还可以增加相当于 的项的插值方式。 §4-6协调位移单元的收敛性和误差估计 1、有限元空间Sh 假设我们已经做到了以下几点 (i) 将区域Ω划分为m个单元(在平面情况下可以是三角元或矩形元)。 (ii) 配置好每个单元的结点,选择好每个结点参数。对于二阶问题总以结点函数值为结点参数,以下我们仅限于讨论二阶问题。 (iii) 选择每个单元内的插值多项式。多项式满足协调性和可微性要求,且多项式的系数要能由单元所包括的结点参数唯一确定。 利用上述三点可得到一组基函数φ1、φ2、…φn。它们满足协调性和可微性要求;每个基函数只在一个结点上为1,在其余结点上为零。 我们把形如u=∑uiφi (其中φi为定义的基函数,ui任意常数)的函数全体称为由基函数φ1、φ2、…φn张成的有限元空间,记作Sh。显然Sh是一个有限维空间,是真解的容许空间(二阶问题的真解的容许空间为H1(Ω))的一个子空间问题。以精确的结点参数值进行插值所得到的分片插值函数πu是Sh的一个元素,也是H1(Ω)的一个元素。有限元解uh的结点参数由δπP=0定出,一般仅是近似值。故一般说来,uh即是Sh的一个元素,也是H1(Ω)的一个元素,但是与πu不同的元素。 2、有限元解 投影关系 设总势能 的势能驻值条件可表示为:对任意δu有 (4-6-1) 二阶问题的真实解u可以表述为:寻找一个u∈H1(Ω), 使得对任何 δu∈H1(Ω)(记作 δu∈H1(Ω))都有 (4-6-2) 有限元解uh则可表述为:寻找一个uh∈Sh,使得对 δu∈Sh 都有 (4-6-3) 注意到Sh是H1(Ω)的一个子空间(记作Sh H1(Ω)),将(4-6-1)与(4-6-2)相减。对 δu∈Sh都有 (4-6-3)表明:真实解u与有限元解uh之差与Sh的任一元素在能量内积的意义下正交。(图4-15)。换句话说:uh是u在Sh上的投影。 从另一个观点看,设v是Sh的一个元素,则v-uh∈Sh,即在能量模意义下,uh是Sh的所有元素中与u“最接近”的一个元素。特别是,对于πu有 (4-6-4) 证明如下: 由于 所以有 (4-6-5) 上式说明了这样的事实:πu(精确的结点参数值构成的插值函数)在能量模的意义下并不比uh(由近似的结点参数值构成的插值函数)在整个求解域内更好(图4-16)。 3、收敛性证明 设真实解u∈HK(Ω),各单元的分片插值多项式完全到k-1次(k-1次称为有限元空间Sh的次数),有限元解为uh。 对于二阶问题,能量模最多同时包括函数及函数的一阶导数。由(4-4-1)和(4-4-2),对每个单元ei都有 (4-6-6) 代入(4-6-4),合并 , (4-6-7) 得到误差的能量模估计 在k≥2的情况下,(4-6-7)表明当h→0时 u→uh,并给出了收敛速度O(hk-1)。利用数学上的Nitsche技巧还可以分别得到函数和一阶导数的误差估计 (4-6-8) (4-6-9) 这两式表明,函数的精度高于导数的精度。 在u∈H1(Ω)的情况下,可选用完全一次多项式作为分片插值多项式。这时(4-6-6)-(4-6-9)不适用,但可通过其他方法证明,当h→0时‖u-uh‖→0,即 uh→u 。收敛性可以保证,但收敛速度可能较慢。 §4-7有限元(协调位移单元)解的特点 1、设u为真实解(u∈H1(Ω)),则总势能为 由对 δu∈H1(Ω)有 当取δu=u时也应该成立,故有 则 总势能 (4-7-1) 类似对于有限元解uh∈Sh也有 (4-7-2) 2、由 展开 由于uh∈Sh,利用投影关系有 (4-7-3) 得到 (4-7-4) 以及 (4-7-5) 协调位移单元以最小势能原理为基础。一般说来:有限元解的总势能高于真实解的总势能;有限元解的变形能小于真实解的变形能;外力在有限元解uh上做的功小于外力在真实解u上做的功。当外载荷仅是一个集中力(假如这种载荷是允许的话),那么力作用点的位移的有限元解一般比真实解小,这就是通常所说的,“有限元解的位移偏小、刚度偏大)的确切涵义。 §4-8结束语 本章借助于一些数学概念对协调位移单元进行了进一步讨论(插值逼近,收敛性,有限元解的特点)。(仅限于讨论二阶问题)。我们已经看到:有限元解的精度不但依赖于插值多项式的次数,而且依赖于广的可微性。还可以看到,对二阶问题,插值多项式必须至少是完全一次多项式,这正是第三章中提到的刚体位移和常应变条件。这一要求必须满足,因而是必要条件。 本章的结论对任何一种协调位移单元都适用,从理论角度讲,我们已经“相当彻底”地解决了问题。但到现在为止,就平面问题而言,我们仅构造出三角元和矩形元两种单元。换句话说,我们构造单元的实际能力相对说来还显得薄弱。这一问题将在下一章得到系统解决。
表中C代表一般常数,如果借助于分数阶索伯列夫空间,还可以对可微性做出更细致的估计。 扇形角φ及值 可能包含的奇异项 可能包含的奇异项
