有限元方法

发表于 分类于 ，

又到了期末考试的时候，大概整理一下有限元的知识吧。

基础理论

有限元总的来说就是一种求偏微分方程数值解的方法。与简单地将时空微分算子变成差分算子不同，我们希望把解限定在一个较为一般的有限维空间中。这样解法就不必对时空区域有很高的要求了，而我们可以通过有限维空间基的线性组合算出比较合适的解。

弱解空间

对于经典解，其实我们很难将解限定在某一个特定的解空间内，而上学期学过的PDE课程恰好提供了适合的工具：弱解。 例如对于微分方程 \[ \begin{cases} -u''=f&,x\in \Omega = (0,1)\\ u(0)=u(1)=0& \end{cases} \] 我们只需要看它的变分形式，找它的弱解 \[ \begin{cases} Find&u\in H_0^1(\Omega),~s.t.\\ &(u',v') = (f,v),\forall v\in H_0^1(\Omega) \end{cases} \] 这样就自然找到了一个比较合适的解空间\(V=H_0^1(\Omega)\)。

有限元空间

正如上面所说，我们希望把解限定在某个有限维空间中。一般情况下，我们希望找的有限元空间正是弱解空间的一个有限维子空间\(V_h\subset V\)。

一旦有了这样的一个有限元空间，我们可以将问题转化为 \[ \begin{cases} Find&u_h\in V_h,~s.t.\\ &(u_h',v') = (f,v),\forall v\in V_h \end{cases} \] 既然是一个有限维的问题，通过找一组基函数，并且将所有函数表示为基函数线性组合的方式我们就能将问题转化为一个线性方程组。例如这里，假如有一组基\(V_h=span\{\phi_1,\cdots,\phi_n\}\) \[ u_h=\sum_{i=1}^n u_i\phi_i \] 则要条件对所有\(V_h\)中的\(v\)都成立等价于对每个基函数成立，即问题等价于 \[ \begin{cases} Find&U\in \mathbb{R}^n,~s.t.\\ &\sum_{j=1}^n u_j(\phi_j',\phi_i')=(f,\phi_i),~\forall i=1,2,\cdots,n \end{cases} \] 由于\(f,\phi_i\)都是已知的，这就转化为了一个线性方程组的问题\(KU=F\)，再通过求出的\(U\)以及取的基函数容易得到有限元空间中的解\(u_h\)。

本例中，一个基础且可行的有限元空间是连续分片线性空间。 \[ V_h=\{v|v\in C^0,~v|_{[x_{j-1},x_j]}\in P^1([x_{j-1},x_j]),~v(0)=v(1)=0\} \] 基取节点基函数即可。

至此，有限元方法的基本理论就结束了。只要找到弱解空间和有限元空间，对于更复杂情形的推广是比较平凡的。接下来的任务包括如何确定算法精度，以及如何对于各种特殊的空间找到合适的有限维子空间。

误差估计

FEM基本框架： \[ \begin{cases} Find&u_h\in V_h,~s.t.\\ &a(u_h,v) = F(v),\forall v\in V_h \end{cases} \] 其中\(a(\cdot,\cdot)\)是一个双线性形式，满足有界性和强制性 \[a(u,v)\le C\Vert u\Vert_V\Vert v\Vert_V\] \[a(v,v)\ge \alpha\Vert v\Vert_V^2\] 其中，记 \[\Vert\cdot\Vert_V^2=a(\cdot,\cdot)\] 这两个性质就是用来保证解的适定性。 我们想要知道有限元解和弱解究竟差了多少，那么也观察一下极其类似弱解形式： \[ \begin{cases} Find&u\in V,~s.t.\\ &a(u,v) = F(v),\forall v\in V \end{cases} \]

如果说这里我们假定有限元空间可以取在弱解空间之中，即\(V_h\subset V\)，则 \[a(u-u_h,v)=0,\forall v\in V_h\] 从而容易得到估计 \[\alpha \Vert u-u_h\Vert_V^2\le a(u-u_h,u-u_h)=a(u-u_h,u-v)\le C\Vert u-u_h\Vert_V\Vert u-v\Vert_V,\forall v\in V_h\] 即 (Cea Lemma) \[\Vert u-u_h\Vert_V\le\frac{C}{\alpha}\min_{v\in V_h}\Vert u-v\Vert_V\] 这意味着在这种范数意义下，\(u_h\)的估计离该空间中最佳估计也就最多差一个常数倍了，所以已经是一个比较好的逼近解了。

当然这个范数本身可能并不是我们想要的，如果想要估计误差的\(L^2\)范数，可能就要动用一种对偶的小技巧(duality)

\[\Vert u-u_h\Vert_{L^2}=(u-u_h,u-u_h):=a(u-u_h,w)\] 通过变分形式\(a(v,w)=(v,u-u_h),\forall v\)来待定\(w\)。

