コンピュータリテラシ 第14回
MATLAB（2b）
コンピュータでの数値計算：やってはいけない計算
今日のテーマ：
コンピュータによる数値情報の扱い
二進数による実数の表現
を見てみよう！
２進の浮動小数点による数の表現
>> f=[1 1];
>> f .* 0.5 .^(1:length(f))
>> sum(f .* 0.5 .^(1:length(f)))
0.11(2)= (1/2)+(1/2)^2
>> f=[0 1 0 1 0 1];
>> f .* 0.5 .^(1:length(f))
>> sum(f.*0.5.^(1:length(f)))
0.010101(2)
コンピュータ内部での実数は、２進の指数表現。
1.1(2) × 2-1
仮数部が、「１．-----」（小数点より左が必ず１になるよ
うな形）に正規化した後に、数の大きさを指数部で表
している。
1.0101(2) × 2-2
IEEE754規格の倍精度（64bits）実数
1bit
０
11bits
52bits
００・・・・・・０
００････････････・・・・・・・・・・・・・・・・・０
仮数部 ： 2進52桁は、10進数の15桁くらいの有効数字
指数部（符号含む）
仮数部の符号
符号は、１のとき負の数
では、いよいよ、
やってはいけない計算、について
自分で確かめる
浮動小数点による数の表現の限界による誤差
オーバーフロー、アンダーフロー
>> realmax
最大の数
>> realmin
最小の数
realmaxやrealminの周辺で起きる
>> (1e150)^2
>> (1e200)^2
無限大
>> 1/Inf
ゼロ
オーバフローとアンダーフローを観察
>> (1e-150)^2
>> (1e-200)^2
ゼロ
>> 1/(1e-200)^2
無限大
実数の丸め誤差による計算結果の狂い
>> x=1；
<= x=1からはじめ
て、 1+xが１と等し
くない間xを半分に
していく。
>> bitex(0.1)
<= 関数bitexは仮数部（有効
桁部分、２進の小数５２桁）の2
進表現を求める。
>> while (1+x ~= 1)
ans =
x=x/2;
1001100110011001100110011001100
110011001100110011010
end
<= 0.1は、無限小数なので、有効桁数で丸め
られることが観察できる。
>> x
x=
1.1102e-016
<= While文が停止した時のｘの値（数学的にはどんなにxを1/2しても1+x
とxは等しくならないはずだが） このwhile文は誤差のおかげで止まる
0.1の二進表現に係わる
計算の狂い
>> x=0.1; disp(x);disp(bibtex(x));
>> while x<1
x = x + 0.1; disp(x);
disp(bitex(x));
0.500000000000000
0000000000000000000000000000000000000000000000000
000
end
(省略）
0.200000000000000
0.900000000000000
10011001100110011001100110011001100110
01100110011010
110011001100110011001100110011001100110011001100110
0
0.300000000000000
1.000000000000000
00110011001100110011001100110011001100
11001100110100
0.400000000000000
10011001100110011001100110011001100110
01100110011010
1111111111111111111111111111111111111111111111111111
1.100000000000000
00011001100110011001100110011001100110011001100110
01
fl.m
• 関数ｆｌは、数値ｘを10進の有効桁数５桁に丸め
る関数である。
使用例：
>> pi
>> format long
>> pi
>> fl(pi)
←有効桁５桁
• 倍精度実数計算ではだいたい１０進数１６桁の
有効桁数であるが、何桁であっても有限桁数
が決まっていればそれによる誤差が生じること
を関数flを使って観察しよう。
有効桁数に伴う情報落ち
>> x1=121730; x2=3.1481;
>> x1+x2
ans =
1.217331481000000e+005
有効桁数１６桁の正しい答え
>> fl(x1+x2)
ans =
121730
赤： 誤差の大きな結果
有効桁数が５桁しかないと、結果は
こうなる。x2の値は無いのと同じ。
有効桁限界あたりの怪しい計算
例： ２数の平均
>> a=8.1651; b=8.1653;
>> (a+b)/2
ans =
8.16520000000000
>> fl(fl(a+b)/2)
ans =
正しい答え
同じ計算を有
効桁数５桁で
すると、誤差
が大きい
平均値であり得ない値！
8.164999999999999
赤： 誤差の大きな結果
>> fl(a+fl(fl(b-a)/2))
ans =
8.165200000000001
有効桁数５桁でも
工夫によって正し
い答え
桁落ち（ごく近い２数の差で有効桁消失）
例： 2次方程式の解
𝑏 2 − 4𝑎𝑐 = 123.98
>> a=1; b=-124;c=1;
>> (-b+sqrt(b^2 - 4*a*c))/(2*a)
ans =
1.239919349593155e+002
>> (-b-sqrt(b^2 - 4*a*c))/(2*a)
ans = 0.008065040684528
>> fl(fl(124+123.98)/2)
ans =
1つの解は大体正しいが、もうひとつの解
は誤差が大きい
1.239900000000000e+002
>> fl(fl(124-123.98)/2)
赤： 誤差の大きな結果
ans = 0.010000000000000
>> fl(1/123.99)
計算の工夫で正しい解
正しい２つの解
ans = 0.008065200000000
緑： 誤差許容範囲の結果
2次方程式の2つの
解をαとβとすると、
β＝c/α*a
を使
用
桁落ち（ごく近い２数の差）
例： 分子の有理化
>> x=0.0012345;
>> sqrt(1+x)-1
𝒙 = 𝟎. 𝟎𝟎𝟏𝟐𝟑𝟒𝟓のとき
𝟏＋𝒙 − 𝟏 を求める
有効桁５桁でも計算の工夫
正しい答え
𝑥
１＋ｘ ＋１
を求める
ans = 6.170596187133892e-04
>> fl(1+x)
>> fl(1.0006+1)
ans = 1.001200000000000
ans = 2.000600000000000
>> fl(sqrt(1.0012))
>> fl(0.0012345/2.0006)
ans = 1.000600000000000
ans =
>> fl(1.0006-1)
ans =
誤差が大きい
6.000000000000000e-004
6.170600000000000e-004
中間結果のオーバフロー
>> a=2.0e200;
>> s=max(abs([a b]));
>> b=1.0e200;
>> c=s*sqrt((a/s)^2 + (b/s)^2)
>> c=a^2 + b^2
c=
c=
Inf
cは非常に大きな値のため
オーバーフロー
2.236067977499790e+200
オーバーフローを避ける計算の工夫
S=max(|a|,|b|)
𝑎
C= 𝑎2 ＋𝑏 2
𝑏
C = s* ( 𝑠 )2 ＋( 𝑠 )2
やってはいけない計算のまとめ
1+ε=1であるようなε>0が存在する(!!)
• 小さい数を足すと小さい数の情報が失われる
ことを『情報落ち』という
• たとえば
1.0000000000000000+0.0000000000123456789 だと、
本当はまったく1にはならないはずが
1.0000000000123456 のようになって、有効桁
数より後ろの部分の情報は勝手に失われる)
実数で数えてはいけない
たとえばプログラミング言語で、ある数を0からはじめて1になるまで0.1を足すとこんなこと
が起きる。
1回目
0.10000000000000000000
2回目 0.20000000000000000000
3回目 0.30000000000000004000
4回目 0.40000000000000000000
5回目 0.50000000000000000000
6回目 0.60000000000000000000
7回目
0.70000000000000000000
8回目
0.79999999999999990000
9回目 0.89999999999999990000
10回目 0.99999999999999990000 ← １０回足しても１にならないから１１回やる
11回目 1.0999999999999999000
・このように、0.1を足していくと微妙に誤差がある。
・このためちょうど1にならず、11回足してしまう。
・逆に「0.333333333333333333」だと3回で1になる。
・教訓: 決まった回数実行するときは実数で数えてはいけない
実数で数えてはいけない
• 根本的な原因は、0.1という数値が2進法だと
正確に表せず、近似値になること。
• もともと、コンピュータでの実数計算はすべて
近似値だと思え!
• 近似値計算による誤差のことを「丸め誤差」と
呼ぶ
数式として等しくても
計算に向く式と向かない式がある
• 例: 「√(1+x) - 1」を計算したいとする。xは0にか
なり近い。
• 方法1: そのまま計算
「 √(1.0 + x) - 1.0」
• 方法2: 元の式の分子と分母に「√(1+x) + 1」を
掛けて整理
「x / (√(1.0 + x) + 1)」
数式として等しくても計算に向く式と向かない式がある
たとえば、x = 0.0000000000002 だとする。
1+x = 1.000000000000002 だから
√(1+x) ～ 1.000000000000001 のはず(展開して確認しよう)
つまり答え ～ 0.000000000000001 のはず
あるプログラミング言語で実際に計算した結果
方法1： 8.8817841970012520000e-16
方法2:
9.9999999999999970000e-16
・方法１（素朴に計算）→有効桁の値がデタラメになる
・方法２（変形してから計算）→有効数字16桁で合っている
・何が起きているのか? たとえば、16桁の精度で表された非常に近い2数を引き算
1.0000000000123456
-) 1.0000000000123333
--------------------0.0000000000000123
・最初は有効桁数16桁あったはずが→計算後は3桁しかない!
数式として等しくても
計算に向く式と向かない式がある
• 非常に近い2数を引き算すると有効桁数が激
減すること→『桁落ち』と呼ぶ
• 桁落ちを起こさないためには、近い2数の引き
算を避ける(この例題での変形のようにして)
級数による計算
• 数値計算では、級数による計算をすることも多
い。
• 例: π = 4 - (4/3) + (4/5) - (4/7) + (4/9) - ...
無限級数をn項求めることで円周率を求められ
るか、あるプログラミング言語で計算してみた
ところ
級数による計算
public class Test {
public static void main(String[] args) {
int n = Integer.parseInt(args[0]);
double pi = 4.0;
double sign = -4.0;
for(int i = 0; i < n; ++i) { pi += sign / (2.0*i + 3.0); sign = -sign; }
System.out.printf("%.20g¥n", pi);
}
}
> java Test 10
： n=10項まで計算
3.2323158094055940000
> java Test 100
： n=100項まで計算
3.1514934010709914000
> java Test 1000
： n=1000項まで計算
3.1425916543395442000
級数による計算
• 級数を計算するとき、無限にはできないので適
当なところで打ち切る
→打ち切った残りは誤差になる
→『打ち切り誤差』
• 打ち切り誤差を小さくするには
→計算する項を増やす(遅くなる)、
→収束の速い級数を使用する(正しい方向、
ただし常に可能とは限らない)
以上のような理由から、
コンピュータによる計算結果が真の値から
狂ってしまうことがあることを知っておく。
MATLABのライセンス
• 東工大生のMATLABの利用：
TAH(Total Academic Headcount)ライセンス
http://wwwnode/39.gsic.titech.ac.jp/
© Copyright 2024 ExpyDoc
