FVM in CFD 學習筆記_第10章_補充專題_多重網格算法

參看資料：
https://ocw.mit.edu/courses/mathematics/18-086-mathematical-methods-for-engineers-ii-spring-2006/readings/am63.pdf
https://computing.llnl.gov/projects/hypre-scalable-linear-solvers-multigrid-methods/CiSE_2006_amg_220851.pdf
user.it.uu.se/~maya/Projects/3phase/AMG_parallel_Falgout.pdf

因爲書中所講的多重網格算法過於簡略，沒法參悟透徹，故重開一篇，再講講如何使用該算法。看了網上不少資料，大多語焉不詳，要麼是大部頭的理論剖析，要麼是寥寥數語不知如何應用，確實讓人稀裏糊塗沒法弄得很透徹。本篇純屬做者我的胡蒙亂猜的，不當錯誤之處在所不免，還請及時指正，省得以訛傳訛，貽誤後人。閒言少敘，開吹：

在傳統的迭代方法中，迭代是在固定網格上進行的，這樣就會使得偏差中的高頻份量（波長較短的部分）衰減較快，而低頻份量（波長較長）的部分衰減較慢，而多重網格的思想是讓粗細不一樣網格層級之間來回計算折騰，從而使得各類高低頻率長短波長的偏差份量都獲得快速地衰減，那麼迭代收斂速度比傳統的迭代方法要快得多了。

多重網格算法分爲幾何多重網格和代數多重網格，幾何多重網格是真的劃分了網格，真的構造了粗細網格層級；而代數多重網格則是直接跳脫出了網格的束縛，僅僅根據得出的係數矩陣 $\bold A$，運用多重網格的思想，來進行虛擬的網格層級處理，從而提升收斂速度。

顯然，幾何多重網格更容易理解，也更簡單一些，而在幾何多重網格中，結構化的網格處理起來則更爲方便一些，那麼我們就先從結構網格開始吧。

1 結構網格上的多重網格算法

因爲要把粗網格和細網格的信息相互傳遞轉化，因此實際上須要3個矩陣來作這些事情（用下標 $h$來標識細網格，用下標 $2h$來標識粗網格）：

• 把矢量從細網格到粗網格傳遞的限制矩陣 $\bold R_h^{2h}$
• 反過來，把矢量從粗網格向細網格延拓的插值矩陣 $\bold I_{2h}^{h}$
• 固然了，還要根據細網格的係數矩陣 $\bold A_h$，來構造粗網格上的係數矩陣 $\bold A_{2h}$，那麼這個矩陣實際上就是 $\bold A_{2h}=\bold R_h^{2h} \bold A_h \bold I_{2h}^{h}$

1.1 插值矩陣 $\bold I_{2h}^{h}$

先來看插值矩陣 $\bold I_{2h}^{h}$，對於1維問題，結構網格，計算域 $0\le x \le1$，邊界值爲0，用 $h=1/8$來剖分細網格，用 $2h=1/4$來剖分粗網格，邊界點的值不參與編碼和計算，則細網格上的變量 $\bold u=[u_1,u_2,u_3,u_4,u_5,u_6,u_7]^T$共7個內部節點，而粗網格上的變量 $\bold v=[v_1, v_2, v_3]^T$共3個內部節點。

如何根據粗網格（ $2h$）上的變量 $\bold v$來插值算得細網格（ $h$）上的變量 $\bold u$呢？不難發現 $u_2,u_4,u_6$是和 $v_1,v_2,v_3$重合的節點上的變量值，因此直接讓它們相等就行了。而對於 $u_1,u_3,u_5,u_7$這些點處的變量值，它們的位置是位於兩個粗單元節點之間的（且是中點），因此用其兩側的粗單元節點值的算術平均來處理就行了（注意，邊界節點的值是0），即 $u_1=v_1/2$， $u_3=v_1/2+v_2/2$， $u_5=v_2/2+v_3/2$， $u_7=v_3/2$。寫成矩陣形式，即
$\begin{bmatrix} u_1 \\ u_2 \\ u_3 \\ u_4 \\ u_5 \\ u_6 \\ u_7 \end{bmatrix} =\frac{1}{2} \begin{bmatrix} 1 & & \\ 2 & & \\ 1 & 1 & \\ & 2 & \\ & 1 & 1 \\ & & 2 \\ & & 1 \\ \end{bmatrix} \begin{bmatrix} v_1 \\ v_2 \\ v_3 \end{bmatrix}$
從而得到了插值矩陣 $\bold I_{2h}^{h}$
$\bold I_{2h}^{h}=\frac{1}{2} \begin{bmatrix} 1 & & \\ 2 & & \\ 1 & 1 & \\ & 2 & \\ & 1 & 1 \\ & & 2 \\ & & 1 \\ \end{bmatrix}$

1.2 限制矩陣 $\bold R_h^{2h}$

