Homotopy Method

解きたい非線形方程式を $\boldsymbol{f}(\boldsymbol{x}) = \boldsymbol{0}$ とします。ホモトピー法は、自明な解を持つ補助方程式 $\boldsymbol{g}(\boldsymbol{x}) = \boldsymbol{0}$ を用意し、連続パラメータ $\lambda \in [0, 1]$ を導入して、以下の一般式（ホモトピー写像）を構成します。

この $H(\boldsymbol{x}, \lambda)$ の零点集合が、$\lambda=0$ の既知の解から $\lambda=1$ の求めるべき解へと繋がる「解曲線」を形成します。補助方程式 $\boldsymbol{g}(\boldsymbol{x})$ の選び方により、以下の代表的な手法に分類されます。

① 不動点ホモトピー (Fixed-Point Homotopy)

$\boldsymbol{g}(\boldsymbol{x}) = \boldsymbol{x} - \boldsymbol{x}_0$ とする最も標準的な形式。$\lambda=0$ において唯一の解 $\boldsymbol{x} = \boldsymbol{x}_0$ を持ち、そこから探索を開始します。

\[ H(\boldsymbol{x}, \lambda) = (1 - \lambda)(\boldsymbol{x} - \boldsymbol{x}_0) + \lambda \boldsymbol{f}(\boldsymbol{x}) = \boldsymbol{0} \]

② ニュートンホモトピー (Newton Homotopy)

$\boldsymbol{g}(\boldsymbol{x}) = \boldsymbol{f}(\boldsymbol{x}) - \boldsymbol{f}(\boldsymbol{x}_0)$ と置く形式。初期値 $\boldsymbol{x}_0$ における残差を強制的にゼロにするようなオフセットを与え、それを徐々に取り除いていく過程に相当します。

\[ H(\boldsymbol{x}, \lambda) = \boldsymbol{f}(\boldsymbol{x}) - (1 - \lambda)\boldsymbol{f}(\boldsymbol{x}_0) = \boldsymbol{0} \]

③ 奇ホモトピー (Odd Homotopy)

$\boldsymbol{g}(\boldsymbol{x})$ に $\boldsymbol{f}(\boldsymbol{x})$ の奇関数部分を利用する形式。解の存在証明や、複数の解を探索する際の位相幾何学的な安定性に優れます。

\[ H(\boldsymbol{x}, \lambda) = \lambda \boldsymbol{f}(\boldsymbol{x}) + (1 - \lambda) \frac{\boldsymbol{f}(\boldsymbol{x}) - \boldsymbol{f}(-\boldsymbol{x})}{2} = \boldsymbol{0} \]

解曲線 $H(\boldsymbol{x}, \lambda) = \boldsymbol{0}$ を追跡する際、$\lambda$ を単純な増分パラメータとして扱うと、曲線が「折り返し点（Turning Point）」を持つ場合にヤコビ行列が特異となり、計算が破綻します。これを回避するため、解曲線の道のり（弧長）$s$ を独立変数とする弧長法 (Arc-length Continuation) を用います。

① 拡大方程式の構成

変数ベクトルを $\boldsymbol{y} = [\boldsymbol{x}^T, \lambda]^T$ とし、$n+1$ 次の拡大方程式を解きます。弧長 $s$ に対する接線ベクトルを $\dot{\boldsymbol{y}} = [d\boldsymbol{x}/ds, d\lambda/ds]^T$ とすると、以下の拘束条件が成り立ちます。

② 予測子・修正子サイクル

• 接線予測 (Predictor): 現ステップの接線ベクトルを $\boldsymbol{J}_H \dot{\boldsymbol{y}} = \boldsymbol{0}$ から算出。オイラー法等で $\boldsymbol{y}_{pred} = \boldsymbol{y}_k + \Delta s \cdot \dot{\boldsymbol{y}}_k$ を求めます。
• 垂直修正 (Corrector): 予測点から接線ベクトルに直交する超平面上で、ニュートン法により解曲線へ引き戻します。このとき、ヤコビ行列に接線条件を加えた $n+1$ 次の拡大ヤコビ行列は、折り返し点においても常に正則となり、安定した収束を保証します。

「円」と「直線」の交点を求める以下の非線形系を考えます。ニュートン法では初期値が円の外側にあると収束が不安定になりますが、ホモトピー法では確実に解へ導かれます。

Step 1: ホモトピーの構築
初期値 $\boldsymbol{x}_0 = [0, 2]^T$（円の外側）を選択し、不動点ホモトピーを適用します。

\[ H(\boldsymbol{x}, \lambda) = (1-\lambda) \begin{bmatrix} x_1 - 0 \\ x_2 - 2 \end{bmatrix} + \lambda \begin{bmatrix} x_1^2 + x_2^2 - 1 \\ x_1 - x_2 \end{bmatrix} = \boldsymbol{0} \]

Step 2: 解曲線の挙動
$\lambda=0$ のとき、解は自明に $(0, 2)$ です。ここから弧長 $s$ を進めると、$\lambda$ の増加に伴い解は円の境界へと近づきます。途中のヤコビ行列 $\partial H / \partial \boldsymbol{x}$ は以下のようになります：

\[ \boldsymbol{J} = (1-\lambda) \begin{bmatrix} 1 & 0 \\ 0 & 1 \end{bmatrix} + \lambda \begin{bmatrix} 2x_1 & 2x_2 \\ 1 & -1 \end{bmatrix} \]

Step 3: 収束
最終的に $\lambda=1$ に到達したとき、方程式は $\boldsymbol{f}(\boldsymbol{x}) = \boldsymbol{0}$ に完全に一致し、解 $\boldsymbol{x} = [1/\sqrt{2}, 1/\sqrt{2}]^T$ を得ます。ニュートン法が「跳ねて」しまうような強い非線形領域も、弧長法によって曲線を正確にトレースすることで回避可能です。

ホモトピー法の実装において最も重要なのは、解曲線の追跡に失敗した際の「適応的ステップ制御」です。曲率 $d^2\boldsymbol{y}/ds^2$ を監視し、急峻な変化がある箇所では $\Delta s$ を極小化し、平坦な箇所では大きく取ることで、計算時間と堅牢性を両立させます。