有限单元法的数学基础 有限单元法的数学基础

有限单元法的数学基础 有限单元法的数学基础
1、引言 、 有限元方法归根结底是一种数值计算方法， 它有严格的数学证明作为其近似的客观性和 合理性的保证。力学问题最终归结为一组微分方程的边值问题或者初值问题抑或是混合问 题。比如弹性静力学最终归结为 L-N 方程的微分提法。在很难或者根本不可能得到所得方 程的理论解的情况下， 究竟用什么样的方法才能得到方程的近似解 （这种近似解已经能够满 足实际工程的需要） ，在这种情况下，二十世纪五六十年代由结构力学家进而由数学家提出 和证明了这种思想方法的合理性。有限元方法产生于力学计算，但是，它本质上并不是力学 的专利。 世间万物的变化过程很多都可以通过微分方程特别是偏微分方程来描述， 也就是说， 微分方程是很多现象和过程的数学结构， 而大多数的微分方程是不能得到理论解的， 这时候 就可以使用有限元方法来求其近似解，因为有限元方法是求解微分方程（组）的数值计算 有限元方法是求解微分方程（ 有限元方法是求解微分方程 方法。它适用于力学的微分方程，也同样适用于其它领域的相应的微分方程的数值求解。 方法 2、有限元方法数学根源 、 对于一个给定的微分方程定解问题， 为了求其近似解， 我们可以使用 Ritz 方法和 Galerkin 方法。下面分别阐述这两种方法，然后讨论有限元方法和他们的关系。 （1） Ritz 法 Ritz 法源于最小势能原理，设 H 是可分的 Hilbert 空间，在 H 中取有限维空间 Sn， 它是由 N 个线性无关向量 φ1 , φ2 , , φN 张成，即：
N S N ≡ ωn , ωn = ∑ Ciφi , (C1 , C2 , CN ) ∈ RN i =1

用 S N 代替 H，在 S N 上求泛函 J(w)的极值,即求 U N ∈ S N ，使得

J (U N ) = min ωN ∈S N J (ω N )
实际上寻求 U N 只需通过解一个线性方程组

1 J (ω ) = ( Dω , ω ) F (ω ) ≥ 0 2
D--------双线性形式 F--------线性泛函

ω N = ∑ Ciφi
i =1

N

J (ω N ) = =

N N N 1 D(∑ Ciφi , ∑ Ciφi )－F (∑ Ciφi ) 2 i =1 i =1 i =1 N 1 N ∑1 D(φi ,φ j )CiC j ∑ F (φi )Ci 2 i, j= i =1

因 此 ， J (ω N ) 是 一 个 以 C1 , C2 , , C N 为 未 知 数 （ 自 变 量 ） 的 二 次 多 项 式

j (C1 , C2 , , C N ) ，如果二次项的系数矩阵
[ D (φi , φ j )]i , j =1,2,, N 是正定的， 那么 j = j (C1 , C2 , , CN ) 在 N+1 维空间是一个
开口向上的椭球抛物面，它有且只有一个极（最）小值点，所谓在 S N 上求 J (ω N ) 的极 值，就是确定 C 1 , C 2 , , C
0 0 0 N

，使得：

j (C 01 , C 0 2 , , C 0 N ) ＝ min C1 ,,CN ∈R j (C 01 , C 0 2 , , C 0 N )
极值条件：

j | 0 0 （ i = 1, , N ） 0 ＝0 Ci C 1 ,C 2 ,,C N
得:
n

∑ D(φ φ )C
i =1 i j 0

0 i

= F (φi ) ( i = 1, , N )
0 0 N

即： C = [C 1 , C 2 , , C KC=F

]T 适合方程组：

F = [ F (φ1 ), , F (φ1 )]T

D(φ1 , φ1 )，D(φ2 , φ1 )， ，D(φN , φ1 ) D(φ1 , φ2 )，D(φ2 , φ2 )， ，D(φN , φ2 ) K= D(φ1 , φN ), D(φ2 , φN ), , D(φN , φN )
。。。。 。。。 总结：我们可以发现，Ritz 法主要体现在两点：第一，用有限维试解函数空间的基 函数生成真解的一个近似函数，这个试验函数的系数是基本未知量；第二，在求解过 程中进行了积分弱化，也就是说，由原来方程的高阶可微要求，通过一次积分和分步 积分，得到了弱解形式；第三，使用势能驻值原理。 （2） Galerkin 法 Galerkin 法源于虚位移原理， 虚位移原理的本质只是一个数学恒等式。 Galerkin 方法也是用 SN 代替 H,求 UN∈SN，使得对任意 WN∈SN 有 任意