若是知道了細網格上的量，又如何將其限制到粗網格上呢？一種簡單的想法是讓相同節點的粗細網格值相互傳遞，即， $[v_1,v_2,v_3]=[u_2,u_4,u_6]$，這樣作雖然看起來蠻好的，但問題在於，沒有充分利用到細網格上面全部節點的信息，即，沒有考慮到 $[u_1,u_3,u_5,u_7]$節點值的影響。經常使用的處理方法是根據以前的插值矩陣 $\bold I_{2h}^{h}$來作反向處理，即讓 $\bold R_h^{2h}=\frac{1}{2}(\bold I_{2h}^{h})^T$，稱爲全加權算子 $\bold R_h^{2h}$，如此一來
$\begin{bmatrix} v_1 \\ v_2 \\ v_3 \end{bmatrix}=\frac{1}{4} \begin{bmatrix} 1 & 2 & 1 & & & & \\ & &1 & 2 & 1 & & \\ & & & &1 & 2 & 1 \end{bmatrix} \begin{bmatrix} u_1 \\ u_2 \\ u_3 \\ u_4 \\ u_5 \\ u_6 \\ u_7 \end{bmatrix}$
寫成份量形式
$v_i = \frac{1}{4}[u_{2i-1}+2u_{2i}+u_{2i+1}],\quad i=1,2,3$
即每一個粗網格節點的值，由與其重合的細網格節點值和其鄰近的細網格節點值加權獲得。全加權算子 $\bold R_h^{2h}$爲
$\bold R_h^{2h}=\frac{1}{4} \begin{bmatrix} 1 & 2 & 1 & & & & \\ & &1 & 2 & 1 & & \\ & & & &1 & 2 & 1 \end{bmatrix}$

1.3 二維問題的插值矩陣和限制矩陣

仍是結構網格，若是問題擴展爲二維，即笛卡爾網格，那麼這時候插值矩陣和限制矩陣是啥樣的呢？

用 $u_{i,j}$和 $v_{i,j}$分別表示密網格和粗網格上的變量，先看從粗網格到細網格轉化變量的插值矩陣，跟一維問題相似，若是細網格節點和粗網格節點重合，那麼直接讓二者存儲的變量相同，即
$u_{2i,2j}=v_{i,j}$
若是細網格節點恰好落在粗網格兩節點之間，就用粗網格這倆節點上變量平均來計算，即
$u_{2i+1,2j}=\frac{1}{2}(v_{i,j}+v_{i+1,j}) \quad u_{2i-1,2j}=\frac{1}{2}(v_{i,j}+v_{i-1,j}) \\ ~\\ u_{2i,2j+1}=\frac{1}{2}(v_{i,j}+v_{i,j+1}) \quad u_{2i,2j-1}=\frac{1}{2}(v_{i,j}+v_{i,j-1})$
若是細網格節點是位於四個粗網格節點之間的，那麼用這四個粗網格節點上變量的平均來計算，即
$u_{2i+1,2j+1}=\frac{1}{4}(v_{i,j}+v_{i+1,j}+v_{i,j+1}+v_{i+1,j+1}) \\ ~\\ u_{2i+1,2j-1}=\frac{1}{4}(v_{i,j}+v_{i+1,j}+v_{i,j-1}+v_{i+1,j-1}) \\ ~\\ u_{2i-1,2j+1}=\frac{1}{4}(v_{i,j}+v_{i-1,j}+v_{i,j+1}+v_{i-1,j+1}) \\ ~\\ u_{2i-1,2j-1}=\frac{1}{4}(v_{i,j}+v_{i-1,j}+v_{i,j-1}+v_{i-1,j-1}) \\$
如此，即可推得插值矩陣 $\bold I_{2h}^{h}$。

二維問題的限制矩陣 $\bold R_h^{2h}=\frac{1}{4}(\bold I_{2h}^{h})^T$，粗網格的每一個節點值由細網格上的9個節點值插值獲得，各節點係數如圖所示。

即
$v_{i,j}=\frac{1}{16}[4u_{2i,2j}+2(u_{2i+1,2j}+u_{2i-1,2j}+u_{2i,2j+1}+u_{2i,2j-1})+(u_{2i+1,2j+1}+u_{2i+1,2j-1}+u_{2i-1,2j+1}+u_{2i-1,2j-1})]$
這樣，也可推得限制矩陣 $\bold R_h^{2h}$。

對於三維問題，也可用相似方法推出插值矩陣和限制矩陣。

1.4 粗網格係數矩陣 $\bold A_2h$

粗網格上的係數矩陣 $\bold A_{2h}$是由細網格上的係數矩陣 $\bold A_h$和變量傳遞關係轉化而來的，故有 $\bold A_{2h}=\bold R_h^{2h} \bold A_h \bold I_{2h}^{h}$，完美。

