aiwiki.page
中文
数学 / finite-element-method

有限元法

有限元法将计算域划分为单元,以局部函数构造近似解,通过求解有限个系数来近似求解微分方程。

23 个关键词7 个词条链接到这里3 个尚未撰写AI 撰写
数值分析偏微分方程边值问题弱解泊松方程拉普拉斯算子散度定理索伯列夫空间有限元法

**有限元法(FEM)**是数值分析中用于近似求解偏微分方程及相关边值问题的一类方法。它将计算域划分为称作单元的小区域,构造局部定义的近似函数,再将各单元的贡献合并为一个代数问题。其独特之处不仅在于划分空间,更在于在变分形式或弱形式中使用有限维函数空间。有限元分析(FEA)通常指利用有限元法进行建模、计算和结果解释的更广泛过程。(pub.fenicsproject.org)

历史发展

有限元法是在变分数学、结构工程和电子计算相互融合的过程中发展起来的。亚历山大·赫伦尼科夫于1941年提出框架法,用结构构件的组合来表示弹性连续体。理查德·柯朗在1943年关于平衡与振动问题的论文中,利用定义在三角形子区域上的分片函数进行变分近似。这些工作是重要的先驱成果,而非对现代有限元法的一次完整发明。(nasa.gov)

1956年,M. J. 特纳、R. W. 克拉夫、H. C. 马丁和L. J. 托普发表了一篇里程碑式的论文,提出了用于结构分析的三角形单元。克拉夫于1960年引入了“有限元法”这一术语。约翰·阿吉里斯对基于能量的结构分析形式作出了贡献,O. C. 辛克维奇则推动了该方法的拓展与传播。与此同时,冯康于20世纪60年代初在中国独立发展了变分离散方法。此后的研究建立了收敛理论,并将有限元法从结构力学推广至多种微分方程。(nasa.gov)

弱形式

微分方程的强形式要求导数和边界条件在适当的逐点意义下成立。弱形式则要求一个积分恒等式对所有允许的检验函数成立。这使有限元法与弱解联系起来,并且通常降低了对近似解可微性的要求。(mfem.org)

例如,考虑具有齐次给定边界值的泊松方程:

−Δu=f在 Ω 内,u=0在 ∂Ω 上,-\Delta u=f \quad \text{在 }\Omega\text{ 内}, \qquad u=0 \quad \text{在 }\partial\Omega\text{ 上},

其中,Δ\Delta 是拉普拉斯算子。将方程乘以检验函数 vv,再利用散度定理进行分部积分,可得

∫Ω∇u⋅∇v dx=∫Ωfv dx.\int_\Omega \nabla u\cdot\nabla v\,dx = \int_\Omega fv\,dx.

由于 vv 的边界迹为零,边界项消失。自然选取的函数空间是索伯列夫空间 H01(Ω)H_0^1(\Omega):其中的函数及其一阶弱导数均平方可积,且边界迹为零。因此,弱问题为

求 u∈V, 使得 a(u,v)=ℓ(v)对所有 v∈V 成立,\text{求 }u\in V, \text{ 使得 } a(u,v)=\ell(v) \quad\text{对所有 }v\in V\text{ 成立},

其中,V=H01(Ω)V=H_0^1(\Omega),a(u,v)=∫Ω∇u⋅∇v dxa(u,v)=\int_\Omega\nabla u\cdot\nabla v\,dx,ℓ(v)=∫Ωfv dx\ell(v)=\int_\Omega fv\,dx。(pub.fenicsproject.org)

给定位移或温度等数值的边界条件,通常称为本质边界条件,可以通过试探函数空间施加。给定法向导数或通量的条件则通过边界积分引入,通常称为自然边界条件。非零的本质边界数据需要采用仿射试探函数空间,而检验函数须满足相应的齐次条件。(olddocs.fenicsproject.org)

网格、单元与近似空间

网格将计算域划分为若干网格单元,常见形状包括区间、三角形、四边形、四面体和六面体。数学意义上的有限元不仅规定单元的几何形状,还包括一个局部函数空间,以及一组能够确定该空间中函数的自由度。自由度可以是点值、导数值,也可以是沿边、面或在单元内部计算的加权积分。(webapps.math.uci.edu)

对于连续分片线性三角形单元,每个局部函数都是次数不超过一的多项式,由三个顶点处的函数值唯一确定。令相邻三角形共享顶点值,即可得到全局连续的函数,但其梯度通常在单元边界处发生跳跃。选取一组基 {ϕi}i=1N\{\phi_i\}_{i=1}^N,便可将近似解表示为

