数値計算の基本について

数値計算を行う上での基本的な注意事項等についてまとめていきます。

桁落ち、情報落ち、丸め誤差、打ち切り誤差、Kahanの誤差補償アルゴリズムについて

数値計算の安定性(陽解法Euler、Runge-Kutta)

桁落ち、情報落ち、丸め誤差、打ち切り誤差、Kahanの誤差補償アルゴリズムについて

桁落ち▼

極めて近い数値同士の減算を行ったり、加算の数値が極めて0に近くなる際に発生。有効桁数が小さくなってしまいます。 

例)準備中

情報落ち▼

極端に桁数の異なる数値同士の加算、減算を行う際に発生。有効桁数によって桁数の小さい数字が計算に反映されなくなります。 

準備中 

丸め誤差▼

有効桁数の大きい数値を、四捨五入や切り捨て、切り上げ等で、有効桁数に丸める際に発生。 

例)準備中

打ち切り誤差▼

循環小数や無理数等、無限に続く数値(2進数において無限に続く数値)を有限桁数で打ち切る際に発生。 

例)準備中

Kahanの誤差補償▼

陽解法オイラーの丸め誤差の影響をKahanのアルゴリズムで抑えることを考えます。 

計算例)詳細準備中

Kahan

数値計算の安定性(陽解法Euler、Runge-Kutta)

陽解法Eulerの安定限界▼

数値計算を行う上で計算刻みは重要です。精度を求めるためには計算刻みを細かくしたいところですが、刻みすぎると計算コストが高くなってしまいます。

一方で計算速度を上げようと、刻みを粗くしていけば採用する計算手法によっては発散してしまう可能性があります。

ちょうど良いバランスを探す必要がありますが、いくつか目安となる条件があります。以下では陽解法Euler、4次のRunge-Kutta法の安定条件についてまとめます。

陽解法Eulerでの計算刻み(Δt)の限界は、系の最小時定数をτminとすると、以下の形で表されます。 

$$2\cdot \tau_{min}>=\Delta t$$

例えば、以下の二階微分方程式を陽解法Eulerで求め、発散するΔtの条件を調べてみます。 

$$\frac{d^2y}{dt^2}+7\frac{dy}{dt}+6y=0$$

ただし、初期条件をy(0)=2,y'(0)=-7とします。 

上式の時定数は、τ=1、1/6となります。 

よって、τmin=1/6となりますので、安定条件はΔt≧2*1/6=1/3となります。 

Δt=1/3前後の時間挙動を見てみます。 

StableCriteria

計算刻みを1/3よりも大きい数値(0.34)にした途端、計算が発散していっていることが確認できます。

 

4次のRunge-Kuttaであれば、このTime Stepの限界値を超えても問題ありません。

 

しかし、実は4次のRunge-Kuttaであっても、安定限界の条件があります。

 

参考書籍は以下 

Philip J Thomas,; ”Simulation of Industrial Processes for Control Engineers”; Butterworth-Heinemann, 1999

Ionut Danaila,Pascal Joly,Marie Postel,;”An Introduction to Scientific Computing Twelve Computational Projects Solved with MATLAB”; Springer, 2007

・また、右のリンクにもう少し詳細な計算を載せています。

当サイトに不具合、ご意見等ございましたらCEsolutionにお知らせください。

Page Top