好比，若是一維問題的控制方程和零邊值條件爲
$-\frac{d^2\phi}{dx^2}=f(x),\quad 0\le x\le 1\\ \phi(x=0)=\phi(x=1)=0$
跟前面講的網格剖分方式同樣，細網格等分8個單元，單元長度爲 $h=1/8$，用中心格式近似2階項，有
$-\frac{d^2\phi}{dx^2}=-\frac{\phi_{i+1}-2\phi_{i}+\phi_{i-1}}{h^2}$
則細網格上係數矩陣 $\bold A_h$爲
$\bold A_h=\frac{1}{h^2}\begin{bmatrix} 2 & -1 \\ -1 & 2 & -1 \\ & -1 & 2 & -1 \\ &&-1 & 2 & -1 \\ &&&-1 & -2 & -1 \\ &&&&-1 & 2 & -1 \\ &&&&&-1 & 2 \\ \end{bmatrix}$
粗網格等分4個單元，則粗網格上係數矩陣爲
$\begin{aligned} \bold A_{2h}&=\bold R_h^{2h} \bold A_h \bold I_{2h}^{h} \\ &=\frac{1}{4} \begin{bmatrix} 1 & 2 & 1 & & & & \\ & &1 & 2 & 1 & & \\ & & & &1 & 2 & 1 \end{bmatrix} \frac{1}{h^2}\begin{bmatrix} 2 & -1 \\ -1 & 2 & -1 \\ & -1 & 2 & -1 \\ &&-1 & 2 & -1 \\ &&&-1 & -2 & -1 \\ &&&&-1 & 2 & -1 \\ &&&&&-1 & 2 \\ \end{bmatrix}\frac{1}{2} \begin{bmatrix} 1 & & \\ 2 & & \\ 1 & 1 & \\ & 2 & \\ & 1 & 1 \\ & & 2 \\ & & 1 \\ \end{bmatrix} \\ &=\frac{1}{(2h)^2} \begin{bmatrix} 2& -1 \\ -1 & 2 & -1 \\ & -1 & 2 \end{bmatrix} \end{aligned}$

1.5 多重網格算法流程

若是僅有兩層網格，細網格-粗網格-細網格這樣作完一個循環，就是兩層網格的V型循環，其算法流程以下：

1. 在細網格上迭代求解 $\bold A_h \bold u = \bold b_h$，用Jacobi或Gauss-Seidel迭代，注意只作3次迭代步，獲得還沒有收斂的解 $\bold u_h$（雖然沒收斂，可是高頻份量衰減了不少）；
2. 計算細網格上面的殘差 $\bold r_h=\bold b_h - \bold A_h \bold u_h$，並將 $\bold r_h$傳遞到粗網格上面去，用 $\bold r_{2h}=\bold R_h^{2h}\bold r_h$獲得粗網格上面的殘差 $\bold r_{2h}$；
3. 計算粗網格上係數矩陣 $\bold A_{2h}=\bold R_h^{2h} \bold A_h \bold I_{2h}^{h}$（這個能夠放到最前面算好存好拿來直接用，不用每次都算的），求解粗網格上方程 $\bold A_{2h} \bold E_{2h}=\bold r_{2h}$，一樣用迭代方法來求解，仍是迭代3步就好，迭代初值選擇爲 $\bold E = \bold 0$，迭代完成後算出了修正值 $\bold E_{2h}$；
4. 將粗網格上的修正值 $\bold E_{2h}$插值回細網格，用 $\bold E_{h}=\bold I_{2h}^{h}\bold E_{2h}$，得到細網格上的修正值 $\bold E_h$；
5. 在細網格上繼續求解 $\bold A_h \bold u = \bold b_h$，仍是用迭代方法迭代3步，注意這時候迭代初值選擇爲 $\bold u_h+\bold E_h$，即包含了粗網格的修正成分。

多重網格通常會構造好多層的粗細網格，而粗細網格遍歷循環方式能夠由V、W、F型，如圖

1.6 例子

若是沒有例子，那麼仍是搞不明白，這時候能夠找個很簡單的小例子來試試看，就用前面提到的一維零邊值問題，控制方程和邊界條件以下
$-\frac{d^2\phi}{dx^2}=f(x),\quad 0\le x\le 1\\ \phi(x=0)=\phi(x=1)=0$
用前面提到的網格剖分方式，細網格等分8單元，粗網格等分4單元，插值矩陣、限制矩陣、細網格係數矩陣、粗網格係數矩陣前面都推出來了。咱們讓方程中的 $f(x)=-\sin(2\pi x)$，那麼這個問題的精確解就很容易獲得是 $\phi(x)=-\sin(2\pi x)/(2\pi)^2$，即細網格上的精確解
$\bold u=[u_1,u_2,u_3,u_4,u_5,u_6,u_7]^T= [-0.0179, -0.0253, -0.0179, -0.0000, 0.0179, 0.0253, 0.0179]^T$
我們用最簡單的V型粗細兩層多重網格算法來算算看，以Gauss-Seidel迭代爲例（能夠寫個小程序來作哈，很簡單的）。