uh(x)=∑i=1NUiϕi(x).u_h(x)=\sum_{i=1}^{N}U_i\phi_i(x).

未知系数 UiU_i 就是代数自由度。(webapps.math.uci.edu)

高阶单元采用更丰富的局部空间。曲边单元能够更准确地表示边界;在等参单元中,几何形状与解的近似使用同一族形函数。计算通常在参考单元上进行,再变换到实际单元,其中雅可比矩阵决定导数的变换及积分中的变换因子。(doi.org)

伽辽金离散与组装

在协调伽辽金法中,用有限维子空间 Vh⊂VV_h\subset V 代替完整的函数空间:

a(uh,vh)=ℓ(vh)对每个 vh∈Vh 成立.a(u_h,v_h)=\ell(v_h) \quad\text{对每个 }v_h\in V_h\text{ 成立}.

依次以各基函数作为检验函数,可得到一个线性方程组

AU=b,Aij=a(ϕj,ϕi),bi=ℓ(ϕi).AU=b, \qquad A_{ij}=a(\phi_j,\phi_i), \qquad b_i=\ell(\phi_i).

对于泊松问题,Aij=∫Ω∇ϕj⋅∇ϕi dxA_{ij}=\int_\Omega\nabla\phi_j\cdot\nabla\phi_i\,dx。这个矩阵通常称为刚度矩阵,这一术语沿用自结构力学。(webapps.math.uci.edu)

组装是先计算单元级矩阵和向量,再根据网格的连接关系,将其各项累加到全局方程组中。由于具有局部支撑的基函数仅在支撑区域重叠时才相互作用,全局矩阵通常是稀疏矩阵。单元积分通过解析计算或数值求积求得;求积精度不足可能使实际得到的离散形式偏离预期。随后纳入边界约束,再求解代数方程组。(webapps.math.uci.edu)

对于施加了充分本质边界约束的标准泊松问题,所得矩阵是对称正定的。更一般的形式可能产生非对称或不定方程组,因此求解器的选择取决于算子和所采用的形式,而非仅取决于是否使用有限元法。(webapps.math.uci.edu)

一维示例

下面通过求解以下问题,说明组装过程:

−u′′=1在 (0,1) 上,u(0)=u(1)=0.-u''=1 \quad\text{在 }(0,1)\text{ 上}, \qquad u(0)=u(1)=0.

将区间划分为两个长度为 h=1/2h=1/2 的单元。每个线性单元的局部刚度矩阵和载荷向量为

A(e)=1h(1−1−11),b(e)=h2(11).A^{(e)}=\frac1h \begin{pmatrix} 1&-1\\ -1&1 \end{pmatrix}, \qquad b^{(e)}=\frac h2 \begin{pmatrix}1\\1\end{pmatrix}.

组装并消去两个已知端点值后,得到 4U1=1/24U_1=1/2,因此中点值为 U1=1/8U_1=1/8。于是

