ピボット選択

ガウスの消去法における前進消去の式をもう一度みてみる．

右辺の第2項にで割るという計算が含まれている．元々対角要素に0を含んでいる，もしくは，前進消去の過程で対角要素に0が出てくるとその時点で計算が止まってしまう．また，0でなくても絶対値が非常に小さい数だと桁落ちが起こり，誤差が増大する原因ともなる．

この問題に対処するために，ピボット(pivot:枢軸)選択(pivotal elimination method)，もしくは，ピボッティング(pivoting) と呼ばれる方法を説明する．まず，連立方程式において式の順番は交換可能であることに注目する．つまり，行目の処理を行う際，その対角要素と同じ列でまだ処理していない要素の値を調べ，その絶対値が最大の行と行を入れ替えた後で前進消去の式を適用する．これがピボット選択である．ちなみに，このように行の入れ替えだけ行うことを部分的ピボッティング(partial pivoting)と呼び，さらに列方向についても調べて入れ替えを適用することを完全ピボッティング(complete pivoting)と呼ぶ．ただし，部分的ピボッティングでは未知数の順番は変わらないが，完全ピボッティングの場合は未知数の順番も変化するので注意．

以下にピボット選択を含めたガウス消去法のコード例を示す．


inline void Pivoting(vector< vector<double> > &A, int n, int k)
{
        int p = k;     double am = fabs(A[k][k]);    for(int i = k+1; i < n; ++i){
        if(fabs(A[i][k]) > am){
            p = i;
            am = fabs(A[i][k]);
        }
    }
        if(k != p) swap(A[k], A[p]);
}
 

int GaussEliminationWithPivoting(vector< vector<double> > &A, int n)
{
            for(int k = 0; k < n-1; ++k){
                Pivoting(A, n, k);
 
        double akk = A[k][k];
        for(int i = k+1; i < n; ++i){
            double aik = A[i][k];
            for(int j = k; j < n+1; ++j){
                A[i][j] = A[i][j]-aik*(A[k][j]/akk);
            }
        }
    }
 
                A[n-1][n] = A[n-1][n]/A[n-1][n-1];
    for(int i = n-2; i >= 0; --i){
        double ax = 0.0;
        for(int j = i+1; j < n; ++j){
            ax += A[i][j]*A[j][n];
        }
        A[i][n] = (A[i][n]-ax)/A[i][i];
    }
 
    return 1;
}

Menu