交互方向陰解法

数値線形代数において、交互方向陰的法（ADI法）は、シルベスター行列方程式を解くために使用される反復法である。これは、システム理論や制御で生じる大規模な行列方程式を解くための一般的な方法であり、^[1]メモリ効率の良い因数分解形式で解を構築するように定式化することができる。^{[2] [3]また、}放物型および楕円型の偏微分方程式を数値的に解くために使用され、熱伝導をモデル化したり、 2次元以上の拡散方程式を解いたりするために使用される古典的な方法である。 ^{[4]これは、}演算子分割法の一例である。^[5]

ADI法は、近似解の列空間と行空間を交互に更新する2段階の反復プロセスです。1回のADI反復は、次の手順で構成されます。^[6]${\displaystyle AX-XB=C}$

1. を解きます。ここで${\displaystyle X^{(j+1/2)}}$${\displaystyle \left(A-\beta _{j+1}I\right)X^{(j+1/2)}=X^{(j)}\left(B-\beta _{j+1}I\right)+C.}$

2. を解きます。ここで です 。${\displaystyle X^{(j+1)}}$${\displaystyle X^{(j+1)}\left(B-\alpha _{j+1}I\right)=\left(A-\alpha _{j+1}I\right)X^{(j+1/2)}-C}$

これらの数値はシフトパラメータと呼ばれ、収束はこれらのパラメータの選択に大きく依存します。^{[7] [8]} ADIの反復を実行するには、初期推定値とシフトパラメータが必要です。 ${\displaystyle (\alpha _{j+1},\beta _{j+1})}$${\displaystyle K}$${\displaystyle X^{(0)}}$${\displaystyle K}$${\displaystyle \{(\alpha _{j},\beta _{j})\}_{j=1}^{K}}$

および の場合、 はBartels-Stewart 法 を用いてで直接解くことができます。^[9]したがって、行列ベクトル乗算と および を含む線形解が安価に適用できる場合にのみ、ADI を使用するのが有益です。 ${\displaystyle A\in \mathbb {C} ^{m\times m}}$${\displaystyle B\in \mathbb {C} ^{n\times n}}$${\displaystyle AX-XB=C}$${\displaystyle {\mathcal {O}}(m^{3}+n^{3})}$${\displaystyle A}$${\displaystyle B}$

方程式は、 のときのみ一意に解を持ちます。ここで、は のスペクトルです。^[1]しかし、ADI 法は、と が十分に分離されており、 とが正規行列である場合に特に優れたパフォーマンスを発揮します。これらの仮定は、たとえば、が正定値である場合のLyapunov 方程式によって満たされます。これらの仮定の下では、およびのいくつかの選択に対して、ほぼ最適なシフトパラメータがわかっています。^{[7] [8]}さらに、事前誤差境界を計算できるため、実装時に残差誤差を監視する必要がなくなります。 ${\displaystyle AX-XB=C}$${\displaystyle \sigma (A)\cap \sigma (B)=\emptyset }$${\displaystyle \sigma (M)}$${\displaystyle M}$${\displaystyle \sigma (A)}$${\displaystyle \sigma (B)}$${\displaystyle A}$${\displaystyle B}$${\displaystyle AX+XA^{*}=C}$${\displaystyle A}$${\displaystyle A}$${\displaystyle B}$

ADI法は、上記の仮定が満たされない場合でも適用できます。最適でないシフトパラメータの使用は収束に悪影響を与える可能性があります。 ^[1]また、収束はまたはの非正規性によっても影響を受けます（時には有利に）。^[10]有理クリロフ部分空間法などのクリロフ部分空間法は、この設定では通常ADIよりも急速に収束することが観察されています。^{[1] [3]これ}により、ハイブリッドADI投影法が開発されました。^[3]${\displaystyle A}$${\displaystyle B}$

適切なシフトパラメータを見つける問題は簡単ではありません。この問題は、ADI誤差方程式を調べることで理解できます。反復後、誤差は次のように表されます。 ${\displaystyle K}$

${\displaystyle X-X^{(K)}=\prod _{j=1}^{K}{\frac {(A-\alpha _{j}I)}{(A-\beta _{j}I)}}\left(X-X^{(0)}\right)\prod _{j=1}^{K}{\frac {(B-\beta _{j}I)}{(B-\alpha _{j}I)}}.}$

を選択すると、相対誤差の境界は次のようになります。 ${\displaystyle X^{(0)}=0}$

${\displaystyle {\frac {\left\|X-X^{(K)}\right\|_{2}}{\|X\|_{2}}}\leq \|r_{K}(A)\|_{2}\|r_{K}(B)^{-1}\|_{2},\quad r_{K}(M)=\prod _{j=1}^{K}{\frac {(M-\alpha _{j}I)}{(M-\beta _{j}I)}}.}$

