山口屋~活動日誌~

私生活で主な出来事をピックアップ

TwoSum TwoProduct Kahan 敗者復活算法 丸め誤差

2022-11-03 20:48:35 | ソフトウェア開発
加算または乗算の丸め誤差を拾うことができるアルゴリズムとして、TwoSum や TwoProduct が知られている。
Donald E. Knuth : The Art of Computer Programming Volume 2: Seminumerical Algorithms, Addison-Wesley (1969)
T. J. Dekker : A Floating-point Technique for Extending the Available Precision, Numerische Mathematik, Vol.18, No.3, pp.224-242 (1971)

エラーフリー変換と呼ばれることもあるが、拾った丸め誤差を上手く使わなければ普通に計算するのと変わらなくなってしまう。なので、エラーフリー変換という表現はあまり好きではない。
TwoSum や TwoProduct の、オーバーフローとアンダーフローへの対策、そこから発展して疑似的な4倍精度演算(倍精度浮動小数点数演算の2倍)を行う手法もあるようだ。
柏木雅英:double-double 演算、double-double 区間演算に関するまとめ(PDF)

級数の和の計算などで、拾った丸め誤差も次の項と一緒に加算する、敗者復活算法とも呼ばれる技法を使うと、精度改善が見られる場合がある。
William Kahan : Further remarks on reducing truncation errors, Communications of the ACM, Vol.8, No.1, pp.40 (January 1965)
S:和>a:加算項という前提で、FastTwoSumで足し込んでいくのである。

浮動小数点数 比較 一致 判定 計算機イプシロン

2022-11-03 20:32:51 | ソフトウェア開発
倍精度浮動小数点数(double型)の値を比較して一致するか判定する際に、 a==b や a!=b では絶対にやってはいけない。計算機イプシロン(C言語では"float.h"のDBL_EPSILONマクロ)を上手く利用し、差が十分小さければ一致と判定する。

ここで、倍精度の浮動小数点数について、改めて振り返っておく。
符号部→1bit。
指数部→11bit。2の-1024~1023乗
仮数部→52bit。0~2の52乗
加算であれば、もとの数に2の-52乗を掛けた数より小さい数は全く反映されない。(情報落ち)

さて、計算機イプシロンであるが、定義が混在している。
→再帰の反復blog:計算機イプシロンのこと
定義1:1.0より大きい最小の浮動小数点数と1.0との差、2^-52(倍精度浮動小数点数の仮数部が52bit)
定義2:「1.0 + eps > 1.0」が真となる最小の浮動小数点数eps、ANSI C
※定義2は丸めモードにより値が異なってくる。
・上への丸め:2^-1074=最小の正の浮動小数点数(倍精度浮動小数点数で、指数部が11bitで-1022乗、仮数部が52bitで-52乗、を表している?)
・下への丸め、0への丸め:2^-52≒2.220446049250313e-16
・最近接丸め(=偶数への丸め):2^-53(≒1.1102230246251565e-16)より大きい最小の数≒1.1102230246251568e-16
実際には
・1.0の大きい方の隣りの数との相対誤差: 2^-52
・1.0の小さい方の隣りの数との相対誤差: 2^-53
のように、
・2^-53 ≦ 隣りの数との相対誤差 ≦ 2^-52
となって、2^-53を計算機イプシロンとして比較の判定に用いた場合は等値演算子で比較するのとほぼ同じになってしまう。