● 前回の講義のまとめ
•
一般に,実数x
を浮動小数点数として表すとき,浮動小数点数として表される数x
との誤差は| x − x | ≤ ǫ 2 | x |
となる. （データ誤差）ここで,
ǫ
は浮動小数点数の体系（何進表現で何桁の浮動小数点数か）と「丸 めモード」に依存する定数である. ただし,この不等式は,オーバフローやアンダーフローが生じな い場合に限る.•
また,上記の不等式は| x | ≤ (1 + ǫ/2) | x |
と等価である.•
一般に,浮動小数点数x, y
に対して,四則演算（加算・減算・乗算・除算）を行なった結果の誤差は| (x · y) − (x ⊙ y) | ≤ ǫ | x · y | , | x ⊙ y | ≤ (1 + ǫ) | x · y | ,
となる. （演算誤差）ここで,
x · y
はx, y
が実数だと考え演算を行なった結果（実数）であり,x ⊙ y
はx, y
が浮動小数点数だと考え,浮動小数点演算を行なった結果（浮動小数点数）である. なお,こ の「誤差評価式」が成り立つためには,「減算」に関しては「１桁の保護桁」が必要となる.•
減算における保護桁の存在の仮定の下に,演算誤差がǫ
となる理由は,おおよそ以下の通りと考える ことができる.–
加減算においては,仮数部の桁移動の際の誤差と,演算結果としてくり上がりが発生した場合の 丸めの誤差のそれぞれでǫ/2
の誤差が発生する可能性がある.–
乗除算においては,演算結果の丸めの誤差と演算結果としてくり上がりが発生した場合の丸めの 誤差のそれぞれでǫ/2
の誤差が発生する可能性がある.•
このように浮動小数点演算の誤差は,| x − x | < ǫ
のような「絶対誤差」による評価ではなく,| x − x | < ǫ | x |
のような「相対誤差」による評価が必要となる.•
上記の最大相対誤差ǫ
は,浮動小数点数の体系での「刻み幅」から決まる値と考えて良く,その「丸 めモード」が「四捨五入」の時には,刻み幅の半分の値となる. なお, IEEE 784浮動小数点数規約に おける倍精度浮動小数点数の体系で,丸めモードが四捨五入の時には,ǫ = 2
−53となる.と考えて良い. そこで,
x
n を計算すると,a
n= x
n= (1 + δ)
nx
n= (1 + δ)
na
n∼ (1 + nδ)a
nとなる. したがって,
| a
n− a
n| ≤ nδ | a
n| ∼ nδe
と考えれば,| e − a
n| ≤ | e − a
n| + | a
n− a
n| ≤ e 2n + nδe
が成り立つ. ここで,
n
−2 以上の項とδ
2 以上の項を無視していることに注意. したがって,e
を(1 + 1/n)
n で近似する際の誤差は,e 2 ( 1
n + nǫ)
程度を見積もる必要がある. 実際に, この値を最小とする
n
はn = √
ǫ ∼ 2
−26∼ 10
−8 程度であり, そのときの値はe √
ǫ ∼ 10
−8程度となる.この計算は,すべての演算誤差を無視しているが,仮に
1 + 1/n
の「データ誤差」が0
であったとし ても, (1 + 1/n)2 にはすでに相対誤差ǫ
の誤差が含まれる可能性がある. よって, 以後の演算にもす べて相対誤差ǫ
の誤差が含まれ,演算がn
回程度行われることを考えれば,データ誤差と演算誤差の 総和はおよそnǫ
（の数倍）程度となり,上の評価の本質には影響はない.• (1 + 1/n)
n の計算における演算誤差は,n
回（非常に大きな回数）の演算を行うことによってデータ誤差や演算誤差が蓄積していることが原因である. このように,浮動小数点演算において, 演算回数 が非常に多くなると,丸め誤差（データ誤差と演算誤差）が極めて大きくなる.
• (1 + 1/n)
nを計算する高速な方法は,高速乗算法（およそlog n
回の乗算）またはpow
関数の利用があるが,これらも初期のデータ誤差が同じ程度蓄積するため,これらの方法を用いても,丸め誤差を小 さくすることはできない.
•
実際にe
を計算する効率的な（かつ,誤差の少ない）方法は,e
xのテイラー展開の最初の数項を計算 する方法である. （倍精度浮動小数点数の範囲では, 18項程度計算すればよい）しかし,この場合に も「高次の項からの計算」と「低次の項からの計算」では計算精度に違いが発生する.● 講義資料
▼ 講義予定
•
円周率π
の計算と浮動小数点数演算の誤差▼ 計算結果
• ℓ
n:
半径1
の円に内接する正n
角形の周長の半分• L
n:
半径1
の円に外接する正n
角形の周長の半分★ 方法１
ℓ
2n= 2n s
1 − p
1 − (ℓ
n/n)
22 ,
L
2n= 2n
2L
n− 1 + p
1 + (L
n/n)
2★ 方法２
ℓ
2n= 2n s
1 − p
1 − (ℓ
n/n)
22 ,
2 L
2n= 1
ℓ
n+ 1 L
n★ 方法３
ℓ
2n=
√ 2ℓ
nq 1 + p
1 − (ℓ
n/n)
2,
2 L
2n= 1 ℓ
n+ 1 L
n★ 方法４
2 L
2n= 1 ℓ + 1
L ,
1e-16 1e-14 1e-12 1e-10 1e-08 1e-06 0.0001 0.01 1
10 100 1000 10000 100000 1e+06 1e+07 1e+08
Inner (1) Outer (1) Outer (2)
1e-16 1e-14 1e-12 1e-10 1e-08 1e-06 0.0001 0.01 1
10 100 1000 10000 100000 1e+06 1e+07 1e+08
Inner (4) Outer (4) Inner (3) Outer (4)
● 実習内容
1.
円周率の値を方法１〜４で計算し,π
の真の値と計算値との相対誤差をグラフに表しなさい● レポート問題
以下の問題は評価対象のレポート問題です. 以下の２問のうちいずれか一方を選択してください. いずれ の問題についても,
A: 10
点B: 6
点C: 0
点と採点します. （中間的な点をつける可能性あり）な お, ２問ともレポートを提出した場合には,高い評価（得点）のレポートを評価対象とします.【レポート問題
01-A
】 倍精度浮動小数点演算を用いて,自然対数の底e
の値を, (1 + 1/n)nで近似計算す るプログラムを書き,e
の真の値（倍精度浮動小数点数として最も精度の良い近似の値）と(1 + 1/n)
n との誤差をグラフに表しなさい. ただし,n
として10
9 程度まで計算すること. (1 + 1/n)nの計算方 法としては高速乗算法を用いること. また,実際の浮動小数点計算における| e − (1 + 1/n)
n|
（計算 誤差を含む）を評価する式を求めなさい. ただし,この評価を求める際には, (1 + 1/n)n の計算はn
回の乗算として計算してもかまわない. （要するに第１回目および２回目の講義内容と実習課題をま とめたもの）【レポート問題
01-B
】 倍精度浮動小数点演算を用いて,半径1
の円に内接および外接する正多角形の周 長で近似することにより,円周率π
の近似値を求めるプログラムを書き,π
の真の値（倍精度浮動小 数点数として最も精度の良い近似の値）と, それぞれの値の誤差をグラフに表しなさい. また,単純 に半角公式を用いたとき,それほどよい近似値を得ることができない理由を説明し, その解決方法を 説明しなさい. （要するに３回目の講義内容と実習課題をまとめたもの）★ 注意事項
•
いずれのレポートも以下のものを提出してください.–
計算を実行するためのプログラム（電子メールで提出）–
計算結果を表示した図（EPS形式またはSVG
形式で提出）–
計算結果の考察（PDF 形式またはレポート用紙で提出）【提出方法】 すべての提出物を添付ファイルとして含んだ電子メールを
[email protected]
宛 に送付する. その際の注意事項は以下の通り•
電子メールのSubject
は「Report-xx-y」とすること. （ただし,xx-y
の部分は問題番号と する.）•
添付ファイルのファイル名は,zzzzzz-xx-y.c
などと,拡張子を適切につけ,拡張子よりも前の 部分は, 「学籍番号-問題番号」の形式とすること.•
電子メール本文にも,学籍番号,学年,氏名を記載すること.ただし,「計算結果の考察」については,レポート用紙に書いたものを提出してもよい.
【締め切り】 電子メールで提出する場合には,２００９年５月２７日（水）１４：００までに届くように
