固定小数点の平方根

別に平方根自体は固定小数点とは何の関係もないけど固定小数点を使うことにより必要性が出てきた。
平方根自体は標準関数のstd::sqrtを利用すればいいのだがこれは浮動小数点用なので固定小数点は使えない。まあ適度にキャストすれば使えなくもないのだが今は訳あって浮動小数点が使えない。というわけで自前実装が必要になった。
とりあえずC++には標準で固定小数点の型が存在しない。ので適当に__int64を32bitシフトしたのを利用することに決めた。普通の四則演算はきわめて簡単で、足し算引き算はそのまんま。掛け算と割り算が微妙に面倒だが適度にシフトすれば問題なし。というわけで難しい平方根が来ました。
一般的にCPUで平方根をどうやって計算しているのかはしらないが、どうやらニュートン法が主流っぽいんじゃないかと思われる(適当)。しかしニュートン法浮動小数点ならともかく固定小数点だと安定した精度を保つのが非常に面倒そうだから却下(ぇー というわけでいろいろ悩んだ結果、一番慣れ親しんだ開平法を利用することにしました。
http://www004.upp.so-net.ne.jp/s_honma/root.htm
私は何をどう間違った教育を受けたのかこの開平方を小学生の頃から使えました。流石に原理までは理解してなかったけど。しかしおかげさまで今でも余裕で使えます。そういうわけでこいつを2進数で書き換えてプログラムに起こしました。

__int64 FixedSqrt( __int64 num )
{
	if( d < 0 ){ throw; }

	__int64 result = 0;
	__int64 rest = 0;
	__int64 sum = 0;
	for( int i = 0; i < 32 + 16; ++i ){
		rest = ( rest << 2 ) | ( ( num >> 62 ) & 3 );
		num <<= 2;
		result <<= 1;
		sum <<= 1;
		if( ( sum | 1 ) <= rest ){
			sum |= 1;
			rest -= sum;
			++sum;
			result |= 1;
		}
	}
	return result;
}

可読性は激悪ですがちゃんと64bitの固定小数点の精度は出ています。速度は、まあ48回もループをまわすのであんまり速くないと思いますが逆に言えば確実に48回ループで出てくるので速度はえらい安定しています。あと、出来上がってから気がついたのですが掛け算と割り算を一切使ってないんですよね。ていうか殆どがシフト演算という・・・なので案外速いかもしれません。まあ面倒なのでベンチ取りませんが。

まあ目的が浮動小数点を使わずに環境依存しない演算がほしかっただけなのでこれで十分かなと(^^;