1. 選擇初始值 $\bold u=\bold 0$，Gauss-Seidel迭代方法迭代3步，求解 $\bold A_h \bold u = \bold b_h$，獲得 $\bold u_h$爲
   $\bold u_h=[ -0.0122, -0.0172, -0.0122, -0.0000, 0.0122, 0.0172, 0.0122]^T$
2. 計算細網格上面的殘差 $\bold r_h=\bold b_h - \bold A_h \bold u_h$，得
   $\bold r_h = [ -0.2500, -0.3536, -0.2500, -0.0000, 0.2500, 0.3536, 0.2500]^T$
   將 $\bold r_h$傳遞到粗網格上面去，用 $\bold r_{2h}=\bold R_h^{2h}\bold r_h$獲得粗網格上面的殘差 $\bold r_{2h}$
   $\bold r_{2h}=[-0.3018, -0.0000, 0.3018]^T$
3. 求解粗網格上方程 $\bold A_{2h} \bold E_{2h}=\bold r_{2h}$，Gauss-Seidel迭代方法迭代3步，注意迭代初值選擇爲 $\bold E = \bold 0$，迭代完成後算出了修正值 $\bold E_{2h}$
   $\bold E_{2h}=[-0.0094, -0.0000, 0.0094]^T$
4. 將粗網格上的修正值 $\bold E_{2h}$插值回細網格，用 $\bold E_{h}=\bold I_{2h}^{h}\bold E_{2h}$，得到細網格上的修正值 $\bold E_h$
   $\bold E_h=[ -0.0047, -0.0094, -0.0047, -0.0000 , 0.0047, 0.0094, 0.0047]^T$
5. 在細網格上繼續求解 $\bold A_h \bold u = \bold b_h$，注意這時候迭代初值選擇爲 $\bold u_h+\bold E_h$
   $\bold u_h=\bold u_h+\bold E_h=[ -0.0169, -0.0267, -0.0169, -0.0000, 0.0169 , 0.0267 , 0.0169]^T$
   顯而易見，粗網格的修正起了做用，這時候的 $\bold u_h$更加接近精確解了，採用Gauss-Seidel迭代方法再迭代3步，得
   $\bold u_h=[ -0.0189， -0.0257， -0.0189， -0.0000， 0.0189， 0.0257， 0.0189]^T$
   其殘差 $\bold r_h=\bold b_h - \bold A_h \bold u_h$的2範數爲0.2165。還能夠繼續進行該細-粗-細的循環直到解收斂。
   注意，若是不用多重網格，只用Gauss-Seidel迭代方法只在細網格上進行迭代，進行了6步迭代後，所獲得的的值 $\bold u_h$爲
   $\bold u_h=[ -0.0165， -0.0233 ， -0.0165 ， -0.0000 ， 0.0165 ， 0.0233 ， 0.0165]^T$
   其殘差 $\bold r_h=\bold b_h - \bold A_h \bold u_h$的2範數爲0.2500，跟多重網格相比（0.2165）仍是有差距的，代表多重網格方法的確是能起到加速收斂的做用。

2 代數多重網格（AMG）算法

把非結構網格上的多重網格算法，跟代數多重網格算法放到一塊兒來說，二者並無區別。由於非結構網格的節點鄰近拓撲關係是未知的、無規則的，而代數多重網格則更是直接刨掉了網格的束縛，然而二者的係數矩陣是相同的，即都是稀疏矩陣，哪裏出現非零元素是無規則的，那麼它倆放一塊講就是了。

注意，係數矩陣 $\bold A$仍然要求是對稱正定矩陣，且對角元素是正值，而非對角元素是非正值（即，要麼是0值，要麼是負值）。雖然這並非代數多重網格的必要條件，然而AMG最初就是針對這種類型的矩陣提出來的，若是矩陣跟這種類型的矩陣誤差較遠的話，那麼極可能AMG並無什麼效果。

跟前面的幾何多重網格算法同樣，須要作下面的事情（注：雖然粗細網格並不像結構網格那樣單元長度爲h和2h，可是這裏仍然用h和2h分別表示細和粗網格層次）

• 選取粗網格節點，通常是在原細網格節點的基礎上，利用矩陣 $\bold A$的係數，根據必定的條件來選取一些網格點做爲粗網格；
• 插值矩陣 $\bold I_{2h}^h$，用來將粗網格上的修正值傳遞到細網格上；
• 限制矩陣 $\bold R_h^{2h}$，用來將細網格上的殘差傳遞到粗網格上；
• 粗網格上的係數矩陣 $\bold A_{2h}$。

2.1 插值矩陣 $\bold I_{2h}^h$

先無論粗網格是如何獲得的，假設已經整好了粗網格，那麼怎麼把粗網格上的值傳遞到細網格上呢？即，怎麼由粗網格節點值插值獲得細網格節點值呢？因爲此時節點關係比較複雜，因此並不能一會兒就把插值矩陣給整出來，而是須要逐個節點去處理纔好的。

