多重网格法

前言

从个人本学期进行的几次实验情况来看,有限元方法耗时最长的几个地方:一是网格划分,需要尽量均匀,甚至最好让外接圆半径小、内切圆半径大,不过我们大部分实验都是一维情况的,仅有的二维实验还给出了网格划分,所以暂且不做考虑;二是数值积分,计算那些\(a(\phi_i, \phi_j)\)\((f,\phi_i)\)的值来组装矩阵,一般来说用比较简单的Gauss积分公式就能够满足精度要求,程序写起来会有点麻烦,但是不会很耗时;最后就是求解线性方程组,在大部分时候这一步成为了限制我的程序性能的枷锁。

这里介绍的多重网格方法就是适用于有限元问题中的一种求解线性方程组的方法,它可以在\(O(n)\)的时间找到一个解,并且达到和原来同样量级的精度!

多重网格法简介

多重网格法顾名思义,要准备一系列的网格。有点类似于多分辨率分析的思路,我们的想法就是把高精度的问题投影到低精度的网格上求解,再将解投影回高精度网格做一些光滑化的操作。

在有限元问题中,要不断提高精度,我们本来就需要准备一系列网格以及有对应的问题 \[ A_k u_k = f_k \] 所以我们只需要再准备投影函数 \[I_k^{k-1},I_{k-1}^k\] 满足 \[ (I_k^{k-1}w,v)_{k-1}=(w,I_{k-1}^kv)_k=(w,v)_k, \forall v\in V_{k-1},w\in V_k \]

算法

先列出一个投影求解的子过程\(MG(k,z_0,g)\),这一步中求解子问题\(A_k z=g\)

  1. \(k=1\)则最底层网络直接求解(通过共轭梯度等常用方程组解法)

  2. 否则,求解过程分为几步:

    1. Presmoothing Step \[ z_l=z_{l-1}+\frac{1}{\Lambda_K}(g-A_kz_{l-1}),1\le l\le m_1 \]

    2. Error Correction Step

      首先将上一步最后迭代得到的结果误差投影到上一层 \[ \bar g=I_k^{k-1}(g-A_k z_{m_1}) \]

      然后上一层中迭代求解一个合适的调整(希望是\(A_k q=g-A_kz_{m_1}\),实际求解\(A_{k-1} q=\bar g\)),这一步中调用子过程p次。

      \[ q_i =MG(k-1,q_{i-1},\bar g) \]

      最终将调整量投影回这一层 \[ z_{m_1+1}=z_{m_1}+I_{k-1}^kq_p \]

    3. Postsmoothing Step

      \[ z_l=z_{l-1}+\frac{1}{\Lambda_K}(g-A_kz_{l-1}),m_1+2\le l\le m_1+m_2+1 \]

    \(MG(k,z_0,g)\)最终结果即为最后的\(z_{m_1+m_2+1}\)

上面的过程中光滑化是用的Richardson迭代,事实上参数的选取还比较麻烦,可能直接用Jacobi迭代或Gauss-Seidel迭代会是更好的可行选择吧。

\(p\)是一个非常重要的参数,一般常用的就是\(p=1\)(V-cycle)或\(p=2\)(W-cycle)。 事实上相邻两层的数据规模比大致是\(2^d\),光滑化过程在稀疏矩阵条件下大致是\(O(n)\)复杂度。那么相当于有条件

\[ T(n)=pT(\frac{n}{2^d})+O(n) \]

\(p<2^d\)条件下为线性的时间复杂度!

方法精度

证明过程还是有点复杂的,但是确实可以说明

\[ \Vert u_k-\tilde u_k\Vert_E \le Ch_k|u|_{H^2(\Omega)} \]

而有限元解误差我们给出的估计也是

\[ \Vert u-u_k\Vert_E \le Ch_k|u|_{H^2(\Omega)} \]

所以这种方法求解并没有损失精度。

o( ̄▽ ̄)d

结尾

其实后面还学了不少内容,包括混合元、间断元等等,但据说不做主要考查要求,而且都是大篇幅的理论推导,与我本来只是想整理一下思路的初衷相悖,所以有限元理论就暂且写到这吧。