迟到了迟到了 ≡_(ゝヽε:)ﾉ
主要参考stable_fluid自己整了个 Eulerian Fluid 的实现 。放出来见笑了 ：Github仓库 。
主要实现：
- Advection 直接就是三阶的RK
- Projection 分别测试了：
- Jacobian
- conjugate gradient
- MGPCG
以前没接触过物理模拟，我也不太确定实现的效果是否准确。
视觉效果上大家发现有什么不对劲还请不吝指正一下_(┐「ε:)_。
这里直接放出来gif对比一下不同projection method的效果 ：
-
Jacobian iteration , steps = 30 , fps = 60.0
emmmm … 比起烟 ，这个更有雾的感觉吧，Jacobian这个帧数拉满了，动画还挺流畅的。 -
Conjugate gradient , max_iter = 20 , fps = 7.1
上了共轭梯度法之后，可以冒出来很多涡了，只需要20次迭代 ， 但帧率低了很多。 -
MGPCG , max_iter = 20 . fps = 0.66
MGPCG搓出来的涡相对的细腻一些 。但是残差降到到个位数时基本就要10个step了。帧率没法直视，这个只能录下来看了。 -
Jacobian iteration , steps = 200 , fps = 48
加大Jacobian 的迭代次数，解出来的效果开始接近共轭梯度法了。
一步步实现下来还是积累了不少问题的：
-
虽然MGPCG的收敛极快，但是Multigrid这个preconditioner的计算量真是太大了。同样是最大迭代20次的情况下，效果上没有比普通的共轭梯度法好多少。从帧率上比较，MGPCG的性能几乎是CG的十分之一 ，而CG的方法又几乎是Jacobian的十分之一…
-
对比example里面mgpcg实现，发现我自己的实现跑起来还要慢上几倍 。详细对比了一下，发现其实在实现流体模拟时，CG方法里面设置的各种中间矩阵，并不能保持稀疏，一个step之后全部都变稠密了…
-
另外跑了mgpcg_advanced的profile , 不带multigrid的时候 ，发现性能瓶颈都集中在各个reduce操作里面。taichi对于这类操作，实际上是编译成了atomic_add ？应该会有性能更好的规约方案（不确定cuda的这个有没有帮助，只要是无锁并行应该都能得到较大提升）。
-
为了一些性能上的方便，example里的MGPCG实现还是跟具体求解的矩阵有很多耦合的地方。怎样单独分离出来一个比较高效的PCG的线性方程组求解器，这个还得再想想。
题外
还有我发现example都比较喜欢直接写出并行for循环来计算矩阵。
事实上我在用structural for来计算稀疏矩阵间，或稀疏和稠密矩阵之间的运算时，尽管两个矩阵的shape匹配，但由于稀疏的存储会导致遍历时跳过未激活的单元，做乘法以外的运算时一不小心就会漏掉一部分，经常算错，所以只好稠密地遍历了，但是也没能充分利用稀疏矩阵的优势。
比如两个稀疏量element-wise相加就必须稠密遍历，但如果是element-wise相乘的话用任意一个ti.grouped(a)或者ti.grouped(b)结果都是正确的。
如果想要给所有类型的运算（乘除、加减等等）都提供尽量稀疏的遍历，则需要针对两边的激活单元进行组合才能保证遍历准确。如果taichi提供了带重载的运算符给ti.var之间运算，也就可以省略很多暴力写for的代码了。