硬い方程式

数学・数値解析において硬い方程式は、常微分方程式の数値解法・偏微分方程式の数値解法において、刻み幅を極めて小さくしない限り、数値的不安定になる微分方程式である。硬さを的確に定義するのが困難であると判明したが、方程式に解の急激な変化を起こせる項が含まれていることは確かである。

1. 導入の例
-------

下記の初期値問題を考える。

この問題は直接に解くことができ、厳密解  が次の公式で与えられる。

公式によって、 も明らかである。

同じ振舞いを持つ数値解を求めよう。様々な数値的方法を用いて得られる数値解は右側の画像に表示される。

よって、オイラー法は上記の硬い方程式に対し数値的不安定である。一方、台形公式は数値的安定である。

他の例として、もっとも有名な硬い方程式の一つは、Robertsonの化学反応を支配する方程式系である。

[0, 40]}} のような短い区間では、上記の方程式系を数値的に積分することに問題はない。しかし区間が極めて大きい場合、多数のコードは方程式系を正しく積分することができなくなる。

2. 常微分方程式
---------

上述の例の示すように、硬い常微分方程式の近似解を計算するとき、数値的に安定な方法を使うべきである。常微分方程式における数値的安定性に複数の定義が存在する。特に、線型方程式に対する安定性と非線型方程式に対する安定性を分けて考える必要がある。

2.1. 硬さの比例
----------

線形常微分方程式系の硬さは簡単に測ることができる。一般的な線型方程式系

を考える。上記方程式に対する 硬さの比例  は、行列 A}} の最大固有値を最小固有値で割った商である。つまり、A}} の固有値を //λsub としたとき、方程式に対する硬さの比例を norm と定義する。//

典型的な硬さの比例は、10sup あたりである。極端な場合に、その数は 10sup まで届けるIserles。

非線型方程式の場合は、代わりに関数のヤコビ行列の固有値を使って比例を同じ公式で計算する。

2.2. 線型安定性
----------

線型常微分方程式に対する安定性は 線型安定性 、あるいは 絶対安定性  という。線型テスト方程式

を考える。この方程式は簡単に解くことができ、厳密解は //y// =}} //esup である。Re //λ// < 0}} が成立するとき、y}} の //t// → ∞}} の極限も 0}} である。理想的に、近似解にもそのような振舞いを期待できる。しかし刻み幅 h}} が一定のとき、すべての方法に対する近似解がそのような振舞いを持つわけではない。それを区別するのが線型安定性である。//

一つの方法による時刻 tsub での近似解を ysub とする。複素数平面上の集合

は方法に対する 線型安定性領域 、あるいは 絶対安定性領域  というIserlesBurden。この集合はすなわち、与えられた方法による近似解が期待通りの振舞いを持つすべての //hλ//}} である。特に、ルンゲ＝クッタ法に対する線型安定性領域は以下の形で与えられる。

ここで、//r//}} は等式 //ysub を成立させる関数であり、時々方法に対する 安定性関数 という。例えば、オイラー法に対応する関数は //r// =}} 1 + //z//}} である。//

一般的に、方法に対する安定性領域が大きいほど、その方法はより安定である。よってもっとも安定な方法に対する安定性領域は左複素数平面すべてを含めるべきである。そのような方法を A-安定  というIserles。A-安定な方法は硬い方程式の場合でも、刻み幅 h}} を精度のみの考慮で選択することができ、よって硬い方程式を解くために適切な方法だと考えられる。しかし、優れる安定性を持つ方法を実装するには通常高い計算コストが所要される。そのため、実践では常にA-安定な方法を使うわけではなく、方程式の性質、精度の要件や計算コストの制限などの条件を共に考えてから適切な方法を選ぶのが必要となる。

2.3. 非線型安定性
-----------

上述の安定性理論に考察されたのは線型方程式のみである。その理論は時折り非線型方程式にも適用できるが、決して正しいわけではない。非線形方程式の研究を完全に一般化するのが困難であるように、すべての方程式に対する安定性を考察するのもほぼ不可能である。現在非線形方程式に対する安定性はほとんど単調性条件 を満足する方程式

のみを考える。ただし、 は標準内積である。この発想は、ダールキストによるものであるDahlquist。また、ルンゲ＝クッタ法と線型多段法に対する安定性の定義は異なる。なぜならば、線型多段法は時刻毎に多数の成分からベクトルを記憶する必要があり、偏差を測るには標準内積と異なる内積を定義しなければならないからである。

上記の方程式に対してルンゲ＝クッタ法の安定性は B-安定性  という。方程式にルンゲ＝クッタ法を適用するときに、異なる初期値 //y//sub と hat に対し不等式

が成立すれば、その方法をB-安定と呼ぶ。ここで、//y//sub と hat は時刻 //t//sub でのそれぞれの近似解である。B-安定な方法は必ずA-安定である。

さらに、ルンゲ＝クッタ法の係数が //bsub かつ行列//

が半正定値であるという条件を満足するとき、その方法を 代数的安定  という。代数的安定な方法は必ずB-安定であるHairer。

線型多段法の安定性は G-安定性  という。G-安定性は、B-安定性と同じアイディアを持つが、上述通り標準内積が通用しないので同じように定義することができない。

k}} 段線型多段法の一般形式は次の公式で与えられる。

ここで、αsub と βsub は定数であり、ベクトル //α// =}}  という。 この方法に対応する one-leg法  は次の公式で与えられるButcher。

ただし、

である。線型多段法は対応するone-leg法と同じ線型安定性を持つため、同じ非線形安定性を持つことも期待できる。しかし上記の方程式に対する安定性を分析するには、線型多段法より対応するone-leg法を用いる方が遥かに簡単である。よって以下の定義はone-leg法に対するものである。

与えられた k}} 次正定値対称行列 G}} に対応する 上の内積を以下のように定義できる：

ここで、である。

one-leg法に対し、行列

を半正定値にする G}} が存在するとき、その方法をG-安定というButcher。この定義は初見でルンゲ＝クッタ法の安定性とは全く違うものと思われるかもしれないが、本質的には同じものである。なぜならば、上記の定義を以て、不等式

を証明できるからである。

また、G-安定性もB-安定性のようにA-安定性より強い条件に見えるかもしれないが、実際にA-安定性とは同値であるDahlquist。すなわち、線型多段法がA-安定であることは対応するone-leg法がG-安定であることの必要十分条件となる。

3. 脚注
-----

3.1. 注
------

3.2. 出典
-------

4. 参考文献
-------

first1 = Richard L. .
last1=Butcher .
first1 = Germund .
first1 = Germund .
first1 = Ernst .
title=A First Course in the Numerical Analysis of Differential Equations .
first1 = Arieh .

5. 外部リンク
--------

[1]
[2] Lawrence F. Shampine and Skip Thompson Scholarpedia, 2:2855.doi:10.4249/scholarpedia.2855