divide on RS/6000

[Bob Alverson]

I've been looking at Markstein's paper "Computation of elementary
functions on the IBM RISC System/6000 processor", and I'm beginning
to wonder if they really use reciprocals for float division.
Can someone try this test case on an RS/6000?

	q = scalb(49, -1074) / scalb(49, -1010);

The correct answer is obviously 2^-64.  However, when I use the
method outlined in the paper, the answer is off by one:

	A = scalb(49, -1074)
	B = scalb(49, -1010)

	Y = 1.0 / B
	Q = A * Y = 2^-64 - 2^-117
	R = A - B * Q = 0 (due to underflow)
	Q' = Q + R * Y = 2^-64 - 2^-117 (off by one!)

Speaking as someone who really wants to do float division using reciprocals,
this is most unpleasant.  If IBM can get it right, so can I.

Ugh.

Bob Alverson (bobatera.com)