例如在之前的问题中\(a(u,v)=(u',v')\)则取对偶问题 \[ \begin{cases} -\Delta w=u-u_h\\ w|_{\partial \Omega}=0 \end{cases} \] 如此有 \[ \begin{aligned} \Vert u-u_h\Vert_{L^2}=a(u-u_h,w)=a(u-u_h,w-w_I)&\le C\Vert u-u_h\Vert_V\Vert w-w_I\Vert_V\\ &\le Ch\Vert u-u_h\Vert_V|w|_{H^2}\\ &\le Ch\Vert u-u_h\Vert_V\Vert u-u_h\Vert_{L^2} \end{aligned} \] 不等式分别由逼近性质和对偶问题解的正则性估计得到。 \[ \Vert u-u_h\Vert_{L^2}\le Ch\Vert u-u_h\Vert_V \] 具有比\(\Vert\cdot\Vert_V\)更优一阶的逼近。

逼近理论

在误差估计中，最基础的估计就是Cea Lemma \[\Vert u-u_h\Vert_V\le\frac{C}{\alpha}\min_{v\in V_h}\Vert u-v\Vert_V\] 它将我们有限元解的误差直接缩小到了空间内最佳逼近的常数倍范围内，但是想要知晓具体的误差量级，还是得要考察一下空间内的最佳逼近误差是在什么量级的。

基本出发点：一般函数可以用多项式函数有效逼近

一般来说为了更轻易地控制收敛阶数，我们都会选取有限元空间为分片多项式空间，那么最佳逼近误差就是希望在每个块内找最佳逼近的多项式。最好的已有理论自然是Taylor多项式逼近 \[ T_y^mu(x)=\sum_{|\alpha|<m}\frac1{\alpha}D^{\alpha}u(y)(x-y)^\alpha \] 然而这里我们讨论的\(u\)一般是在一个Sobelev空间中的，\(D^\alpha u\)可能并不能在逐点情况下有意义。这里我们对其做一个卷积光滑化

令 \[ \psi(x)= \begin{cases} e^{-(1-(|x-x_0|/\rho)^2)^{-1}}&,|x-x_0|<\rho\\ 0,&,|x-x_0|\ge \rho \end{cases} \] \[ c=\int_{\mathbb{R}^n}\psi(x)dx \] 则\(\phi(x)=(1/c)\psi(x)\)满足

  1. \(\phi\in C_0^\infty(\mathbb{R}^n)\)
  2. \(supp~\phi=\overline{B_\rho(x_0)}\)
  3. \(\int_{\mathbb{R}^n}\phi=1\)
  4. \(\max |\phi|\le C\cdot \rho^{-n}\)

这其实就是截断函数的一个简单构造。

令\(B=B_\rho(x_0)\)，光滑化的Taylor多项式为 \[ \begin{aligned} Q^mu(x)&=\int_BT_y^mu(x)\phi(y)dy\\ &=\int_B\sum_{|\alpha|<m}\frac1{\alpha}D^{\alpha}u(y)(x-y)^\alpha\phi(y)dy \end{aligned} \]

这样对于\(u\in W_p^{m-1}\)就是有意义的了。容易看出对于\(x\)，这是一个至多\(m-1\)次多项式。

关于\(x\)整理一下可得 \[ Q^mu(x)=\sum_{|\alpha|<m}x^{\lambda}\int_B{\psi_\lambda}(y)u(y)dy \] 其中 \[ supp(\psi_\lambda)\subset\overline{B} \] \[ \psi_\lambda\in C_0^\infty(\mathbb{R}^n) \] 于是自然有估计 \[ \Vert Q^mu\Vert_{W_{\infty}^k(\Omega)}\le C_{m,n,\rho,\Omega}\Vert u\Vert_{L^1(B)} \]

经过一些比较复杂的理论推导，可以推出以下比较重要的结论

（以下记\(d=diam(\Omega),\rho_{max}\)为区域内包含的圆的最大半径，\(\gamma=\frac{d}{\rho_{max}}\)）

余项估计

对于 \(u\in W_p^m(\Omega)\) \[ \Vert R^mu\Vert_{L^\infty(\Omega)}\le C_{m,n,\gamma,p}d^{m-n/p}|u|_{W_p^m(\Omega)} \]

Sobelev不等式

当区域\(\Omega\)关于球\(B\)为星形域时，若满足条件之一 1. \(1<p<\infty,m>n/p\) 2. \(p=1,m\ge n\)

则\(u\)连续且 \[ \Vert u\Vert_{L^\infty(\Omega)}\le C_{m,n,\gamma,d}\Vert u\Vert_{W_p^m(\Omega)} \]

Bramble-Hilbert引理

当\(\Omega\)是star-shaped区域，\(u\in W_p^m(\Omega),p\ge 1\)时 \[ |u-Q^mu|_{W_p^k(\Omega)}\le C_{m,n,\gamma}d^{m-k}|u|_{W_p^m(\Omega)} \]

结尾

准备工作大致结束了，更具体的问题下次再写吧。