数値計算の基本/numerical calculation basic

数値計算の基本について

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

桁落ち▼

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

例）$\sqrt{1000001}-\sqrt{1000000}$ のように近い値同士をそのまま引くと、有効桁が失われます。こうした場合は有理化して $\frac{1}{\sqrt{1000001}+\sqrt{1000000}}$ の形に変形すると安定して計算できます。

情報落ち▼

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

例）倍精度（有効桁およそ16桁）で $10^{16}+1-10^{16}$ を計算すると、結果が 1 ではなく 0 になることがあります。大きい数と小さい数を同時に扱うときは、加算順序の工夫が重要です。

丸め誤差▼

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

例）$\pi=3.14159265...$ を小数第3位までで 3.142 と丸めると、真値との差 $|3.142-\pi|$ が丸め誤差です。反復計算ではこの誤差が累積する場合があります。

打ち切り誤差▼

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

例）$1/3=0.3333...$ を 0.3333 で打ち切ると、差分 $1/3-0.3333$ が打ち切り誤差です。2進数でも 0.1 は有限桁で表せないため、同種の誤差が生じます。

Kahanの誤差補償▼

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

計算例）例えば 0.1 を100万回足し込むと、通常の単純加算では丸め誤差が蓄積し、理論値 100000 からわずかにずれます。Kahan の補償和を使うと、このずれを大きく抑えられます。

陽解法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前後の時間挙動を見てみます。　

計算刻みを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にお知らせください。