工作规划(1) 基于LU分解的千万阶稀疏矩阵求解器

以下是投给某个研究院的一份工作计划。在稀疏求解器上花了六、七年功夫，也是对这些年的一个小结吧。

目前，千万阶矩阵求解以迭代法（含多重网格法）为主。迭代法算法简单，易于实现。处理对称正定矩阵及某些特殊类型矩阵有较好的效果。但总体上很难做到通用，稳定。精度较低，效率上先天不足。直接法（基于LU分解）的精度高，通用性好。尽管直接法同样有很多困难，但这些困难是可以克服的；其潜力巨大，有望得到广泛的应用。

1 迭代法的不足
尤其对于非对称矩阵，以GMRES，BICG为代表的迭代法需要高质量的预处理器，而目前还没有通用而成熟的预处理算法，在理论上更是难以分析。以ILU为代表的预处理算法越来越复杂，代码长度，运行时间甚至超过迭代法本身，效果仍依赖于具体问题。
2直接法的优势
1）稳定，通用，精度高
基于恰当的主元技术及预处理技术，目前直接法的表现相当稳定，精度也比迭代法高很多。
2）CPU效率可高一个数量级
目前单颗CPU的峰值浮点性能已有超过100G，即每秒1000亿次运算。但只有为数不多的代码能实现这个峰值，其中就有BLAS-3。通过算法优化，稀疏LU分解的绝大部分运算可直接调用BLAS-3, 而迭代法一般只能调用BLAS-2。两者的CPU效率相差一个数量级。
3）更适合于反复求解的场合
稀疏矩阵通常需要反复求解，即矩阵不变，而右端项不断变化。例如非线性方程组求解，拟牛顿迭代以，有限元中的多种荷载等。直接法只需数值分解一次，然后针对不同的右端项回代。而回代的速度很快，一般远小于迭代法的时间。因此，当反复求解的次数较多，直接法所需的时间会更少。

3 直接法的难点及分析
目前，直接法难以推广的因素很多，主要如下：
1） 算法复杂，要求很高的编程水平。目前主要的稀疏LU分解算法有multifrontal, frontal, supernodal等，都需要一定的图论与组合数学基础，主流求解器的核心代码都在万行以上。这些求解器都采用很多参数，针对具体问题如何取值也有一些难度。目前成熟的直接法求解器较少，有一些已商业化。总体上还有很多提高的空间。
2）需要大量内存
这个是限制直接法的主要瓶颈所在。就百万阶矩阵实测来看，需要10G左右内存，千万阶所需内存估计在100G以上。估计PC，工作站的内存要数年后才能达到这个数量级。因此需要外存模式，这将增加额外的时间，并使算法更加复杂。
3) 计算量巨大
千万阶稀疏LU分解的计算量可认为是k*1.014。看上去很大，但仍是可以克服的。目前已有CPU浮点性能达到1.011次/秒，考虑并行效率，稀疏分解本身的效率，在多CPU上可达到1.012~1.013次/秒。这样总的求解时间预期为数分钟到数个小时，还是可以接受的。

4 GSS的成功实践
GSS已可在PC上求解百万阶矩阵，基于这些算法和模型，开发千万阶求解器是可行的
1）GSS运用多种矩阵排序算法，包括一些独创的算法，可有效降低内存需求。
2）GSS的并行模块采用多线程，应可推广到多机，多CPU模型。
3）GSS采用先进的符号分析模型，可实现新的并行回代算法，要优于现有的算法。
4）千万阶矩阵求解需要采用外存模式，而GSS的out-core算法稳定，高效。

5具体工作目标
1）矩阵排序算法研究，将千万阶矩阵的内存控制在100G甚至更少。
2）基于多CPU（含大型机，网格，工作站等）的千万阶并行求解模块。
3）基于多CPU（含大型机，网格，工作站等）的千万阶并行回代模块。
4）争取推出可在现有工作站上运行的千万阶求解器。