假設已經找好了粗網格節點（即從本來的細網格中篩選好了粗網格節點），即，把最初的細網格節點 $\{1,2,3,...,n\}$分爲了 $C \cup F$，粗網格節點集合 $C$（從本來的細網格節點中篩選出來的做爲粗網格的那些節點們，它們也與原來的某些細網格節點重合，也屬於細網格節點），和僅僅是細網格的節點集合 $F$（即從細網格中篩掉的沒資格作粗網格的那些節點們，大家這些渣渣就不要摻合粗網格的計算了）。那麼插值要解決的問題就是，知道了在粗網格節點 $C$上的修正值 $e_i$，如何將它們傳遞到僅僅是細網格的那些節點 $F$上去呢？

若是細網格上的節點 $i$是和粗網格節點重合的（ $i \in C$），即以前選出來作粗網格的那些細網格節點，那麼直接用 $e_i$就行了，即 $e_i=e_i$。

若是細網格上的節點 $i$是隻屬於細網格的（ $i \in F$），即以前選粗網格節點時被篩掉的那些細網格節點，那就須要插值來計算了，即 $e_i=\displaystyle \sum_{j \in C_i}\omega_{ij}e_j$，注意插值的時候只能用到粗網格節點的值，並且只能用到節點 $i$的鄰近節點的值，還得是那些對節點 $i$的影響較大的節點的值（人輕言微，影響太小的話我們直接就把它們咔咔了，另行處理），一句話，只能用到那些對節點 $i$影響較大的鄰近粗網格節點的值，也就是加和範圍中的集合 $C_i$。

