NDSolveの"Projection"法
はじめに
微分方程式系がある種の構造を持つ場合,数値積分法がその構造を保存すると利点が多い.ある場合は,解が制約される微分方程式を解くと役に立つことがある.投影法は,数値積分法で時間刻みを取り,実際の解が発展していく多様体上に近似解を投影することで作用する.
NDSolveには,適切であろうと思われる微分代数ソルバが含まれている.これについては「微分代数方程式の数値解法」でより詳細に記述してある.
方程式の形式は,DAEソルバによって要求される形式に還元できないことがある.さらに,いわゆる指標低減法では,微分方程式系が持ち得るシンプレクティック性等の構造特性を破壊する可能性がある([HW96]と[HLW02]を参照のこと).この例は,DAEのドキュメントに記載してある.
このような場合,制約条件が確実に保存できるように,微分方程式系を解いてから投影法を使うことができる場合が多い.これが"Projection"メソッドの考え方である.
微分方程式系が について可逆であれば,対称投影法が有利である([H00]を参照のこと).対称投影法は通常の投影よりも一般にコストがかかるため,NDSolveにはまだ実装されていない.
不変量
定義:非定数関数 は,すべての について が成り立つとき,(1)の不変量と呼ばれる.
不変量と同義語の「第1積分」,「保存量」,「運動定数」等もよく使われる.
多様体
微分方程式(3)があるなら, はすべての について ということである.これは不変性よりも弱い仮定であり,は弱い不変量といわれる([HLW02]を参照のこと).
投影アルゴリズム
が1ステップの数値積分の解を表すとする.制約条件付きの最小化問題は,以下の系を導く([AP91],[HW96],[HLW02]を参照).
作業を保存するために,は と近似される.(4)の最初の式を2つ目の式に代入すると, について以下のように簡略化されたニュートン法が導かれる.
高次の積分法を使うときに追加されるコストは,投影ステップでニュートン法の反復回数を少なくすることで相殺できる.
"Projection"メソッドの終了基準には,オプション"IterationSafetyFactor"がNDSolveによって使われる作業精度の1ULP(Unit in the Last Place)と結合される.
例題
一次不変量
従って,この例では"Projection"メソッドを使う必要はない.
数値法の中には,正確に二次不変量を保存するものもある([C87]の例を参照のこと).このようなメソッドのひとつに,陰的中点公式(あるいは1段のガウス公式,陰的ルンゲ・クッタ法)がある.
調和振動
摂動ケプラー(Kepler)問題
実験により,投影過程では使用できる不変量をすべて使うことが大切であることが分かる([HLW02]を参照).以下を使って得られた解を考えてみる.
両方の不変量に投影した解のみが,正しい質的動作をしていることが観察できる.これと比較するために,効果的なシンプレクティックソルバを使った結果が,「NDSolveのSymplecticPartitionedRungeKuttaメソッド」に記載されている.
ロトカ・ヴォルテラ系
ロトカ・ヴォルテラ系の制約射影の例は,「ロトカ・ヴォルテラ(Lotka–Volterra)方程式の数値解法」に記載されている.
オイラー方程式
オイラー方程式の制約射影の例は,「剛体ソルバ」に記載されている.
オプションの要約
オプション名
|
デフォルト値
| |
"Invariants" | None | 微分方程式系の不変量を指定する |
"IterationSafetyFactor" | 不変量の反復解法で使う安全率を指定する | |
MaxIterations | Automatic | 不変量の反復解法の最大反復回数を指定する |
Method | "StiffnessSwitching" | 微分方程式系の数値積分に使うメソッドを指定する |