共役勾配法

**共役勾配法**（きょうやくこうばいほう、、CG法とも呼ばれる）は対称正定値行列を係数とする連立一次方程式を解くためのアルゴリズムである。反復法として利用され、コレスキー分解のような直接法では大きすぎて取り扱えない、大規模な疎行列を解くために利用される。そのような問題は偏微分方程式などを数値的に解く際に常に現れる。

共役勾配法は、エネルギー最小化などの最適化問題を解くために用いることもできる。

は、共役勾配法の非対称問題への拡張である。

また、非線形問題を解くために、さまざまな非線形共役勾配法が提案されている。

## 詳説
対称正定値行列**A**を係数とする*n*元連立一次方程式

の解を**x***とする。

### 直接法としての共役勾配法
非零ベクトル**u**、**v**が
 \mathbf{A} \mathbf{v} = \mathbf{0}}}
を満たすとき、**u**、**v**は**A**に関して共役であるという。**A**は対称正定値なので、左辺から内積
 \mathbf{u}, \mathbf{v}\rangle = \langle \mathbf{A} \mathbf{u}, \mathbf{v}\rangle = \langle \mathbf{u}, \mathbf{A}\mathbf{v} \rangle = \mathbf{u}^{\mathrm{T}} \mathbf{A} \mathbf{v}}}
を定義することができる。この内積に関して2つのベクトルが直交するなら、それらのベクトルは互いに共役である。
この関係は対称で、**u**が**v**に対して共役なら、**v**も**u**に対して共役である（この場合の「共役」は複素共役と無関係であることに注意）。

{**p***k*}を*n*個の互いに共役なベクトル列とする。**p***k*は基底**R***n*を構成するので、**Ax** = **b**の解**x*****をこの基底で展開すると、

と書ける。ただし係数は
 \mathbf{A}\mathbf{x}_* = \sum^{n}_{i=1} \alpha_i\mathbf{p}_k^{\mathrm{T}} \mathbf{A} \mathbf{p}_i= \mathbf{p}_k^{\mathrm{T}} \mathbf{b}.
 \alpha_k = \frac{\mathbf{p}_k^{\mathrm{T}} \mathbf{b}}{\mathbf{p}_k^{\mathrm{T}} \mathbf{A} \mathbf{p}_k} = \frac{\langle \mathbf{p}_k, \mathbf{b}\rangle}{\,\,\,\langle \mathbf{p}_k,  \mathbf{p}_k\rangle_\mathbf{A}} = \frac{\langle \mathbf{p}_k, \mathbf{b}\rangle}{\,\,\,\|\mathbf{p}_k\|_\mathbf{A}^2}. }}
で与えられる。

この結果は、上で定義した内積を考えるのが最も分かりやすいと思われる。

以上から、**Ax** = **b**を解くための方法が得られる。すなわち、まず*n*個の共役な方向を見つけ、それから係数α*k*を計算すればよい。

### 反復法としての共役勾配法
共役なベクトル列**p***k*を注意深く選ぶことにより、一部のベクトルから**x***の良い近似を得られる可能性がある。そこで、共役勾配法を反復法として利用することを考える。こうすることで、*n*が非常に大きく、直接法では解くのに時間がかかりすぎるような問題にも適用することができる。

**x***の初期値を**x**0 = 0 とする。**x***が二次形式
 \mathbf{A}\mathbf{x} - \mathbf{b}^{\mathrm{T}} \mathbf{x} , \quad \mathbf{x}\in\mathbf{R}^n. }}
を最小化する一意な解であることに注意し、最初の基底ベクトル**p**1を**x** = **x**0での*f*の勾配**Ax**0&minus;**b**=**&minus;b**となるように取る。このとき、基底の他のベクトルは勾配に共役である。そこで、この方法を**共役勾配法**と呼ぶ。

**r***k*を*k*ステップ目での残差

