アンダーソン加速

数学において、アンダーソン加速法（アンダーソン混合法とも呼ばれる）は、固定小数点反復の収束速度を加速する方法である。ドナルド・G・アンダーソン^{[1]によって導入されたこの手法は、}計算科学の分野で頻繁に生じる固定小数点方程式の解を求めるために使用できる。 ${\displaystyle f(x)=x}$

関数 が与えられたとき、方程式 の解である の固定点を求める問題を考えてみましょう。この問題に対する古典的なアプローチは、固定点反復スキーム^[2]を使用することです。つまり、解の初期推定値が与えられた場合、ある収束基準が満たされるまでシーケンスを計算します。しかし、このようなスキームの収束は一般に保証されません。さらに、収束速度は通常線形であるため、関数の評価に計算コストがかかる場合は遅くなりすぎる可能性があります。 ^[2]アンダーソン加速は、固定点シーケンスの収束を加速する方法です。^[2]${\displaystyle f:\mathbb {R} ^{n}\to \mathbb {R} ^{n}}$${\displaystyle f}$${\displaystyle f(x)=x}$${\displaystyle x_{0}}$${\displaystyle x_{i+1}=f(x_{i})}$${\displaystyle f}$

残差 を定義し、 およびと表記する（ここで は前の段落の反復シーケンスに対応する）。初期推定値と整数パラメータが与えられた場合、この方法は次のように定式化できる：^{[3] [注 1]}${\displaystyle g(x)=f(x)-x}$${\displaystyle f_{k}=f(x_{k})}$${\displaystyle g_{k}=g(x_{k})}$${\displaystyle x_{k}}$${\displaystyle x_{0}}$${\displaystyle m\geq 1}$

${\displaystyle x_{1}=f(x_{0})}$

${\displaystyle \forall k=1,2,\dots }$

${\displaystyle m_{k}=\min\{m,k\}}$

${\displaystyle G_{k}={\begin{bmatrix}g_{k-m_{k}}&\dots &g_{k}\end{bmatrix}}}$

${\displaystyle \alpha _{k}=\operatorname {argmin} _{\alpha \in A_{k}}\|G_{k}\alpha \|_{2},\quad {\text{where}}\;A_{k}=\{\alpha =(\alpha _{0},\dots ,\alpha _{m_{k}})\in \mathbb {R} ^{m_{k}+1}:\sum _{i=0}^{m_{k}}\alpha _{i}=1\}}$

${\displaystyle x_{k+1}=\sum _{i=0}^{m_{k}}(\alpha _{k})_{i}f_{k-m_{k}+i}}$

ここで、行列とベクトルの積は、は の 番目の要素です。この方法の反復を終了するには、従来の停止基準を使用できます。たとえば、が所定の許容値を下回ったとき、または残差が所定の許容値を下回ったときに反復を停止できます。 ^[2]${\displaystyle G_{k}\alpha =\sum _{i=0}^{m_{k}}(\alpha )_{i}g_{k-m_{k}+i}}$${\displaystyle (\alpha )_{i}}$${\displaystyle i}$${\displaystyle \alpha}$${\displaystyle \|x_{k+1}-x_{k}\|}$${\displaystyle g(x_{k})}$

標準的な固定小数点反復法と比較すると、この方法はより速く収束し、より堅牢であり、場合によっては固定小数点シーケンスの発散を回避できることがわかっています。^{[3] [4]}

解については、 であることが分かっており、これは と言うことと同等です。したがって、この問題は を最小化したい最適化問題として言い換えることができます。 ${\displaystyle x^{*}}$${\displaystyle f(x^{*})=x^{*}}$${\displaystyle g(x^{*})={\vec {0}}}$${\displaystyle \|g(x)\|_{2}}$