ここで は演算子ノルムです。シフトパラメータの理想的なセットは、量 を最小化する有理関数を定義します。と が正規行列であり、固有分解と を持つ場合、 ${\displaystyle \|\cdot \|_{2}}$${\displaystyle \{(\alpha _{j},\beta _{j})\}_{j=1}^{K}}$ ${\displaystyle r_{K}}$${\displaystyle \|r_{K}(A)\|_{2}\|r_{K}(B)^{-1}\|_{2}}$${\displaystyle A}$${\displaystyle B}$ ${\displaystyle A=V_{A}\Lambda _{A}V_{A}^{*}}$${\displaystyle B=V_{B}\Lambda _{B}V_{B}^{*}}$

${\displaystyle \|r_{K}(A)\|_{2}\|r_{K}(B)^{-1}\|_{2}=\|r_{K}(\Lambda _{A})\|_{2}\|r_{K}(\Lambda _{B})^{-1}\|_{2}}$。

ほぼ最適なシフトパラメータは、および （ただし、 および は実数直線上の互いに素な区間）の場合など、特定のケースでは既知である。^{[7] [8]}たとえば、リャプノフ方程式 は、が正定値の場合にこれらの仮定を満たします。この場合、シフトパラメータは楕円積分 を使用して閉じた形式で表すことができ、数値的に簡単に計算できます。 ${\displaystyle \Lambda _{A}\subset [a,b]}$${\displaystyle \Lambda _{B}\subset [c,d]}$${\displaystyle [a,b]}$${\displaystyle [c,d]}$ ${\displaystyle AX+XA^{*}=C}$${\displaystyle A}$

より一般的には、閉じた互いに素な集合と（ただし 、 と）が既知である場合、最適なシフトパラメータ選択問題は、次の値を達成する極値有理関数を見つけることによって近似的に解かれる。 ${\displaystyle E}$${\displaystyle F}$${\displaystyle \Lambda _{A}\subset E}$${\displaystyle \Lambda _{B}\subset F}$

${\displaystyle Z_{K}(E,F):=\inf _{r}{\frac {\sup _{z\in E}|r(z)|}{\inf _{z\in F}|r(z)|}},}$

ここで、下限は次数 のすべての有理関数に対して取られる。^[8]この近似問題は、ポテンシャル理論におけるいくつかの結果と関連しており、^{[12] [13]} 、 = [a, b] および ^{[14]に対して}ゾロタレフによって1877 年に解決された。また、 と が複素平面上の互いに素な円板である場合の解も知られている。 ^[15]${\displaystyle (K,K)}$${\displaystyle E}$${\displaystyle F=-E.}$${\displaystyle E}$${\displaystyle F}$

とについてあまり知られていない場合、またはまたは が非正規行列である場合、最適に近いシフトパラメータを見つけることができない可能性があります。 このような設定では、適切なシフトパラメータを生成するためのさまざまな戦略を使用できます。 これらには、ポテンシャル理論の漸近結果に基づく戦略、^[16 ]、行列、、、のリッツ値を使用して貪欲なアプローチを定式化する戦略、^[17]、および収束許容値が満たされるまでシフトパラメータの同じ小さなコレクションを再利用する巡回法が含まれます。^{[17] [10]}すべての反復で同じシフトパラメータが使用される場合、ADI はスミス法と呼ばれるアルゴリズムと同等です。^[18]${\displaystyle \sigma (A)}$${\displaystyle \sigma (B)}$${\displaystyle A}$${\displaystyle B}$${\displaystyle A}$${\displaystyle A^{-1}}$${\displaystyle B}$${\displaystyle B^{-1}}$

多くのアプリケーションでは、と は非常に大きな疎行列であり、（、 ）として因数分解できます。^[1] このような設定では、潜在的に密な行列を明示的に保存することは現実的ではありません。 ADI の変形である因数分解 ADI ^{[3] [2]}を使用すると、 （ ）を計算することができます。因数分解 ADI の有効性は、 が低ランク行列で十分に近似されるかどうかに依存します。 これは、とに関するさまざまな仮定の下で真であることが知られています。^{[10] [8]}${\displaystyle A}$${\displaystyle B}$${\displaystyle C}$${\displaystyle C=C_{1}C_{2}^{*}}$${\displaystyle C_{1}\in \mathbb {C} ^{m\times r},C_{2}\in \mathbb {C} ^{n\times r}}$${\displaystyle r=1,2}$${\displaystyle X}$${\displaystyle ZY^{*}}$${\displaystyle X\approx ZY^{*}}$${\displaystyle X}$${\displaystyle A}$${\displaystyle B}$

歴史的に、ADI法は有限差分法を用いて正方領域上の2次元拡散方程式を解くために開発されました。^[4]行列方程式のADIとは異なり、放物型方程式のADIでは、各反復で現れるシフトがタイムステップ、拡散係数、グリッド間隔などのパラメータによって決まるため、シフトパラメータの選択は必要ありません。行列方程式のADIとの関連は、定常状態のシステムに対するADI反復の作用を考えるとわかります。