とする。**r***k*は**x** = **x***k*での*f*の負の勾配であることに注意されたい。最急降下法は**r***k*の方向に進む解法である。**p***k*は互いに共役でなければならないので、**r***k*に最も近い方向を共役性を満たすように取る。これは
 \mathbf{A} \mathbf{r}_{k+1}}{\mathbf{p}_k^{\mathrm{T}}\mathbf{A} \mathbf{p}_k} \mathbf{p}_k }}
のように表すことができる（記事冒頭の図を参照）。

### アルゴリズム
以上の方法を簡素化することにより、**A**が実対称正定値である場合に**Ax** = **b**を解くための以下のアルゴリズムを得る。初期ベクトル**x**0は近似解もしくは0とする。
 r_0 = b - A x_0
 p_0 = r_0
 
 **for** (k = 0; ; k++) 
     \alpha_k = \frac{r_k^{\rm T} p_k}{p_k^{\rm T} A p_k}
     x_{k+1} = x_k + \alpha_k p_k
     r_{k+1} = r_k - \alpha_k A p_k
 
     **if** r_{k + 1} が十分に小さい **then**
         **break**
 
     \beta_k = \frac{r_{k+1}^{\rm T} r_{k+1}}{r_k^{\rm T} r_k}
     p_{k+1} = r_{k+1} + \beta_k p_k
 結果は x_{k+1}

### Octaveでの共役勾配法の記述例
Gnu Octaveで書くと以下のようになる。

```
function [x] = conjgrad(A,b,x0)

    r = b - A*x0;
    p = r;
    z = A*p;
    a = (r'*p)/(p'*z);
    x = x0 + a*p;
    B = 0;

    for i = 1:size(A)(1);
       r = r - a*z;
       if( norm(r) < 1e-10 )
            break;
       endif
       B = (r'*z)/(p'*z);
       p = -r + B*p;
       z = A*p;
       a = (r'*p)/(p'*z);
       x = x + a*p;
    end
 end
```


### 前処理
前処理行列とは、**A**と同値な**P**-1**A** (**P**T)-1の条件数が**A**より小さく、**Ax**=**b**より**P**-1**A**  (**P**T)-1 **x**&prime; =**P**-1**b**&prime;の方が容易に解けるような正定値行列 **P**.**P**Tを指す。前処理行列の生成には、ヤコビ法、ガウス・ザイデル法、対称SOR法などが用いられる。

最も単純な前処理行列は、**A**の対角要素のみからなる対角行列である。これはヤコビ前処理または対角スケーリングとして知られている。対角行列は逆行列の計算が容易かつメモリも消費しない点で、入門用として優れた方法である。より洗練された方法では、&kappa;(**A**)の減少による収束の高速化と**P**-1の計算に要する時間とのトレードオフを考えることになる。

## 正規方程式に対する共役勾配法
任意の実行列**A**に対して**A**T**A**は対称（半）正定値となるので、係数行列を**A**T**A**、右辺を**A**T**b**とする正規方程式を解くことにより、共役勾配法を任意の*n*×*m*行列に対して適用することができる（CGNR法）。

反復法としては、**A**T**A**を明示的に保持する必要がなく、行列ベクトル積、転置行列ベクトル積を計算すればよいので、*A*が疎行列である場合にはCGNR法は特に有効である。ただし、条件数κ(**A**T**A**)がκ(**A**2)に等しいことから収束は遅くなる傾向があり、前処理行列を使用するCGLS (Conjugate Gradient Least Squares)、LSQRなどの解法が提案されている。LSQRは**A**が悪条件である場合に最も数値的に安定な解法である。

## 脚注

### 出典

## 参考文献

*
*
* 
* 
* 

## 関連項目
*双共役勾配法 (BICG)
*非線形共役勾配法

## 外部リンク
*Conjugate Gradient Method by Nadir Soualem.
*Preconditioned conjugate gradient method by Nadir Soualem.
*An Introduction to the Conjugate Gradient Method Without the Agonizing Pain by Jonathan Richard Shewchuk.
*Iterative methods for sparse linear systems by Yousef Saad
*LSQR: Sparse Equations and Least Squares by Christopher Paige and Michael Saunders.