2018年1月2日火曜日

[Python] 数値計算における誤差

桁落ち
桁落ちとは,値のほぼ等しい数値同士の減算で,有効数字が失われることです.桁落ちを避けるためには,分子の有理化を行います.桁落ちによる誤差は発見が難しいこともあるので,値のほぼ等しい数値同士の減算は避けるべきです.
通常の計算と分子の有理化を行った場合の計算結果例がわかるプログラム例(error1.py)を以下に示します.

# Import module
import math

# Main execution unit
# Case 1: x = 1e15
print("Case 1: x = 1e15")
x = 1e15
res1 = math.sqrt(x + 1) - math.sqrt(x)          # Normal calculation method
res2 = 1 / (math.sqrt(x + 1) + math.sqrt(x))    # Calculation method rationalizing molecules

# Output
print("Normal                 : ", res1)
print("Rationalizing molecules: ", res2)
print()

# Case 2: x = 1e16
print("Case 1: x = 1e16")
x = 1e16
res1 = math.sqrt(x + 1) - math.sqrt(x)          # Normal calculation method
res2 = 1 / (math.sqrt(x + 1) + math.sqrt(x))    # Calculation method rationalizing molecules

# Output
print("Normal                 : ", res1)
print("Rationalizing molecules: ", res2)
print()

実行結果:
$ python error1.py
Case 1: x = 1e15
Normal                 :  1.862645149230957e-08
Rationalizing molecules:  1.5811388300841893e-08

Case 1: x = 1e16
Normal                 :  0.0
Rationalizing molecules:  5e-09


上記のプログラムにおいて,"1e15"は 10^(15) を,"1e16"は 10^(16) を意味します.

丸め誤差
丸め誤差は,実数を有限の桁数の2進数で表現するために生じる誤差です.一般に用いる10進数,無理数,循環小数などを計算機で扱う場合には必ず丸め誤差が生じます.また,10進数で有限小数の場合でも,2進数では循環小数となるような数値を扱う際には丸め誤差が生じることがあります.

丸め誤差が生じる場合の計算結果例がわかるプログラム例(error2.py)を以下に示します.

# Import module
# Decimal value of 0.1
print(0.1)

# Add 0.1 to 1000000 times
x = 0.0
for i in range (1000000):
    x = x + 0.1             # 0.1 is a circulating fraction in binary number

# Output
print(x)

実行結果:
$ python error2.py
0.1
100000.00000133288

実行結果では,0.1を100万回加えると10万にならなければならないところ,0.00000133288の丸め誤差が生じています.
丸め誤差は,有効桁数の2進法を数値表現として用いるコンピュータで実数を表現する限りは,不可避の誤差です.従って,計算過程において丸め誤差が大きく影響を及ぼすようなアルゴリズムは採用するべきではありません.

Pythonには数値計算で問題となる誤差を適切に扱うためのモジュールが用意されています.
decimalモジュールは,2進浮動小数点の誤差を正確に管理するための,10進演算モジュールです.decimalモジュールを利用した計算では,数値は以下のように記述します.

# Import module
from decimal import *

# Main execution unit
# The decimal 0.1 value
print(Decimal("0.1"))

# Add 1 million times Decimal ("0.1")
x = Decimal("0.0")
for i in range(1000000):
    x = x + Decimal("0.1")      # Decimal ("0.1") is different from 0.1

# Output
print(x)

上記の Decimal("0.1") は10進法の0.1であり,error2.py で扱った2進不動小数点表現の0.1とは別の数値です.error2.py をdecimalモジュールを使って書き直したプログラム例を以下に示します.

実行結果:
$ python decimalex.py
0.1
100000.0

浮動小数点数の丸め誤差を低減させるには,約分や通分といった分数計算をそのまま実行するfractionsというモジュールを利用することもできます.
fractionsモジュールを利用すると分数計算において,少数表現による丸め誤差を生じさせずに,分数として計算を進めることができます.

以下に,fractionモジュールを利用したプログラム例を示します.

# Import module
from fractions import Fraction

# Main execution unit
# Fraction calculation
print(Fraction(5, 10), Fraction(3, 15))
print(Fraction(1, 3) + Fraction(1, 7))
print(Fraction(5, 3) + Fraction(6, 7) * Fraction(3, 2))

fractionsモジュールを用いると,分数は以下のように表されます.
    1/5   →  Fraction(510)
    3/15 →  Fraction(315))

実行結果:
$ python fracex.py
1/2 1/5
10/21
62/21

プログラム中では,約分,通分,および分数の乗算などの分数計算が実行され,数学的に正しい値が得られています.

情報落ち
情報落ちは,絶対値の大きく異なる数値同士の演算において,絶対値の小さな数値が演算結果に反映されない現象です.

情報落ちが生じる場合の計算結果例がわかるプログラム例(error3.py)を以下に示します.

# Initial setting
x = 1e10
y = 1e-8
temp = 0.0

# Add y (1e-8) to x (1e10) 10000000 times
for i in range(10000000):
    x = x + y

# Output
print(x)

# First add y (1e-8) 10000000 times
for i in range(10000000):
    temp += y

# Add the result to x (1e10)
x = 1e10
x += temp

# Output
print(x)

実行結果:
$ python error3.py
10000000000.0
10000000000.1

上記の実行例から分かるように,絶対値の小さな値の数値を絶対値の大きな数値に先に加えると,加算の結果が反映されていません.一方で,はじめに絶対値の小さな値の数値同士を加えてから,後で絶対値の大きな数値に加えると,加算結果が最終結果に反映されています.このように,計算順序を変更する(絶対値の小さな数値から順に加えるなど)ことで,情報落ちを防ぐことができる場合があります.

0 件のコメント :

コメントを投稿