NDSolveの合成・分割法
はじめに
微分方程式系を部分系に分割し，それぞれを適切な積分法を用いて解くと便利なことがある．それぞれの解を合成することで，体積等の動的特性を保存することができることがある．分割と合成に関する詳細は[
MQ02,
HLW02]に，
NDSolveに関連する特定の問題は[
SS05,
SS06]に記載されている．
定義
y (0)=y0
nであるときの初期値問題
y' (t)=f (y (t))について考える．
合成
合成は，数値積分法の次数を上げるための便利な手法である．
外挿で使われるAitken-Nevilleアルゴリズムとは対照的に，合成ではシンプレクティック法等の基礎的な積分法の幾何的特性を保存することができる．
を，指定された実数
1, ..., sを持ち，刻み幅
ih を取る基本的積分法とする．
すると，
s 段階の合成法
f, h は以下のように与えられる．
同じ基本的メソッド
= (i), i=1, ..., s.を含む合成法
f, h に興味があることがある．
面白い特殊ケースとして，対称合成
i=s-i+1（ここで
i=1, ..., 
s/2
）がある．
- 1次法（随伴行列*を持つメソッド）の対称合成
分割
s 段階の分割法は，加法的に分割された
f の合成法の一般化である．
大切なのは，
fではなく
fiを含む問題を解く際に，計算上の利点がある場合が多いということである．
ここで
f1, ..., fs は必ずしも異なるものでなくてもよい．
これでそれぞれの基本的積分法が，問題の一部だけを解くようになったが，適切な合成により有利な特性を持つ数値スキームが生まれる可能性がまだある．
べクトル場
fi が積分可能ならば，厳密解もしくは流れ
fi, h を数値積分法の代りに使うことができる．
分割法では，流れと数値法を混合して使うこともできる．
その例にLie-Trotterの分割[
T59]がある．
f=f1+f2を
1=2=1で分割する．すると
が1次の積分法を生成する．
計算上では，集団特性を使って流れを合成すると利点が多い．
実装
分割法と合成法を実装するために，新しい
NDSolveフレームワークにいくつかの変更を加えることが必要であった．
- メソッドにより，任意数の部分メソッドが呼び出せるようにする．
- べクトル場全体ではなく部分体を数値的に評価するための関数を順に渡していく機能．
- 流れを計算するLocallyExactメソッド．部分系を分析的に解いて（局所的な）解を数値的に求める．
- 初期化の繰返しを避けるために同一メソッドのデータをキャッシュに格納．同一の部分体を数値的に評価するためのデータもキャッシュに格納される．
簡略化された入力文法により，省略されたべクトル場とメソッドを循環的に挿入することができる．これらは明確に定義しなければならない．
{f1, f2, f1, f2}は
{f1, f2}で入力できる．
{f1, f2, f3, f2, f1}は
{f1, f2, f3}で入力できない．
{f1, f2, f3, f2, f1}は
{f1, f2, f3, f1, f2}に対応するからである．
ネストされたメソッド
以下の例では合成を使って，低次の分割法から高次の分割法を構築する．
簡略化
上記の例で，流れの集合特性を使い，
Splittingメソッドを直接呼び出すと，より効率的な積分法を得ることができる．
例題
次の例題は，ストラング(Strang)の分割[
S68]，[
M68]として知られる2次の対称分割を使う．分割係数は，方程式の構造から自動的に決定される．
メソッド SymplecticPartitionedRungeKuttaについてシンプレクティック法のリープフロッグ公式として知られるメソッドを定義する． |
シンプレクティック分割
シンプレクティックリープフロッグ
SymplecticPartitionedRungeKuttaは長期力学を使った，分離可能なハミルトン系
H (p, q)=T (p)+V (q)を解くための効率的なメソッドである．
Splittingの方が汎用のメソッドであるが，分割されたシンプレクティック法を構築するために使うことができる（しかし，
SymplecticPartitionedRungeKuttaよりもやや効率性に劣る）．
分離可能なハミルトン系 H (p, q)=p2/2+q2/2によって支配される線形微分方程式系から派生する調和振動を例に取る．
| Out[5]= |  |
|
ハミルトンべクトル場を，運動量および位置を支配する独立の成分に分割する．これは，方程式の右辺をゼロに設定することで実行できる． |
この重み付き（1次）オイラー積分ステップの合成は，シンプレクティック（2次）リープフロッグ公式に当たる．
| Out[14]= |  |
|
ExplicitEulerメソッドは1度しか指定できなかった．2回目と3回目は循環的に埋められるからである．
Out[15]//InputForm= |
| |  |
|
シンプレクティックリープフロッグ積分法に対応する組込みの積分法を実行する．
| Out[16]= |  |
|
積分ステップの最後の結果は，分割法の結果と同一である．
Out[17]//InputForm= |
| |  |
|
シンプレクティックリープフロッグの合成
これは，基本的な積分法としてシンプレクティック法のリープフロッグ公式を使い，Ruth-Yoshidaの対称合成[ Y90]を使って4次のシンプレクティック積分法を構築する．
| Out[20]= |  |
|
Out[21]//InputForm= |
| |  |
|
これはRuth，Yoshidaの4次メソッドの係数を使って組込みのシンプレクティック積分法を実行する．
| Out[23]= |  |
|
積分ステップの最後の結果は，合成法の結果と同一である．
Out[24]//InputForm= |
| |  |
|
複合メソッド
閉形式の解はべクトル場全体に存在しないことがよくあるが，ベクトル場の一部について微分方程式系を解析的に解くことが可能な場合がある．
DSolveで解が見付かる場合は，局所的に解を求めるのに直接数値評価を使うことができる．
この考えはメソッド
LocallyExactに実装されている．
調和振動テスト例題
以下の例題は，調和振動方程式の分割成分の厳密な流れの解が，各分割成分にオイラー法を適用したものと同じであることを調べる．
Out[34]//InputForm= |
| |  |
Out[38]//InputForm= |
| |  |
|
数値・記号複合分割法（ABC flow）
3次元流れを体積保存するためのモデルで，広く研究されているABC flow (Arnold，Beltrami，Childress flow)を考えてみる．
| Out[39]= |  |
|
体積保存積分法を直接適用すると，一般に対称性は保存しない．陰的中点公式等の対称保存積分法では，体積を保存しない．
右辺の要素を零に設定することで系の分割を定義する．
| Out[45]= |  |
| Out[46]= |  |
|
Y1系は
DSolveで厳密に解くことができるので，
LocallyExactメソッドを使うことができる．
しかし，
Y2は解けないので，分割法の目標の属性を得るためには，適切な数値積分法を使う必要がある．
組込みの ImplicitRungeKuttaメソッドについて陰的中点公式を計算するメソッドを定義する． |
2次の体積保存・反転対称群積分法を定義する[ MQ02]．
| Out[48]= |  |
|
ロトカ・ボルテラ(Lotka-Volterra)方程式
この系のさまざまな数値積分法が
「ロトカ・ボルテラ方程式の数値解法」の中で比較してある．
オイラー方程式
オイラー方程式のさまざまな数値積分法が，
「剛体ソルバ」の中で比較してある．
非自律系のべクトル場
強制平面非自律微分方程式系のDuffing振動子を考える．
| Out[49]= |  |
|
べクトル場の時間成分の分割はあいまいであるため，メソッドはエラーメッセージを発する．
| Out[52]= |  |
|
Z[T] Tとなるような新しいダミー変数 Z[T]を導入し，分割された微分方程式系でそれをどのように分配するかを指定することにより，方程式を拡張することができる． |
これは，あらゆる有限時間に対して  1+2=- （ここで 1と 2はリアプノフ指数）を満足する幾何的分割法を定義する[ MQ02]．
| Out[59]= |  |
|
| Out[60]= |  |
|
オプションの要約
Compositionのデフォルトの係数オプションは，
Methodオプションを使って指定するメソッドの特性によって，
SymmetricCompositionCoefficientsと
SymmetricCompositionSymmetricMethodCoefficientsのどちらかを自動的に選ぶことを試みる．
Compositionメソッドのオプション
| | |
| Coefficients | {} | 分割法で使用する係数を指定する |
| DifferenceOrder | Automatic | メソッドの局所的な精度の次数を指定する |
| Equations | {} | 方程式を分割する方法を指定する |
| Method | None | 数値積分で使用する基礎メソッドを指定する |
Splittingメソッドのオプション