D (U N , ω N ) = ( f , ω N )

且满足

U N = ∑ Ciφi
i =1

n

仍 然 只 需 要 解 一 个 代 数 方 程 组 ， 事 实 上 ， 设 ：

ω N = ∑ Ciφi
i =1

n

U N = ∑ C 0iφi , ω N = ∑ Ciφi
i =1 i =1

n

n

则：
N N

0=

i , j =1

∑ D(φi ,φ j )C 0 j Ci ∑ Ci F (φi )
i =1

＝

∑ [∑ D(φ , φ )C
i =1 j =1 i j

N

N

0 j

∑ F (φi )]Ci
i =1

N

由于 Ci = (i = 1, 2, , N ) 是任意的，必须有：

∑ D(φ ,φ )C
j =1 i j

N

0
j

= F (φi ), ( i = 1, 2, N )

因此得到与 Ritz 方法相同的代数方程组，此时 U 就是在 S N 取的 Galerkin 意义 下的广义解。 总结：我们可以发现，Galerkin 法主要体现在两点：第一，用有限维试解函数空 Galerkin 间的基函数生成真解的一个近似函数，这个试验函数的系数是基本未知量；第二，由 高阶可微变换成弱解积分形式；第三，使用虚位移原理。 （3） 有限元法 有限元方法和 Ritz 法、Galerkin 法有什么样的关系呢？ 通过上面的分析可以看出，除去由于方法不同可能带来的数值运算差别之外，仅从 原理上来说，Ritz 法和 Galerkin 法位移的区别就在于前者使用了势能驻值原理，而后者 使用了虚位移原理。而有限元法属于 Ritz—Galerkin 法的范畴，但又是改进，关键在于有 限元法选取分块多项式函数类作为子空间 SN。也就是说有限元方法可以从 Ritz 法导出， 也可以从 Galerkin 法导出 （即既可以基于势能驻值原理导出也可以基于虚位移原理导出） ， 但是唯一的区别是它的试解函数空间选取不同。它采用分片插值函数空间，这正好体现 了有限元剖分，而上面两种方法虽然也是近似，但是没有分片的概念。而这正是有限元 的核心。

具体说，就是先把Ω区间剖分为一系列单元，然后在每个单元内做多项式插值，而使他 们在单元的公共点、 面上满足一定的连续性条件， 边、 以保证这种函数组成的有限维空间 S N 是 H 的子空间。 所以，有限元的基础是变分原理和分片多项式插值。 下面，以二维 Poisson 方程的第一边值问题为例来说明有限元方法的这种思想。 Poisson 方程的边值问题：

△u = f ( x, y )
满足边值条件：

in

u = α ( x, y )
即

in



第一边值条件

△u = f ( x, y ) u = 0

in in
0



求广义解，就是求 U ∈ H 1 () ，使得：

D (u , ω ) = ( f , ω ), ω ∈ H 01 ()
其中： D (u , ω ) =

∫∫ (u w
x

x

+ u y wy )dxdy

( f , ω ) = ∫∫ f ω dxdy


对Ω进行三角剖分，记三角形单元为 e1 , e2 , , eNE ，Ω的内部节点为

P , P2 , , PNP ，边界节点为 PNP +1 , PNP + 2 , , PNP + M ，所有单元 ei 的最大边长是 h，对每一个 1
单元 e=△ PPj Pm 做线性插值，使它在三个顶点各取已知值 ωi , ω j , ωm ,当节点满足 上时， i 给定值为 0，这样在Ω上就构成了一个多面体函数 ωh ( x, y ) ，当 (ω1 , ω2 , , ωNP ) 遍历空间

R NP 所有向量时， 函数 ωh ( x, y ) 就构成了一个 N 维函数空间 SN 。 显然，SN 是一个线性空间，

且是 H 1 () （ C 0 () 在 意义下完备化） 。
0 1

SN 作 为 H 01 () 的 NP 维 线 性 子 空 间 ， 它 在 NP 个 基 函 数

φ1 ( x, y ), φ2 ( x, y ), , φNP ( x, y ) ，使得 ωn ( x, y ) 可以表示成：
NP

