油藏数值模拟中大型稀疏矩阵算法综述.docx

大型稀疏矩阵数值解法综述
前言
线性方程组的求解方法基本上可以分为两种：直接法和迭代法。直接法是以Gauss消去法为基础的方法，实际上这种方法是初等代数中解多元一次方程组的直接推广。其基础思路是：对原方程组经过一定的运算处理后，逐个消去部分变量，最后得到一个与原方程等价的、便于逐步求解的方程组，以解除各个变量的值。如果不考虑计算式可能产生的舍入误差，则可以认为直接法是一种精确的方法，它可以一次性求得原线性方程组的解。与直接法相比，迭代法是一种近似的解法，它的基本思想是：先估计一组变量的数值，作为原方程组的第一次近似，也称作迭代初值，然后利用某种迭代处理，逐次修改这组数值，得到解的第二次、第三次以至第k次近似值。在方程组满足一定条件的前提下，这些结果可以逐次逼近原方程组的真实解。进行有限次重复计算后，如果所得到的近似解满足所规定的误差范围，则认为最后一次迭代的近似值就是原方程的解。
随着石油行业的发展，对油藏的认识也不断地深入，通过应用大型计算机，油藏数值模拟的规模也越来越大，另外随着一些新的方法如：、全隐式、局部加密网格方法的使用，以及混合驱、热力驱、化学驱等油藏模型的出现，使方程组矩阵的性质更加复杂多变，甚至带有严重的病态性质。开发能够快速求解各种复杂的矩阵线性方程组的新方法就成为油藏数值模拟进一步发展的一个重要方向。
通过对IMPES、全隐式的油藏数值模拟常用方法的研究发现，油藏数值模拟中所研究的矩阵只是在少数元素上为非零值，大部分元素值为零值，数学上称之为稀疏矩阵。当研究对象为二维油藏数值模拟时，稀疏方程组的系数矩阵为对称格式，当研究对象为三维油藏数值模拟时，稀疏方程组的系数矩阵则变为非对称格式，直接法在求解这种稀疏矩阵时，会使原矩阵新添加大量的非零元素，增加计算量和存储量，对于求解大规模稀疏矩阵，速度会变得非常缓慢，甚至无法求解，因此直接法已经逐渐被迭代法所取代，不再作为求解大型油藏数值模拟线性方程组的方法。
迭代法包括古典迭代法和Krylov子空间迭代法。古典迭代法包括Jacobi.、GaussSeidel、SOR、ssor等方法，因其收敛速度很慢，目前也已很少用于求解大型稀疏线性方程组，而是用于与一些效率更高的方法相结合使用。目前，Krylov子空间方法是求解大型稀疏矩阵线性方程组最流行和最有效的方法之一，也是当前研究的热点，其主要思想是为各迭代步递归地构造残差向量，即第n步的残差向量通过系数矩阵A的某个多项式与第一个残差向量相乘得到：

通常，迭代多项式的选取应使用所构造的残差向量在某种内积意义下相互正交，从而保证某种极小性（极小残差性），达到快速收敛的目的。Krylov子空间方法具有存储量少，计算量少且易于并行等优点，非常适合于求解大型稀疏线性方程组，结合预条件的Krylov子空间迭代法是目前求解大型稀疏线性方程组的最主要方法，
1Krylov子空间概念
给定初值，求解稀疏线性方程组：

其中
设为m维子空间，一般投影方法是从m维仿射子空间中寻找近
似解使残差向量满足Petrov—Galerkin条件

其中为另一个m维子空间，为迭代初值。如果是Krylov子空间，
则上述投影方法就称作Krylov子空间方法。Krylov子空间定义为

理想的Krylov子空间方法求解稀疏线性方程组的过程具有以下特征：
极小剩余性或极小误差性（保证收敛速度快）；
每一步迭代的计算量少，存储量小，以保证计算高效性。
2Krylov子空间方法分类
选择不同的和就得到了不同的Krylov子空间算法。
2.1基于正交投影方法
取，这一类Krylov子空间方法称为正交投影方法。共轭梯度（CG）方法是其中最重要的方法，此时要求矩阵A是对称正定矩阵。CG该方法要求在子空间中的能量范数到达极小，且具有短的迭代计算公式。此后发展得到的与各种预条件相结合的方法成为求解大型对称正定稀疏线性方程组最主要的方法。
2.1.1共轭梯度法（CG）
具体流程为：
1.任选初值，计算初值残差和内积（），并置,置j=0。
2.计算参数更新向量与残向量，若满足精度要求，则停机。
3.计算，置j=j+1,转（2）
精度测试
维数为9*9矩阵：
通过直接法所求的的矩阵真实解和通过CG法迭代求的的近似解误差值的数量级在10-12。
维数为36*36矩阵：
通过直接法所求的的矩阵真实解和通过CG法迭代求的的近似解误差值的数量级在10-8—10-12。
维数为64*64矩阵：
通过直接法所求的的矩阵真实解和通过CG法迭代求的的近似解误差值的数量级在10-8—10-9。
维数为169*169矩阵：
通过直接法所求的的矩阵真实解和通过CG法迭代求的的近似解误差值的数量级在10-9—10-10。
2.1.2预条件共轭梯度法
共轭梯度法在解决大型线性方程组本身是有不足之处的，当系数