以前32ビットのx86で64ビット乗算をするコードを書いていましたが、今回は除算です。除算は乗算と違って割る数に対して分配法則が成り立ちませんから、上位32ビットと下位32ビットに分けることができません。あれこれずっと頭の片隅で考えていましたが、ようやく閃きました。
そうです。x86にも64ビット(以上)の演算装置があるではありませんか。その名はFPU、浮動小数点数演算装置です(以下x87と呼びます)。誤差はどうするのかと思うかもしれません。しかし心配無用です。x87では64ビット整数(符号有り無し両方)を誤差を起こさずに表現できます。
ご存知の方もいらっしゃるでしょうが、x87内部では数値を拡張倍精度と呼ばれる80ビットで扱います。その内訳は符号1ビット、指数15ビット、仮数64ビットです。一般に二進法の浮動小数点数では、(仮数部がよっぽど少なくない限り)指数部のビット数に収まる大きさの整数は誤差を生じずに表現できます(どこでこれを知ったのかは忘れてしまいましたが)。x87拡張倍精度では仮数部が64ビットあるので64ビット整数は正確に表現できると言うわけです。そして、結果も仮数部に収まる限り、正確に表現できる数同士の加減乗除もまた誤差を起こさず結果を表現できます。なお、除算の場合には結果が整数にならないことがありますが、整数除算の場合は小数点以下を切り捨てるので、浮動小数点数演算でも除算後に切捨てを行えばよいのです。これは64ビット整数型へのキャストの際に行われます。
もっとも、残念ながらABではこの80ビットの浮動小数点数を扱うことはできません。代わりにBorland C++ 5.5.1で次のようなコードを書いて、実際に結果が一致することを試してみました。なおlong doubleが80ビットでないといけません(Visual C++はだめです)。
#include <iostream>
typedef __int64 int64_t;
int main()
{
int64_t x, y;
std::cin >> x >> y;
long double z = static_cast<long double>(x) / y; //64ビット整数除算
std::cout << static_cast<int64_t>(z) << 'n';
std::cout << x / y << std::endl; //80ビット浮動小数点数除算
}
ABで書くとこのような雰囲気になります(80ビット拡張倍精度を表す型を仮にLongDoubleとしています)。
#console
Dim x As Int64, y As Int64
Input x, y
Dim z As LongDouble
z = x As LongDouble / z
Print z As Int64
Print x \\ y
End
速度の比較などは行っていません。こういう方法もあるということを思いついたまでです。
スポンサード リンク |
引き戻し法という(除算の)アルゴリズムをご存知でしょうか。
詳しい仕組みはここでは説明しません。というか、僕ができません。
大雑把に言えば、人間が行なう筆算をそのままアルゴリズムとして手順化した感じのものです。このアルゴリズムによって除算を 64bit の減算とシフトに置き換えることができます。
このアルゴリズムを適用すれば、分配法則など数学的なトリックを用いずに直接的に商を算出できます。ちなみに副産物として剰余も出ます。
試行錯誤的にループするので、ハードウェアによる実装(=FPU)にはかなわないでしょうけど・・・。
どうもです。筆算のような方法を使うと言う話は聞いたことがありますが、きちんとした名前があるのですね。今度図書館で探してみます。