熱伝導方程式を数値的に解く従来の方法は、クランク・ニコルソン法です。この方法では、多次元の非常に複雑な方程式のセットが生成され、解くのにコストがかかります。ADI法の利点は、各ステップで解く必要のある方程式の構造がより単純で、三角行列アルゴリズムを使用して効率的に解くことができることです。

2次元の線形拡散方程式を考えてみましょう。

${\displaystyle {\partial u \over \partial t}=\left({\partial ^{2}u \over \partial x^{2}}+{\partial ^{2}u \over \partial y^{2}}\right)=(u_{xx}+u_{yy})}$

暗黙的なクランク・ニコルソン法では、次の差分方程式が生成されます。

${\displaystyle {u_{ij}^{n+1}-u_{ij}^{n} \over \Delta t}={1 \over 2(\Delta x)^{2}}\left(\delta _{x}^{2}+\delta _{y}^{2}\right)\left(u_{ij}^{n+1}+u_{ij}^{n}\right)}$

どこ：

${\displaystyle \Delta x=\Delta y}$

はp番目の座標 の中心となる2階差分演算子である。${\displaystyle \delta _{p}^{2}}$

${\displaystyle \delta _{p}^{2}u_{ij}=u_{ij+e_{p}}-2u_{ij}+u_{ij-e_{p}}}$

またはの場合はそれぞれまたは です(格子点 の省略形)。 ${\displaystyle e_{p}=(1,0)}$${\displaystyle (0,1)}$${\displaystyle p=x}$${\displaystyle y}$${\displaystyle ij}$${\displaystyle (i,j)}$

安定性解析を実行した後、この方法は任意の に対して安定していることが示されます。 ${\displaystyle \Delta t}$

クランク・ニコルソン法の欠点は、上記の式の行列が通常かなり大きなバンド幅でバンド化されることです。これにより、線形方程式の直接解法は非常にコストがかかります (ただし、不完全コレスキー分解を前提とする共役勾配法の使用など、効率的な近似解法は存在します)。

ADI法の考え方は、差分方程式を2つに分割し、1つはx微分を暗黙的に取り、もう1つはy微分を暗黙的に取ります。

${\displaystyle {u_{ij}^{n+1/2}-u_{ij}^{n} \over \Delta t/2}={\left(\delta _{x}^{2}u_{ij}^{n+1/2}+\delta _{y}^{2}u_{ij}^{n}\right) \over \Delta x^{2}}}$

${\displaystyle {u_{ij}^{n+1}-u_{ij}^{n+1/2} \over \Delta t/2}={\left(\delta _{x}^{2}u_{ij}^{n+1/2}+\delta _{y}^{2}u_{ij}^{n+1}\right) \over \Delta y^{2}}}$

関係する方程式系は対称かつ三重対角（帯域幅 3 でバンド化）であり、通常は三重対角行列アルゴリズムを使用して解かれます。

この方法は無条件に安定しており、時間と空間において2次であることが示されています。^[19] ダグラス法^[20]やf因子法^[21]など、3次元以上に使用できる、より洗練されたADI法があります。

ADI法を演算子分割法として用いることは一般化できる。つまり、一般的な発展方程式を考えることができる。

${\displaystyle {\dot {u}}=F_{1}u+F_{2}u,}$

ここで、およびはバナッハ空間上で定義された（非線形の可能性のある）演算子です。^{[22] [23]}上記の拡散の例では、およびが存在します。 ${\displaystyle F_{1}}$${\displaystyle F_{2}}$${\displaystyle F_{1}={\partial ^{2} \over \partial x^{2}}}$${\displaystyle F_{2}={\partial ^{2} \over \partial y^{2}}}$

従来のADI法を、左辺にのみ類似の演算子を持ち、右辺には演算子がない基本ADI法に簡略化することができます。これは、方程式の両辺に演算子が含まれる従来のほとんどの陰的解法とは異なり、右辺に（簡約される）演算子がなくなったADI法の基本（基本）スキームと見なすことができます^{[24] [25]。FADI} 法は、従来のADI法の精度を低下させることなく、より単純で簡潔で効率的な更新方程式をもたらします。

ピースマン・ラッチフォード、ダグラス・ガン、ディアコノフ、ビーム・ウォーミング、クランク・ニコルソンなどによる多くの古典的な陰解法は、演算子のない右辺を持つ基本的な陰解法に簡略化できる。^{[25]基本的な形式では、2次の時間精度のFADI法は、計算電磁気学における3次元マクスウェル方程式[ 26] [27]}など、2次の時間精度にアップグレードできる基本的な局所1次元（FLOD）法と密接に関連している。2次元および3次元の熱伝導および拡散方程式の場合、FADI法とFLOD法はどちらも、従来の方法と比較して、より単純で効率的で安定した方法で実装できる。^{[28] [29]}