3つの方法、3つの領域
テスト方程式 dy/dx = λy に対して、3 つの明示的な ODE メソッドは複素数 hλ 平面内の以下の安定性領域を持っています：
オイラー法（1 次）：安定性領域はディスク |1 + hλ| ≤ 1 です。これは (-1, 0) を中心とした半径 1 の円です。実数の負の hλ は [-2, 0] に含まれる必要があります。
ルンゲ-クッタ 2（中点法）（2 次）：安定性領域は |1 + hλ + (hλ)²/2| ≤ 1 です。オイラーのディスクより大きいですが、依然として有限です。
ルンゲ-クッタ 4（4 次）：安定性領域は |1 + hλ + (hλ)²/2 + (hλ)³/6 + (hλ)⁴/24| ≤ 1 を満たします。実数の負の hλ は約 -2.785 まで拡張されます。この領域はオイラーのものより実質的に大きいです。
後退オイラー（陰解法）：安定性領域は、ディスク |1 - hλ|⁻¹ > 1（同等に |1/(1-hλ)| ≤ 1）以外の複素平面全体です。λ が左半平面にある場合（Re(λ) < 0）、これは無条件に安定しています — h に対する安定性からの制約はありません。
増幅関数
任意のルンゲ-クッタ法について、ステップごとの増幅係数 R(hλ) は e^(hλ) への多項式近似です：
- オイラー：R(z) = 1 + z（次数 1 で切り捨て）
- RK2：R(z) = 1 + z + z²/2（次数 2 で切り捨て）
- RK4：R(z) = 1 + z + z²/2 + z³/6 + z⁴/24（次数 4 で切り捨て）
安定性領域は {z ∈ ℂ : |R(z)| ≤ 1} です。真の解の増幅：|e^z| = e^(Re(z))。Re(z) < 0（安定な ODE）の場合、真の解は減衰します。数値メソッドは |R(z)| ≤ 1 の場合に安定です — 減衰挙動と一致しています。
純虚数固有値：振動系
多くの物理系は純虚数固有値を持っています：λ = iω（減衰なしの振動）。ばね-質量系、軌道力学、振り子の動力学。
λ = iω の場合：hλ = ihω は虚軸上にあります。
虚軸上のオイラーの安定性： |1 + ihω|² = 1 + (hω)² > 1 は任意の h > 0 に対して成立します。オイラーは純虚数固有値に対する任意のステップサイズで不安定です。計算された「振動」は限界なく成長します。
虚軸上の RK4 の安定性： 安定性領域は虚軸上で約 |hω| ≤ 2.83 まで拡張されます。十分に小さい h では、RK4 は減衰なしの振動を処理できます。オイラーはできません。
この幾何学は、オイラーが保存的システム（ばね-質量、軌道、波動方程式）で小さい h でも失敗する一方、RK4 はそれらを適切に処理する理由です。
硬い問題の幾何学
硬い ODE システムは非常に異なる大きさの固有値を持っています。硬さ比：κ = max|Re(λᵢ)| / min|Re(λᵢ)| >> 1。
なぜ硬さは明示的なソルバーに対して高くつくのか：
安定性は h·max|λᵢ| ≤ C を要求します（ここで C はメソッドに依存します）。最も負の固有値が上限を設定します。
遅い動力学の精度は h·min|λᵢ| ≥ ε を要求します（最も遅いモードを十分に解決するため）。
κ が大きい場合、これら 2 つの要件は小さな h を強制します：速いモードの安定性に対して十分に小さく、遅いモードをサンプリングするのに十分に大きい。ステップ数は κ でスケールします。
固有値スペクトルの幾何学的描像： ヤコビアン ∂f/∂y の固有値は複素平面内の点の集合を形成します。明示的なソルバーの安定性領域はすべての点 h·λᵢ を含む必要があります。固有値が -1 から -1000 まで及ぶ場合、安定性領域は実軸に沿って 1000 の範囲をカバーする必要があり、非常に小さな h が必要です。
陰解法： 後退オイラーの安定性領域は左半平面全体をカバーします。Re(λ) < 0 のすべての固有値は、h に関係なく安定性領域内に自動的に含まれます。h の制約は安定性ではなく精度からのみ来ます。
硬さ比 & コスト
高速反応（時間スケール 10⁻⁶ s）と遅い反応（時間スケール 1 s）を持つ化学反応ネットワークを考えてください。
硬さ比：κ = 10⁶ / 1 = 10⁶。
RK4（安定性限界 h·|λ_max| ≤ 2.785）を使用：h_max = 2.785 / 10⁶ ≈ 2.8 × 10⁻⁶ s。
10 s の反応時間にわたって積分するには：ステップ = 10 / (2.8 × 10⁻⁶) ≈ 3.6 × 10⁶。
後退オイラー（無条件に安定）を使用：h は遅い反応の精度のために選択できます。h = 10⁻² s（1 s スケールにわたって 100 サンプル）。ステップ = 10 / 10⁻² = 1000。
コスト比：明示的な 360 万ステップ vs 陰的な 1000 ステップ — 3600 の係数。各陰的ステップは線形システムを解く必要があります（ステップごとのコストはより高い）。ただし、非常に硬い問題の場合、総コストはずっと低いです。
n 次元チューブが期待と異なる理由
2D では、曲線 C の周りの半径 ε の「チューブ」は、C から距離 ε 以内の点の集合です。断面は半径 ε の円です。チューブの体積はその長さに比例して成長します。
n 次元では、第 9 章の現象のため、チューブの幾何学は根本的に変わります：
n 次元のコーナーパラドックス： n 次元空間では、n 次元超立方体の体積のほぼすべてはコーナーにあります — 中央領域にはありません。n が増加するにつれて、中心から距離 ε 以内の体積の割合は、任意の固定 ε に対してゼロになります。
ODE 解チューブに適用：
2D では：真の解がチューブの中心を通過する場合、ほとんどの近い点は曲線に近いです。小さな摂動は真の解の近くに留まります。
高次元では：チューブのバウンディングボックス内のほとんどの点は、実際には真の解曲線から遠いです。チューブの「体積」はコーナーによって支配されます — 複数の次元で中心から遠い領域。
シミュレーションへの結果： 28 つの結合した ODE（ハミング海軍傍受問題）を使用して、各次元でサイズ ε の摂動は、真の解から ε√28 ≈ 5.3ε の総変位を生じることができます。チューブは、1 つの次元での最大変位ではなく、すべての次元にわたる L2 ノルムの観点から理解する必要があります。
高次元での安定性： 各成分が独立して減衰するシステム（各固有値は負の実部を持つ）でも、成分のエラーが L2 ノルムで追加されるため、大きな結合変位を示すことができます。28 次元チューブは単に 28 の独立した 1 次元チューブではありません — 幾何学はそれらをカップリングします。
幾何学からデザインへ
第 18-20 章の幾何学的洞察は、数値シミュレーションの設計原則のセットとして一緒に来ます：
ステップサイズ選択： h はすべての固有値に対して安定性領域内に h·λ を配置する必要があります。硬いシステムの場合、陰的メソッドは安定性制約を削除し、精度要件のみを残します。
高次元でのエラー累積： 全体的エラーは n 次元空間内のベクトルです。そのノルムはコンポーネントごとのエラーの √n 倍として成長します。高次元シミュレーションはより厳しいコンポーネントごとの精度要件が必要です。
フィードバックとしての安定化装置： シミュレーションがフィードバックを組み込む場合（計算された出力は後続の入力に影響を与える、たとえば誘導システム）、収束フィードバックはエラーをダンプします。シミュレーションはフィードバックループ内の量に対して不正確な入力を許容できます。
信号としての不安定性： 発散方向場を持つ問題の場合、不安定性を利用できます：発散の方向は初期条件エラーに関する情報を運び、補正調整を可能にします。