ωn ( x, y ) = ∑ ωiφi ( x, y )
i =1

有限元方法就是求 U1 , U 2 , , U NP ，使得：

U n ( x, y ) = ∑ U iφi ( x, y )
i =1

NP

满足： D (un , wn ) ( f , wn ) = 0, wn ∈ S N

设： ωn =
NP

∑ ω φ ( x, y) ，以 U ,ω 的表达式代入得：
i =1 i i

NP

n

n

i , j =1

∑ ωiu j D(φi ,φ j ) ∑ ωi ( f , φi ) = 0
i =1

NP

把 D (φi , φ j ) ,( f , φi )分别写成：

D(φi , φ j ) = ∑ ∫∫ (
n =1 en

NE

φi φ j φi φ j + )dxdy x x y y

( f , φi ) = ∑ ∫∫ f φi dxdy
n =1 en

NE

简记： Dn (φi , φ j ) =

∫∫ ( x
en

φi φ j φi φ j + )dxdy x y y

( f , φi ) n = ∫∫ f φi dxdy
en

设单元 en =△ Pk PPm ，因此在单元 en 上只有三个基函数 φk , φl , φm 不为零，其它所以基函 l 数在 en 上的值都为零。从而有：

i , j =1

∑ ω u D(φ ,φ ) = ∑ ω u ∑ D (φ , φ )
i j i j i , j =1 i j n =1 n i j

NP

NP

NE

= ∑ ∑ ωiu j Dn (φi , φ j )
n =1 i , j =1 NE

NE NP

=∑

n =1 i , j = k ,l , m

∑

ωi u j Dn (φi , φ j )

ωk D(φk , φk )，D(φk , φl )，D(φk , φm ) ui = ∑ ωl D (φl , φk )，D (φl , φl ) ，D (φl , φm ) u j n =1 ωm D(φm , φk ), D(φm , φl ) , D(φm , φm ) um
NE

由于在单元 en 上:

φk = N k , φl = N l , φm = N m
得：

D(un , ωn ) =

i , j =1

∑ ω u D(φ , φ ) = ∑ {δ *} [ K ] {δ }
T

NP

NE

i

j

i

j

n =1

en

en

en

其中：

{δ *}e

n

ωk = ωl ωm

{δ }e
位移向量

n

uk = ul um

虚伪移向量 单元上的刚度矩阵

[ K ]e

n

Dn ( N k , N k )，Dn ( N k , N l )，Dn ( N k , N m ) = Dn ( N l , N k )，Dn ( N l , N l ) ，Dn ( N l , N m ) D ( N , N ), D ( N , N ), D ( N , N ) n m k n m l n m m

同理： ( f , ωn ) =

∑ ωi ( f , φi ) = ∑ ωi ∑ (φ , φi )n =∑
i =1 i =1 n =1
T

NP

NP

NE

NE

n =1 i = k ,l , m

∑

ωi ( f , φi ) n

ωk = ∑ ωl n =1 ωm
NE

( f , φk )n ( f , φl ) n ( f , φm ) n

= ∑ {δ *}
n =1

NE

T

{F }e
en

n

{ F }e
NE T

:单元荷载向量

n

最终：

∑ {δ *} [ K ] {δ } ∑ {δ *} {F }
n =1 en en en n =1 en

NE

T

en

=0

意义：

1．

∫∫ (


φi φ j φi φ j + )dxdy x x y y

∫∫ φα fdxdy


Ω上的计算积分，转换到离散单元上计算，再叠加。 2．稀疏性：不相邻的单元， (φi , φ j ) ＝0 （4） 有限元法的其他解释 进一步上升和总结，可以将 Ritz 法、Galerkin 理解为加权残差法的特例。具体见《土 工计算机分析》 （龚晓楠 等 1999） 。 3、补充 、 实际上，有限元方法的产生应该是一件很自然的事情。因为，对于一个一般的定义于连 续区域上（物质的或空间的）微分方程，将其积分是一件很自然的事情，对其利用分部积分 也是一件很自然的事情， 然后对积分区域剖分也是很自然的， 因为根据高等数学的知识我们 就已经知道积分对于连续区域的可加性了。这些在探索其理论解的过程中都应该是常见的， 唯一的“不自然”的地方就是用插值函数来代替精确函数。因此，有限元方法最先由力学家 而不是数学家发现却是“不自然”的。