参考文献：A Fast Unsmoothed Aggregation Algebraic Multigrid Framework for the Large-Scale Simulation of Incompressible Flow
考虑求解方程-\nabla^2 p=b，或写成矩阵形式Ap=b，加负号是为了让矩阵正定。向量p的每一项对应某个n_x\times n_y\times n_z里面的一个格子。我们采用转换积分的方式：
\int_V\nabla^2 p\mathrm{d}V=\oint_S\nabla p\cdot\bm{n}\mathrm{d}S.
这里面V是一个格子，S是格子的六个面，实际就是在六个面上把p的相应偏导乘以面积加起来。假设格子的长度是h，那么相应的系数矩阵就是（忽略各种边界条件）：
A_{\alpha,\beta}=\begin{cases}
6h,\text{ }\alpha=\beta,\\
-h,\text{ }\alpha\neq \beta\text{ but adjacent,}\\
0,\text{ otherwise.}
\end{cases}
这个h是这么来的：想象p_\alpha和p_\beta在\alpha对应方程里面的形式，它会产生一项
\frac{-p_\beta+p_\alpha}{h}\cdot h^2=h(-p_\beta+p_\alpha).
系数就是h。然后，我们用R和P矩阵把它变换到粗网格上，格子的长度为2h。这里面取
R_{\alpha,\gamma}=\begin{cases}
1/8,\text{ }\text{if coarse voxel }\alpha\text{ covers fine voxel }\gamma,\\
0,\text{ }\text{otherwise.}
\end{cases}
我们用A^h代表细网格，A^H代表粗网格，那么A^H=RA^hP，定义P=8R^T。计算这个矩阵乘法（以下使用爱因斯坦求和约定）：
A^H_{\alpha,\beta}=R_{\alpha,\gamma}A^h_{\gamma,\eta}P_{\eta,\beta}=8R_{\alpha,\gamma}A^h_{\gamma,\eta}R_{\beta,\eta}.
只有当\alpha覆盖\gamma，\beta覆盖\eta时两个R项非零。也就是说，在\alpha和\beta的八个儿子里面各选一格，计算相应的乘积。假设\alpha,\beta相邻，那么只有边界上四对儿子的A^h非零，合起来就是4/8\cdot(-h)\cdot 8/8=-(2h)/4。如果\alpha=\beta，那么有8种\gamma=\eta的情况和24种\gamma和\eta相邻的情况，合起来是8/8\cdot(6h)\cdot 8/8+24/8\cdot(-h)\cdot 8/8=3h=6\cdot(2h)/4，是\alpha\neq \beta的六倍，符合泊松方程的性质。
随后考虑restriction和prolongation操作。Restriction应当为：
r^H=Rr^h.
即
r^H_{\alpha}=R_{\alpha,\gamma}r^h_{\gamma}.
这实际上就是对\alpha的八个儿子取平均。求解方程A^Hp^H=r^H得到p^H，然后做prolongation：
p^h_{\text{new}}=p^h_{\text{old}}+\lambda P p^H.
这里面\lambda是松弛系数，取1或者2。计算
(P p^H)_{\gamma}=P_{\gamma,\alpha}p^H_{\alpha}=8R_{\alpha,\gamma}p^H_\alpha.
这就是直接取\gamma父亲的值。但是有一个问题是，现在我们得到的A^H=-\nabla^2/4，即粗网格上求解的是
-\nabla^2p/4=r^H.
而我们希望在粗网格上求解-\nabla^2。这时候可以把右端r^H乘以4，求解
-\nabla^2p = 4r.
这个额外的scaling是采用积分形式导致的。假如求解的是逐点的微分，那么A^h里面的相邻项是-1/h，而A^H里面的相邻项是-1/2h，正好符合-\nabla^2的形式。