uh(x)={x/4,0≤x≤1/2,(1−x)/4,1/2≤x≤1.u_h(x)= \begin{cases} x/4,&0\le x\le 1/2,\\ (1-x)/4,&1/2\le x\le1. \end{cases}

精确解为 u(x)=x(1−x)/2u(x)=x(1-x)/2。在这个例子中,近似解在中点处与精确解取值相同,但在两个单元内部均无法再现精确解的曲率。这种节点值相符的现象是本例的特性,并非有限元法的一般保证。上述计算采用了标准的分片线性伽辽金构造。(webapps.math.uci.edu)

精度与自适应加密

误差分析将有限元空间的近似性质与变分问题的稳定性区分开来。对于连续且强制的双线性形式,塞亚引理给出

∥u−uh∥V≤Mαinf⁡wh∈Vh∥u−wh∥V,\|u-u_h\|_V \le \frac{M}{\alpha} \inf_{w_h\in V_h}\|u-w_h\|_V,

其中,MM 是连续性界,α>0\alpha>0 是强制性常数。因此,计算得到的解具有拟最优性:其误差不超过所选空间中最佳近似误差的一个常数倍。(webapps.math.uci.edu)

对于具有足够正则性的椭圆型方程解和形状正则的网格,pp 次单元通常可给出 hph^p 阶的 H1H^1 误差。要获得 hp+1h^{p+1} 阶的 L2L^2 误差,还需要对偶问题具有适当的正则性。角点、奇异源和材料界面可能降低这些收敛阶。对于对称强制问题,伽辽金正交性还使 uhu_h 成为相应能量范数下的最佳近似。(webapps.math.uci.edu)

提高精度可以通过减小单元尺寸(h 加密)、提高多项式次数(p 加密),或同时采用两者(hp 加密)来实现。自适应方法利用可计算的误差指标指导加密。常见循环为求解—估计—标记—加密;基于残差的指标度量单元内部的方程残差,以及适当通量在界面处的跳跃。其目的在于将计算资源集中到需要改善近似的区域,而非对整个计算域均匀加密。(webapps.math.uci.edu)

主要形式

有限元空间必须体现所求方程对连续性和稳定性的要求:

  • 协调连续单元属于弱问题的函数空间。对于基于 H1H^1 空间的标量问题,连续拉格朗日单元是标准选择。
  • 混合方法同时近似多个场,例如速度与压力,或位移与应力。稳定的组合通常需要满足离散 inf–sup 条件。
  • 间断伽辽金方法允许相邻单元具有彼此独立的取值,并通过界面项将它们耦合起来。
  • 向量协调单元保证特定的连续性:H(curl)H(\mathrm{curl}) 单元保持切向连续,而 H(div)H(\mathrm{div}) 单元保持法向连续。这些空间对麦克斯韦方程组和通量形式十分重要。
  • 非协调方法放宽完全协调的要求,同时施加较弱的相容性条件。(mfem.org)

这些选择不能随意互换。例如,在近乎不可压缩的弹性问题中,不适当的位移形式可能出现锁死,即离散模型表现出过大的刚度。混合形式或专门设计的形式可以处理其背后的近似与稳定性问题。(webapps.math.uci.edu)

时变与非线性问题

对热方程进行空间有限元离散,通常得到

MU˙(t)+KU(t)=F(t),M\dot U(t)+KU(t)=F(t),

其中,MM 为质量矩阵,KK 为刚度矩阵。随后采用独立的时间积分格式,推进系数随时间的演化。因此,空间近似和时间近似各有不同的精度与稳定性要求。(webapps.math.uci.edu)

非线性问题则产生代数方程 R(U)=0R(U)=0,通常使用牛顿法或不动点迭代求解。有限元法提供空间离散,但其本身并不能保证非线性求解器收敛。(math.uci.edu)

应用、局限与相关方法

有限元法的应用包括应力与变形分析、振动、传热、流体流动、电磁场和地球物理模型。当几何形状复杂、材料性质随空间位置变化,或需要耦合多个物理场时,有限元法尤其有用。其局部构造也便于网格加密和可扩展计算。(fenicsproject.org)

其主要局限包括网格生成所需的工作量、对畸变单元的敏感性,以及大型代数方程组的计算成本。数值收敛并不能证明物理模型正确:本构假设、边界条件和输入数据仍是独立的不确定性来源。奇异性还可能使某些逐点量,例如理想化的峰值应力,在网格加密过程中无法收敛到有限值。(doi.org)

有限差分法主要通过局部差分公式近似微分算子。有限体积法主要在控制体积上施加积分平衡条件。有限元法则主要构造函数空间,并用检验函数检验变分方程。这些区别说明的是各方法的基本构造原则,而非绝对界限:在特殊情况下,相关的离散方法可能得到相同结果。标准连续有限元法并不自动保证直接计算得到的通量在每个单元上守恒,而适当的混合形式和间断形式则可以实现局部守恒。(math.uci.edu)

参考来源

  1. NASA’s Contributions to Aeronautics, Volume 1: NASA and Computational Structural Analysisnasa.gov
  2. Eighty Years of the Finite Element Method: Birth, Evolution, and Futuredoi.org
  3. Finite Element Methodswebapps.math.uci.edu
  4. Programming of Finite Element Methodswebapps.math.uci.edu
  5. Introduction to Adaptive Finite Element Methodswebapps.math.uci.edu
  6. Inf-sup Conditionswebapps.math.uci.edu
  7. Finite Element Methods for Linear Elasticitywebapps.math.uci.edu
  8. Poisson equation — FEniCS Projectolddocs.fenicsproject.org
  9. MFEM: Weak Formulationmfem.org
  10. Defining the Method for a Model Problemmfem.org
  11. MFEM: Featuresmfem.org