固定小数点反復法で を選択して から に直接進む代わりに、係数ベクトル、および が最後の点を含む行列である線形結合 として選択する中間点を考えてみましょう。が最小化されるように選択します。 の要素の合計は 1 になるため、1 次近似を行うことができ、問題は を最小化するを見つけることになります。 を見つけた後、原理的には を計算できます。 ${\displaystyle x_{k}}$${\displaystyle x_{k+1}}$${\displaystyle x_{k+1}=f(x_{k})}$${\displaystyle x'_{k+1}}$${\displaystyle x'_{k+1}=X_{k}\alpha _{k}}$${\displaystyle \alpha _{k}\in A_{k}}$${\displaystyle X_{k}={\begin{bmatrix}x_{k-m_{k}}&\dots &x_{k}\end{bmatrix}}}$${\displaystyle m_{k}+1}$${\displaystyle x'_{k+1}}$${\displaystyle \|g(x'_{k+1})\|_{2}}$${\displaystyle \alpha _{k}}$${\displaystyle g(X_{k}\alpha _{k})=g\left(\sum _{i=0}^{m_{k}}(\alpha _{k})_{i}x_{k-m_{k}+i}\right)\approx \sum _{i=0}^{m_{k}}(\alpha _{k})_{i}g(x_{k-m_{k}+i})=G_{k}\alpha _{k}}$${\displaystyle \alpha }$${\displaystyle \|G_{k}\alpha \|_{2}}$${\displaystyle \alpha _{k}}$${\displaystyle x'_{k+1}}$

しかし、は点を に近づけるように設計されているため、はよりも に近いと考えられます。そのため、ではなくを選択するのが合理的です。さらに、 の要素の合計は 1 になるため、1次の近似値 を作成できます。したがって、 を選択します。 ${\displaystyle f}$${\displaystyle x^{*}}$${\displaystyle f(x'_{k+1})}$${\displaystyle x^{*}}$${\displaystyle x'_{k+1}}$${\displaystyle x_{k+1}=f(x'_{k+1})}$${\displaystyle x_{k+1}=x'_{k+1}}$${\displaystyle \alpha _{k}}$${\displaystyle f(x'_{k+1})=f\left(\sum _{i=0}^{m_{k}}(\alpha _{k})_{i}x_{k-m_{k}+i}\right)\approx \sum _{i=0}^{m_{k}}(\alpha _{k})_{i}f(x_{k-m_{k}+i})=\sum _{i=0}^{m_{k}}(\alpha _{k})_{i}f_{k-m_{k}+i}}$

${\displaystyle x_{k+1}=\sum _{i=0}^{m_{k}}(\alpha _{k})_{i}f_{k-m_{k}+i}}$。

アルゴリズムの各反復において、制約付き最適化問題（ ）を解く必要がある。この問題は、いくつかの同等の定式化に書き直すことができ、[ ^3]より便利な実装につながる可能性のある異なる解法を生み出すことができる。 ${\displaystyle \operatorname {argmin} \|G_{k}\alpha \|_{2}}$${\displaystyle \alpha \in A_{k}}$

• 行列とを定義し、を解き、を設定する。^{[3] [4]}${\displaystyle {\mathcal {G}}_{k}={\begin{bmatrix}g_{k-m_{k}+1}-g_{k-m_{k}}&\dots &g_{k}-g_{k-1}\end{bmatrix}}}$${\displaystyle {\mathcal {X}}_{k}={\begin{bmatrix}x_{k-m_{k}+1}-x_{k-m_{k}}&\dots &x_{k}-x_{k-1}\end{bmatrix}}}$${\displaystyle \gamma _{k}=\operatorname {argmin} _{\gamma \in \mathbb {R} ^{m_{k}}}\|g_{k}-{\mathcal {G}}_{k}\gamma \|_{2}}$${\displaystyle x_{k+1}=x_{k}+g_{k}-({\mathcal {X}}_{k}+{\mathcal {G}}_{k})\gamma _{k}}$
• を解き、 を設定します。^[1]${\displaystyle \theta _{k}=\{(\theta _{k})_{i}\}_{i=1}^{m_{k}}=\operatorname {argmin} _{\theta \in \mathbb {R} ^{m_{k}}}\left\|g_{k}+\sum _{i=1}^{m_{k}}\theta _{i}(g_{k-i}-g_{k})\right\|_{2}}$${\displaystyle x_{k+1}=x_{k}+g_{k}+\sum _{j=1}^{m_{k}}(\theta _{k})_{j}(x_{k-j}-x_{k}+g_{k-j}-g_{k})}$

どちらの選択肢でも、最適化問題は制約のない線形最小二乗問題の形式をとり、 QR分解^[3]や特異値分解^[4]などの標準的な方法で解くことができます。また、最適化問題におけるランク不足や条件付けの問題に対処するための正則化技術も含まれる可能性があります。正規方程式を解くことで最小二乗問題を解くことは、数値不安定性の可能性や一般的に計算コストが高いため、一般的には推奨されません。^[4]

方法の停滞（つまり、同じ値、 での後続の反復）は、最小二乗問題の特異性のために、方法の破綻を引き起こします。同様に、停滞に近い状態（）は、最小二乗問題の条件付けを悪化させます。さらに、以下で説明するように、パラメータの選択は、最小二乗問題の条件付けを決定する際に関連している可能性があります。^[3]${\displaystyle x_{k+1}=x_{k}}$${\displaystyle x_{k+1}\approx x_{k}}$${\displaystyle m}$

このアルゴリズムは、可変緩和パラメータ（または混合パラメータ）を導入することで修正できる。^{[1] [3] [4]}各ステップで、新しい反復を次のように計算する。の選択は、この方法の収束特性にとって重要であり、原則として、反復ごとに変化する可能性があるが、一定になるように選択されることが多い。^[4]${\displaystyle \beta _{k}>0}$$${\displaystyle x_{k+1}=(1-\beta _{k})\sum _{i=0}^{m_{k}}(\alpha _{k})_{i}x_{k-m_{k}+i}+\beta _{k}\sum _{i=0}^{m_{k}}(\alpha _{k})_{i}f(x_{k-m_{k}+i})\;.}$$${\displaystyle \beta _{k}}$${\displaystyle \beta _{k}}$

パラメータ は、新しい反復 を計算するために以前の反復からのどれだけの情報を使用するかを決定します。一方で、 が小さすぎるように選択されると、使用される情報が少なすぎて、収束が望ましくないほど遅くなる可能性があります。一方、が大きすぎると、古い反復からの情報が後続の反復に保持されすぎて、収束が再び遅くなる可能性があります。^[3]さらに、 の選択は最適化問題のサイズに影響します。 の値が大きすぎると、最小二乗問題の条件付けとその解決コストが悪化する可能性があります。^[3]一般に、解決する特定の問題によって、パラメータの最適な選択が決まります。^[3]${\displaystyle m}$${\displaystyle x_{k+1}}$${\displaystyle m}$${\displaystyle m}$${\displaystyle m}$${\displaystyle m}$${\displaystyle m}$

上で説明したアルゴリズムに関して、各反復における の選択は変更できます。1 つの可能性は、各反復で を選択することです(切り捨てなしのアンダーソン加速と呼ばれることもあります)。^[3]この方法では、すべての新しい反復は、以前に計算されたすべての反復を使用して計算されます。より洗練された手法は、最小二乗問題に対して十分に小さい条件を維持するように を選択することに基づいています。^[3]${\displaystyle m_{k}}$${\displaystyle m_{k}=k}$${\displaystyle k}$${\displaystyle x_{k+1}}$${\displaystyle m_{k}}$

ニュートン法をの解に適用して、2 次収束するの固定点を計算することができます。しかし、この方法では の正確な導関数を評価する必要があり、非常にコストがかかります。^[4]差分法で導関数を近似することも代替案として考えられますが、反復ごとに を複数回評価する必要があり、これもまた非常にコストがかかります。アンダーソン加速法では、反復ごとに関数を 1 回評価するだけでよく、導関数の評価は必要ありません。一方、アンダーソン加速された固定点列の収束は、一般に線形のままです。^[5]${\displaystyle f(x)-x=0}$${\displaystyle f(x)}$${\displaystyle f(x)}$${\displaystyle f(x)}$${\displaystyle f(x)}$

アンダーソン加速方式と非線形方程式を解く他の方法との類似点を指摘する著者は数名いる。特に、

• Eyert ^[6]とFangとSaad ^[4]は、よく知られているセカント法を非線形方程式の解法として一般化した準ニュートン法とマルチセカント法のクラス内でアルゴリズムを解釈しました。また、この方式がBroydenクラスの方法としてどのように見えるかを示しました。^[7]${\displaystyle g(x)=0}$
• ウォーカーとニ^{[3] [8]は}、アンダーソン加速法は線形問題（つまり、ある正方行列に対する解を求める問題）の場合にGMRES法と同等であり、したがって非線形の場合にGMRES法を一般化したものと考えることができることを示した。同様の結果はワシオとオースターリーによっても得られている。^[9]${\displaystyle A\mathbf {x} =\mathbf {x} }$${\displaystyle A}$

さらに、他の著者らによって、同等またはほぼ同等の方法がいくつか独立して開発されているが^{[9] [10] [11] [12] [13]}、ほとんどの場合、固定小数点方程式の一般的な方法としてではなく、興味のある特定のアプリケーションの文脈で開発されている。

以下は、関数 の固定小数点を求めるためのアンダーソン加速方式のMATLAB言語での実装例です。次の点に注意してください。 ${\displaystyle f(x)=\sin(x)+\arctan(x)}$

• 最適化問題はQR分解を使用して解決されました。${\displaystyle \gamma _{k}=\operatorname {argmin} _{\gamma \in \mathbb {R} ^{m_{k}}}\|g_{k}-{\mathcal {G}}_{k}\gamma \|_{2}}$
• QR分解の計算は最適ではない。実際、各反復で行列に1列が追加され、場合によっては1列が削除される。この事実を利用して、より少ない計算量でQR分解を効率的に更新することができる。^[14]${\displaystyle {\mathcal {G}}_{k}}$
• 反復ベクトル全体が必要ない場合は、最新の反復と残差のみを保存することで、アルゴリズムのメモリ効率を高めることができます。${\displaystyle x_{k}}$
• コードはベクトル値のケースにそのまま一般化されます。${\displaystyle f(x)}$

f = @( x ) sin ( x ) + atan ( x ); % 固定点を計算する関数。x0 = 1 ; % 初期推定値。      
   

k_max = 100 ; % 最大反復回数。tol_res = 1e-6 ; % 残差の許容値。m = 3 ; % パラメータ m。   
   
   

x = [ x0 , f ( x0 )]; % 反復 x のベクトル。g = f ( x ) - x ; % 残差のベクトル。    
     

G_k = g ( 2 ) - g ( 1 ); % 残差の増分行列。X_k = x ( 2 ) - x ( 1 ); % xの増分行列。     
     

k = 2 ; while k < k_max && abs ( g ( k )) > tol_res m_k = min ( k , m ); % 最適化問題を QR 分解で解きます。[ Q , R ] = qr ( G_k ); gamma_k = R \ ( Q ' * g ( k )); % 新しい反復と新しい残差を計算します。x ( k + 1 ) = x ( k ) + g ( k ) - ( X_k + G_k ) * gamma_k ; g ( k + 1 ) = f ( x ( k + 1 )) - x ( k + 1 ); % 増分行列を新しい要素で更新します。X_k = [ X_k , x ( k + 1 ) - x ( k )]; G_k = [ G_k , g ( k + 1 ) - g ( k )]; n = size ( X_k , 2 ); n > m_kの場合X_k = X_k (:, n - m_k + 1 :終了); G_k = G_k (:, n - m_k + 1 :終了);終了k = k + 1 ;終了  
       
       
 
    
       
          
 
    
                
              
 
    
           
           
 
       
       
               
               
    
 
        


% 結果を出力: 9 回の反復後に計算された固定小数点 2.013444 
fprintf ( "%d 回の反復後に計算された固定小数点 %f\n" , x ( end ), k );  

• 固定小数点
• 固定小数点反復
• 根を求めるアルゴリズム