3D Now!命令を使ったlogルーチン

浮動小数の格納形式の所で触れたように小数xが与えられ時, メモリ上にはx=(-1)^s・2^b(a, bは整数, 1<=a<2)と格納されています. すると, y:=log(|x|)=b・log2+log aとなります. 従って1<=a<2の範囲だけのlogを求めることに専念すればよいことになります.

高校の数学で|z|<1の時,

という式を習いました. (細かい事を無視して)この両辺を積分すると

という式が成り立つことになります. ただ, この式をすぐ使うには問題がいくつかあります. 1つは|z|<1という制限, もう1つは収束が遅いということです. そのため, 次のトリックを使います. 上の式でzを-zに変えてみると
log(1+z)=z-z²/2+z³/3-...
となります. 下の式から上の式を引いてみます.

という式が出来ました.

欲しいlogの範囲は1<=a<2です. ところが上の公式では|z|<1という制限がありますからなんとか変換する必要があります. そこで

という変換を考えてみます. すると

となります. つまり,

となるわけです. これは都合がいいです. 何故かと言うと(*)の式は|z|が0に近い時すぐ収束するからです.

(1+z) / (1-z) = a/√2 ですから

さて問題は最後のカッコ内でどこまで和をとるかということです. 少な過ぎては計算誤差が大きくなるし, 多過ぎては無駄な時間がかかります. ここで先程のzの範囲を思い出してみます. |z| < 0.17157...でした. ということは z⁹/9≒1.43*10^-8です. 単精度小数の有効範囲を考えると z⁷/7 までで計算を止めても問題ないと思われます. という訳で z⁷/7 までで打ち切ることに決めました.

結局次の様にするとlog xの近似計算が出来ることが分かりました.

と変形します. aからzを求める部分で2回の足し算と1回の割り算, log xを求めるところで12回程度の足し算 / 掛け算をしています. あまりペアリングに向きそうにない計算ですから1命令あたり1.5clkとして (2+12)*1.5+6(割り算の分)=27clk. 後, aやbを求めたりして大体40clkぐらいを目標にします.

float xが与えられた時 a と b は

a= ( (*(int *)&x) & 0x007FFFFF ) | 0x3F800000;
b=(( (*(int *)&x) & 0x7F800000 ) >> 23 ) - 127;

で求まります. 3D Now!命令でペアリングが可能な命令は

• A群:pfadd, pfsub(r), pfacc, pfcmp?, pf{min|max}, pi2fd, pf2id, pfrcp, pfrsqrt
• B群:pfmul, pfrcpit?, pfrsqit1

としたときA群とB群です.

mov{q|d}命令はどちらの群とも, またmov{q|d}命令同士でもペアリング出来ます. 整数演算命令と違って足し算や引き算同士はペアリング出来ないことに注意します. なるべく掛け算命令の間にA群の命令を挟む様に組んでみました. femmsを除くと大体38clkで動作します. fdivより速い(^^. x=0の時は考慮していません. 一応AMDが提供していたAMDLOG.ASMのlog関数より速かったです.

いかにも3D Now!らしいプログラム, というわけではありませんが, やはり割り算が4〜6clkで出来るのが大きいです. fdivを使うとここだけで40clkぐらいかかりますからね.