微分方程式の数値解法
「入門：微分方程式の数値解法」で説明した
NDSolveを使うと，微分方程式を数値解析的に解くことができる．また，単独の微分方程式でも連立の微分方程式でも
NDSolveで解ける．さらに，広い範囲の常微分方程式と偏微分方程式の一部を扱える．常微分方程式系の場合，未知の関数
yi がいくつあっても構わないが，すべて単一の独立変数
x のみに依存しているものとする．偏微分方程式の場合は，複数の独立変数を含んでいる．
NDSolveは微分方程式と代数方程式が混ざった「微分代数方程式」も扱うことができる．
| NDSolve[{eqn1,eqn2,...},y,{x,xmin,xmax}] |
| x がxmin からxmax の範囲で関数y の数値解を求める |
| NDSolve[{eqn1,eqn2,...},{y1,y2,...},{x,xmin,xmax}] |
| 複数の関数yi の数値解を求める |
微分方程式の数値解法
NDSolveで求まる関数
yiの解は，
InterpolatingFunctionオブジェクトを使って表される．
InterpolatingFunctionオブジェクトは，独立変数
x が
xmin から
xmax で指定される範囲にあるとき
yi の近似を与える．
NDSolveでは解を求めるために反復手法が使われる．それは
x のある値を始点とし，一連の刻みを取りながら
xmin から
xmax の範囲をくまなく探索していく．
探索を開始するには，
NDSolveの式の中に
yiとその導関数に関する適切な初期条件や境界条件が備わっている必要がある．これらの条件を与えることで，
x の特定点における
yi[x]，さらに導関数
yi'[x]の値を決定していく．少なくとも常微分方程式において，この条件は
x のどの点でもよい．それは
NDSolveが
xmin から
xmax の全範囲をくまなく自動的に探すからである．
yについて解く． y[0]の初期条件を与え， xは 0から 2の範囲で探索する．
| Out[1]= |  |
|
xの探索範囲は 0から 2のままだが，初期値は y[3]について与えている．
| Out[2]= |  |
|
| Out[3]= |  |
|
NDSolveを使う際は，
yiの解を完全に決定するのに十分な初期条件，または，境界条件を与えておく必要がある．
DSolveを使って微分方程式の厳密解を求めるときには，初期条件が足りなくてもどうにかなる．これは，指定しなかった初期条件のため発生する自由度を明示化するため，積分定数
C[i]が自動的に挿入されるからである．ところが，
NDSolveでは解は数値として与えなければならないので，数値では余分に発生する自由度が表せない．このため，解を決定するのに十分な初期条件，または，境界条件を明示的に与えなければならない．
n 階の微分方程式を解くとき，通常，最大で
(n-1)次までの導関数について初期条件を与えるか，または，
n 点における境界条件を与える必要がある．
3階の微分方程式を解く．3階なので，最大で2次までの導関数を与えなければいけない．
| Out[4]= |  |
|
| Out[5]= |  |
|
別の3階微分方程式を解く．境界条件を3点について与える．
| Out[6]= |  |
|
Mathematica では，線形結合した関数値や微分値なら何でも境界条件として使える．
| Out[7]= |  |
|
多くの場合，与える初期条件は必ず同一の
x 値，例えば
x0 を含む必要がある．このため，
xminから
xmaxの範囲は別に与えなくてもよくなる．
x の範囲を
{x, x1}として指定すると，
Mathematica では自動的に
x0から
x1の範囲で有効な解が生成される．
| Out[8]= |  |
|
初期条件とする方程式はどんな種類のものでもよい．場合によっては，これらの方程式は複数の解を持つことがある．そのような場合は，それに対応した形で
NDSolveが複数の解を生成する．
| Out[9]= |  |
|
| Out[10]= |  |
|
NDSolveを使って，結合された微分方程式の系を解くことも可能である．
| Out[11]= |  |
|
| Out[12]= |  |
|
xと yの両方を使いパラメトリックプロットを生成する．
| Out[13]= |  |
|
微分方程式に未知の関数を使うとき，関数の名前として必ずしも単一のシンボルを与える必要はない．未知の関数が複数個あるとき等は，例えば
y[i]といった関数名にするとよい．
| Out[14]= |  |
|
| Out[15]= |  |
|
| Out[16]= |  |
|
NDSolveは値がリストや配列である関数も扱える．
y[0]
{v1, v2, ..., vn}のような初期条件を与えると，
NDSolveは
y を値が長さ
n のリストである関数だと仮定する．
これは4つの微分方程式が組み合された系の解を求める．
| Out[17]= |  |
|
| Out[18]= |  |
|
NDSolveの特別オプション
NDSolveには方程式を解く多くのメソッドが備わっているが，基本的にはすべてのメソッドが独立変数
x 内で一連のステップを取ることで解を求める．この手順は適応型であり，刻み幅は常に適当であるように自動調整される．値が特定の領域で急激に変化するときは，
NDSolveが解を見逃さないように刻み幅が徐々に狭められるか，あるいはメソッドが変更される．
| Out[19]= |  |
|
NDSolveの解の探索において， x=0の近傍では刻み幅が狭められ，折れ部でも正確に追跡ができている．
| Out[20]= |  |
|
NDSolveの探索手順は適応型であるから，変数
x の変化率が大きく異なるいくつかの要素を持つ，いわゆる「扱いにくい」微分方程式を解くこともできる．
これらの方程式において， yは zより急激に変化する．
| Out[21]= |  |
|
| Out[22]= |  |
|
NDSolveで使われる手順では，解を確実に追跡できるまで刻み幅を徐々に狭めている．ただし，この手順では，真の解に特異点が存在すると問題が起る．つまり，
NDSolveで刻み幅の細分化が限りなく行われてしまい，探索が終了しなくなってしまう．この問題を避けるため，
NDSolveに許す最大ステップ数をオプション
MaxStepsで指定しておく．常微分方程式の場合，このオプションはデフォルトで
MaxSteps->10000になっている．
| Out[23]= |  |
|
| Out[24]= |  |
|
解が滑らかなほとんどの関数では，
MaxStepsのデフォルト値で十分間に合う．解が複雑な構造を取るようなときは，
MaxStepsを大きく取っておくとよいだろう．また，
MaxSteps->Infinityと設定しておくとステップの上限なしで探索ができる．
ローレンツ(Lorenz)方程式の解をここまで求めるためには， MaxStepsのデフォルトによる制限を除いておく必要がある．
| Out[25]= |  |
|
求まった解を3Dパラメトリック形式でプロットする．
| Out[26]= |  |
|
連立した微分方程式を解くとき，
NDSolveではすべての方程式に適応するように刻み幅が自動調整される．しかし，場合によっては，
NDSolveで最初に取るステップ幅が大きすぎて解の重要な特徴を見逃しかねない．このような問題を避けるには，オプション
StartingStepSizeを使い第1ステップの幅を特別に指定しておく．
NDSolveに与える方程式すべてが導関数を含んでいる必要はない．単に代数的なだけでもよい．
NDSolveを用いて多くの「微分代数方程式」を解くことができる．
| Out[27]= |  |
|
| Out[28]= |  |
|
| NDSolve[{eqn1,eqn2,...},u,{t,tmin,tmax},{x,xmin,xmax},...] |
| 偏微分方程式系をu について解く |
| NDSolve[{eqn1,eqn2,...},{u1,u2,...},{t,tmin,tmax},{x,xmin,xmax},...] |
| 偏微分方程式系を複数の関数ui について解く |
偏微分方程式を数値的に解く
これは波動方程式の数値解を求める．結果は2次元の補間関数である．
| Out[29]= |  |
|
| Out[30]= |  |
|
| Out[31]= |  |
|
| Out[32]= |  |
|
| Out[33]= |  |
|
これはこの方程式の2+1次元におけるバージョンである．
| Out[34]= |  |
|
| Out[35]= |  |
|
| Out[36]= |  |
|
