実数演算の落とし穴。
とあるblogを見て、ある問題の原因がすぐわからなくて調べたことから。
C言語のmath.hに付属するfmod関数があります。
これは小数点以下の余りを算出するときに使えます。
しかし常に期待した結果とは限りません。
#include <iostream> #include <math.h> using namespace std; int main(){ double x = 5.3, y = 0.1; cout << fmod(x,y) << endl; }
上の出力結果は僕の環境では0.1となります。
本来の計算結果ならば余りは出ないので、0を返すべきです。
どうしてこうなるのでしょうか。
これには浮動小数点演算の精度とfmod関数の仕様が絡みます。
浮動小数点演算の仕組み。
wikipedia - 浮動小数点数
http://ja.wikipedia.org/wiki/%E6%B5%AE%E5%8B%95%E5%B0%8F%E6%95%B0%E7%82%B9%E6%95%B0
浮動小数点数の精度の項目が参考になります。
float型は32bit、double型は64bitのbit長を持ちます。
しかし基本的な仕組みは同じで、符号部、指数部、仮数部で構成されています。
符号部がついていたら負の値として、以下の計算式になります。
基数^(指数部-バイアス値)*(1+仮数部)
基数は大体2、バイアス値はfloatが127、doubleが1023とされています。
wikipediaの例にある、float表現での以下のbit長の算出をして見ます。
0 01111100 01000000000000000000000
3つの区間に区切りました。前から符号部、指数部、仮数部です。
符号部の1bitは0なので、負の値でないことがわかります。
指数部の8bitは01111100なので、符号無しで計算すると124です。
仮数部の23bitは0100...0000.で、これは2進小数表現であり、
上位bitから、0.5, 0.25, 0.125, 0.0625, ...となり、bitがある位置の総和となります。
ちなみに、以下の式でfloatでの仮数部の最大値を求める事ができます。
23
Σ1/2^n
n=1
今回は0100...0000であれば、0.25しかないので、仮数部は0.25となります。
以下を式に当てはめると、
2^(124-127)*(1+0.25)
=(2^-3)*1.25
=0.125*1.25
=0.15625
となります。
浮動小数点は2の指数部(-127〜128)乗と仮数部(1.0000...〜1.9999...)の値の掛け算で数値を表現しています。
この表現方法の強みは、整数型では表現できないような大きな値でも、同じbit数で表現できる点です（正確には近似値）。
ここで注意したいのは、仮数部の算出方法によってどうしても表現できない値があります。
floatであれば、仮数部は23bitなので、1/2^24より小さい値は表現できませんし、
計算方法を見てわかるとおり、1/2+1/4+1/8+...としているため、
1/3の値を仮数部では正確に表現することはできません。
これらの事らが誤差を生む事に繋がります。
実際よく出力で見る実数の値は、実際の表現を切り上げた近似値に過ぎないのです。
fmod関数
次にfmod関数は以下の実装がされているようです（実際はエラー処理が入ります）
double fmod(double x, double y){ int n = x/y; return x - n*y; }
ここではxをyで除算しています。
上に上げたとおり、表現に限界があります。今回挙げた例では、5.3/0.1を行いますが、
算出される値は53でなければなりません。
実際にx/yの出力結果をざっと見ると53と出ると思います。
しかし、小数点以下20位あたりまで出力するようにしてみると真実が見えてきます。
http://codepad.org/aELjR2GR
5.3/0.1=52.99999999999999289457という表現になっています。
これをint型にするわけなので、整数部分の52となります。
なのでfmod関数では、5.3-52*0.1 = 0.1となってしまいます。
これはバグではなく、コンピュータ上の仕様と捕らえるしかないでしょう。
実数型はdouble型であっても精度に難がある型です。
若干の誤差が許容できる場面でない限り使ってはいけません。
Jpegの変換にも組み込み系での実効速度を稼ぐために、
質が落ちてでも整数演算で行うという手段をとっているライブラリもあります。
それだけ浮動小数点演算は計算コストがかかるのです。
整数演算にすれば、ビットシフトが容易に使えるので何倍も早くなります。