也就是說，插值矩陣所提供的插值算法應使得每一個細網格上的節點 $i$的修正值爲
$(\bold I_{2h}^h \bold e)_i = \left\{ \begin{array}{rcl} e_i & & if~i \in C\\ \displaystyle \sum_{j \in C_i}\omega_{ij}e_j&& if ~ i \in F \end{array} \right.$
問題來了，這個 $j \in C_i$究竟是哪些節點呢？ $\omega_{ij}$又該如何計算呢？

對於某個僅僅是細網格的節點 $i \in F$，能夠把細網格上關於它的係數矩陣 $\bold A$中的那一行挑出來，並用其來反映各節點修正值的關係，即
$a_{ii} e_i \approx - \sum_{j \in N_i} a_{ij} e_j$
其中 $N_i$表示在細網格中節點 $i$的鄰近節點，那麼這些節點要麼屬於 $C$（粗網格節點），要麼屬於 $F$（僅爲細網格節點），再考慮到節點之間影響程度的強弱，能夠 $N_i$分紅三類來處理：

• 對節點 $i$影響較強的鄰近粗網格節點，即 $C_i$，用於作節點 $i$插值的粗網格節點集合；
• 對節點 $i$影響較強的鄰近細網格節點，用 $D_i^s$表示（s即strong）；
• 對節點 $i$影響較弱的鄰近網格節點，用 $D_i^w$表示（w即weak），這些節點既有多是來自粗網格 $C$的，也有多是來自僅細網格 $F$的；

那啥叫影響強，啥叫影響弱呢？能夠簡單地用係數矩陣 $\bold A$的係數來定義兩個節點的影響程度：

---

給定一個閾值因子 $0\lt\theta\le1$，若
$-a_{ij} \ge \theta \displaystyle \max_{k\ne i}\{-a_{ik}\}$
則說 $j$節點的值 $u_j$對 $i$節點的值 $u_i$的影響較大，或者反過來， $i$節點的值 $u_i$受（被） $j$節點的值 $u_j$的影響較大， $u_i$對 $u_j$的依賴較大。

---

用更通俗易懂的話來說一下，係數 $-a_{ij}$（注意 $j \ne i$）（前面說過，係數矩陣的非對角元素是0或者負值，因此加了負號就是正值了）反映的是節點 $j$對節點 $i$的影響程度的大小，那麼從第 $i$行全部的係數中（除掉對角元 $a_{ii}$）選擇最大的那個 $\displaystyle\max_{k\ne i}\{-a_{ik}\}$，它反映的就是對節點 $i$影響最強的那個節點 $j$的影響程度，再乘上定義好的係數 $\theta$，做爲臨界值 $\theta \displaystyle \max_{k\ne i}\{-a_{ik}\}$，再把那些影響程度大過此臨界值的節點，即 $-a_{ij} \ge \theta \displaystyle \max_{k\ne i}\{-a_{ik}\}$的節點，選擇爲對 $i$節點影響較大的節點，就理所固然了。顯然，這個影響較強的節點跟粗網格節點沒有必然聯繫，即，這些節點有的可能會被選成粗網格節點，而有的則不會被選成粗網格節點。

繼續回來看 $a_{ii} e_i \approx - \displaystyle\sum_{j \in N_i} a_{ij} e_j$，因爲把 $N_i$分紅了三類 $C_i,D_i^s,D_i^w$，因此能夠把它寫成
$\begin{aligned} a_{ii} e_i & \approx - \sum_{j \in N_i} a_{ij} e_j \\ & = - \sum_{j \in C_i} a_{ij} e_j- \sum_{j \in D_i^s} a_{ij} e_j- \sum_{j \in D_i^w} a_{ij} e_j \end{aligned}$
跟插值公式 $e_i=\displaystyle \sum_{j \in C_i}\omega_{ij}e_j$相對照，不難發現，右端第1項就不用管了，直接把 $-a_{ij}/a_{ii}$給到 $\omega_{ij}$就行了。而右端的第2和第3項是須要處理的，先看第3項，因爲它們對於 $i$節點的影響較弱，因此直接把它們加到對角元素上就好，即把 $e_j$替換成 $e_i$，或者說，把它們挪到左端去，加和到對角元上，以下：
$\left(a_{ii} + \sum_{j \in D_i^w} a_{ij}\right) e_i \approx - \sum_{j \in C_i} a_{ij} e_j- \sum_{j \in D_i^s} a_{ij} e_j$
可能有人會好奇，爲啥不直接把它們給扔掉呢？由於要保持權係數加和是1的準則啊，本來的 $a_{ii}=\displaystyle \sum_{j \in N_i} a_{ij}$，即 $1=\displaystyle \sum_{j \in N_i} a_{ij}/a_{ii}$，若是直接把這些影響小的項拋掉，就不知足這個原則了，即獲得的加權係數 $\omega_{ij}$加起來也不是 $1$了。

接下來處理第2項，因爲對 $i$節點影響較大，這些節點的值必須考慮，可是這些節點的值又不知道，因此稍微麻煩些……把這些 $D_i^s$中節點 $j$的值 $e_j$想辦法用 $C_i$中的節點值來表示，從 $j$節點的離散方程中尋找答案，它們一樣能夠寫成和 $i$節點相同的形式
$\begin{aligned} a_{jj} e_j & \approx - \sum_{k \in N_j} a_{jk} e_k \end{aligned}$
儘管 $j$節點有不少鄰近節點 $N_j$，但咱們僅考慮哪些既屬於 $C_i$又屬於 $N_j$的那些節點，即，這些節點一方面是節點 $j$的鄰近節點，另外一方面它們又是對節點 $i$影響較大的鄰近粗網格節點。用這些節點值的加權插值獲取 $j$節點的值，加權係數用 $a_{jk}$來獲取，即
$e_j \approx \frac{\displaystyle\sum_{k\in C_i}a_{jk}e_k}{\displaystyle\sum_{k\in C_i}a_{jk}}$
這樣子，就能夠獲得粗網格向細網格插值的式子了。
$\left(a_{ii} + \sum_{j \in D_i^w} a_{ij}\right) e_i \approx - \sum_{j \in C_i} a_{ij} e_j- \sum_{j \in D_i^s} a_{ij} \left( \frac{\displaystyle\sum_{k\in C_i}a_{jk}e_k}{\displaystyle\sum_{k\in C_i}a_{jk}} \right)$
與插值式子 $e_i=\displaystyle \sum_{j \in C_i}\omega_{ij}e_j$比照便可獲得插值係數 $\omega_{ij}$的值。

光說不練假把式，須要一個小例子來搞搞明白，假設二維問題，結構網格（只是爲了說明問題），均分爲 $n\times n$網格，每一個內部節點的離散框架爲
$\begin{bmatrix} -\displaystyle\frac{1}{2} & -2 & -\displaystyle\frac{1}{2} \\ -1 & \displaystyle\frac{29}{4} & -1 \\ -\displaystyle\frac{1}{8} & -2 & -\displaystyle\frac{1}{8} \end{bmatrix}$
對於每一個內部節點 $i$，我們用到了它的東、西、南、北、東北、西北、東南、西南共8個鄰近網格點來構造離散框架。假設粗糙網格直接選擇 $i$節點的東、西、南、北四個節點來構造（ $i$節點不是粗網格點），選擇閾值係數 $\theta=0.2$來定義影響較強的節點，則根據 $-a_{ij} \ge \theta \displaystyle \max_{k\ne i}\{-a_{ik}\}=0.2*2=0.4$，可知東、西、南、北這4個節點屬於對 $i$節點影響較強的粗網格點 $C_i$，東北、西北節點屬於對 $i$節點影響較強的細網格節點 $D_i^s$，而東南、西南節點則屬於對 $i$節點影響較弱的細網格節點 $D_i^w$（這個例子中沒有對 $i$節點影響較弱的粗網格節點，否則也放在 $D_i^w$中）。如圖（圖中對網格節點的編碼是先 $x$再 $y$的順序）

圖中的空心圓圈標識粗網格節點，實心圓圈標識細網格節點，實線表示對 $i$節點影響較強的粗網格節點 $C_i$，虛線表示對 $i$節點影響較強的細網格節點 $D_i^s$，點線表示對 $i$節點影響較弱的細網格節點 $D_i^w$。

那麼寫出 $i$節點的離散格式
$\begin{aligned} \frac{29}{4}e_i&=\underbrace{2e_{i+n}+2e_{i-n}+e_{i+1}+e_{i-1}}_{C_i}\\ &+\underbrace{\frac{1}{2}e_{i+n-1}+\frac{1}{2}e_{i+n+1}}_{D_i^s}\\ &+\underbrace{\frac{1}{8}e_{i-n-1}+\frac{1}{8}e_{i-n+1}}_{D_i^w} \end{aligned}$
$D_i^w$中的節點爲東南 $i-n+1$、西南 $i-n-1$節點，對 $i$節點影響較弱，直接將其用 $e_i$代替，歸併到左端項中，即
$\begin{aligned} \left( \frac{29}{4}-\frac{1}{8}-\frac{1}{8} \right) e_i&=\underbrace{2e_{i+n}+2e_{i-n}+e_{i+1}+e_{i-1}}_{C_i}\\ &+\underbrace{\frac{1}{2}e_{i+n-1}+\frac{1}{2}e_{i+n+1}}_{D_i^s} \end{aligned}$
$D_i^s$中的節點爲東北 $i+n+1$、西北 $i+n-1$節點，它們對 $i$節點影響較強。東北 $i+n+1$節點的鄰近節點中屬於 $i$節點的 $C_i$集合的（i節點的東、西、南、北節點），只有i節點的北 $i+n$、東 $i+1$兩節點（分別是 $i+n+1$節點的西、南節點），它們在東北節點離散方程中對應的係數分別爲 $-1$和 $-2$；而西北 $i+n-1$節點的鄰近節點中屬於 $C_i$的，只有 $i$節點的北 $i+n$、西 $i-1$節點（ $i+n-1$節點的東、南節點），它們在西北節點離散方程中對應的係數分別爲 $-1$和 $-2$，故可將 $D_i^s$中的節點，即東北 $i+n+1$、西北 $i+n-1$節點值表示爲：
$e_{i+n+1} = \frac{-1(e_{i+n+1})_W-2(e_{i+n+1})_S}{-1-2}=\frac{-1e_{i+n}-2e_{i+1}}{-3}$
$e_{i+n-1} = \frac{-1(e_{i+n-1})_E-2(e_{i+n-1})_S}{-1-2}=\frac{-1e_{i+n}-2e_{i-1}}{-3}$
將其代入上式，可得
$\begin{aligned} \left( \frac{29}{4}-\frac{1}{8}-\frac{1}{8} \right) e_i&=\underbrace{2e_{i+n}+2e_{i-n}+e_{i+1}+e_{i-1}}_{C_i}+\underbrace{\frac{1}{2}e_{i+n-1}+\frac{1}{2}e_{i+n+1}}_{D_i^s} \Rightarrow \\ 7 e_i&=\underbrace{2e_{i+n}+2e_{i-n}+e_{i+1}+e_{i-1}}_{C_i}+\underbrace{\frac{1}{2}\frac{-1e_{i+n}-2e_{i-1}}{-3}+\frac{1}{2}\frac{-1e_{i+n}-2e_{i+1}}{-3}}_{C_i} \Rightarrow \\ 7 e_i&=\frac{7}{3}e_{i+n}+2e_{i-n}+\frac{4}{3}e_{i+1}+\frac{4}{3}e_{i-1} \end{aligned}$
最終獲得了插值公式
$e_i =\frac{7}{21}e_{i+n}+\frac{6}{21}e_{i-n}+\frac{4}{21}e_{i+1}+\frac{4}{21}e_{i-1}$
可見，其知足係數加和爲1的準則，且僅用到了那些對 $i$節點影響較大的粗網格節點 $C_i$的值，完美哈。

這時候，有人可能會想到了，對於那些 $D_i^s$中的節點（對 $i$節點影響較大的細網格節點），若是在它們的鄰近網格節點中找不到屬於 $C_i$（對 $i$節點影響較大的粗網格節點）的，那上面這個方法豈不是要玩完了麼？沒錯，就是這樣的，可是不會出現這個問題，由於在篩選粗網格節點的時候根據必定的準則來篩選，讓選出來的粗網格節點不會出現上面這種尷尬的無從處理的狀況不就行了麼，下面我們就看看怎麼去篩選粗網格節點。

2.2 粗網格構造（篩選）

粗網格節點是從本來的細網格節點中篩選出來的，那麼哪些節點纔有資格選中做爲粗節點呢？先把那些對細網格中 $i$節點影響較強的節點集合定義爲 $S_i$（那些個對 $i$節點愛的死去活來的粉絲們），再把 $i$節點影響強烈的節點集合定義爲 $S_i^T$（那些個 $i$節點愛的死去活來的明星們），注意這倆集合可未必同樣啊，就像蘿蔔愛青菜、青菜愛大豆同樣，未必是相互的。那麼在粗網格節點的選擇中，要知足下述兩個準則H1和H2

• H1 對於細網格中的每一個節點 $i$，對 $i$節點影響較強的節點 $j \in S_i$，要麼選中做爲對 $i$節點影響較強的粗網格節點集合 $C_i$（如前所述，這個集合將用於對 $i$節點的插值計算），要麼（這時候它屬於對 $i$節點影響較強的僅細網格節點 $F_i$中的 $D_i^s$）其至少被一個 $C_i$中的節點強烈影響着（這樣子的話，就不會出現上節最後提到的那種尷尬的狀況了）。
• H2 粗網格節點集 $C$應該是原網格集合中知足這個特性的最大子集合，即， $C$中任何一個節點都不會強烈依賴於 $C$中的其它節點。

H2仍是蠻好理解的，它意思是說，選擇出來的粗網格節點呢，應該是這樣子的，每個節點都不受其它節點的強烈影響，反過來講，任何一個粗網格節點對其餘的粗網格節點並不會產生強烈的影響。爲啥要知足這個呢？由於呢，粗網格節點之間不影響的話（影響較弱），粗網格節點就只好影響細網格節點了，這不正是咱們想要的麼？正好能夠把粗網格上的偏差修正很好地傳遞到細網格節點上去呢。

H1相對比較難理解一些，它的意思是說，在細網格中，某個節點 $j$若是對節點 $i$的影響較強，那麼（1）這個節點 $j$要麼被選擇成爲粗網格節點，若是這樣子的話，它就成了對 $i$節點影響較強的粗網格節點集合 $C_i$中的一員，在計算中作粗網格向細網格插值的時候，主要就是用到這個粗網格點上的數據了；（2）若是不幸地，這個網格節點 $j$沒被選擇成爲粗網格節點，這時候它其實是僅細網格節點（屬於集合 $F$），且是對 $i$節點影響較強的細網格節點（屬於集合 $D_i^s$），那麼它至少得有一個能夠依賴的屬於 $C_i$的節點才行，即，至少得有一個對 $i$節點影響較強的粗網格節點，其同時對這個 $j$節點的影響也較強，如此一來，才能在上節講的插值過程當中，把這個 $j\in D_i^s$節點值表示成 $k\in C_i$的形式啊，纔不會出現上節課末尾所講的那種尷尬的狀況呢。

不難想象，要同時知足H1和H2是蠻難的一件事情。鑑於插值過程的思路是須要知足H1的，因此H1是不管如何都必需要嚴格知足的，否則沒辦法搞插值了，而H2則是個指導性的準則，未必須要嚴格知足，由於即使是多上那麼幾個節點，也不會對插值形成什麼影響，只不過是增長了一些計算量罷了。

粗網格的篩選也並不是一蹴而就的，通常要分兩步來走。第一步是對網格點作個初篩（着色，coloring）的工做，即大致上把最初的細網格劃分爲粗網格 $C$和僅細網格 $F$，第二步則更爲細化，將交換某些 $F$和 $C$中的節點，以便讓H1準則嚴格知足。

2.2.1 第1步 初篩

一開始，給每一個細網格節點 $i$賦個值，用來衡量其成爲粗網格節點 $C$的可能性的大小，一個最簡單的辦法是計算那些 $i$節點影響較強的節點的數目 $\lambda_i$，這些個 $i$節點影響較強的節點集合是 $S_i^T$（前面定義過了）。舉個栗子，若是某節點對於其它7個節點影響較強，那麼其成爲粗網格的可能行就是7，而若是某節點對其餘2個節點影響較強，其成爲粗網格的可能性就是2。顯然這個值越大，其越有可能成爲粗網格節點集合 $C$中的一員。一旦肯定好了每一個節點的可能性，選擇具備最大可能性 $\lambda_i$的節點做爲第1個粗網格節點 $C$。

接下來，咱們選好的這個粗節點 $i$，其強烈影響着周圍的節點，那麼這些被其強烈影響的節點就不該該再被選入粗節點了（由於插值的時候它們能夠用粗網格節點值來獲取了），它們應該被放入僅細節點 $F$中，即把 $S_i^T$中的節點所有放入 $F$中。

隨後，再看看其它節點中，有沒有那些對這些新的僅細節點 $F$影響較強的節點，它們有可能做爲新的粗網格節點呢（即 $F$中的節點一方面受 $i$節點的強烈影響，一方面又受到這些新節點的強烈影響，因此用它們來給 $F$中的節點作插值真是再好不過了）。因此，對於每一個新的僅細網格節點 $F$（同時也是屬於 $S_i^T$， $i$節點影響較強的那些節點），即 $j \in S_i^T$節點，若是有其它節點 $k$強烈影響着 $j$節點的話，即 $k \in S_j$，把 $k$節點的權值 $\lambda_k$加1，表示其成爲粗網格的但願又更進了一步。簡單來講，對於每一個 $j \in S_i^T$，將 $k \in S_j$的權值 $\lambda_k$增長1。

緊接着，繼續選第2個節點了，固然是選擇具備最大可能性 $\lambda$的節點了，做爲新的粗網格點 $i$放入 $C$中；以後，再把這個點強烈影響的節點集合 $S_i^T$