🔒 Locked

64 bit 32.32 fixed point divide

Started by Structure Mar 30, 2005 at 11:57 AM 8 replies 6.2k views

Original Post

Mar 30, 2005 at 11:57 AM

Hi, ive been experimenting with making a fixed point 64bit 32.32 data type for use in places where I need the numerical precision but cant sacrifice the range (isn’t physics fun!!! Grr) i dont have a fpu to play with so no why dont you use doubles etc please.. Ive managed to make a working 64bit multiply without overflow like this



CFixed32L operator*(const CFixed32L& other) const {
   		
		
	U32 ah = static_cast<U32>(m_v >> 32);
	U32 al = static_cast<U32>(m_v);
		
	U32 bh = static_cast<U32>(other.m_v >> 32);
	U32 bl = static_cast<U32>(other.m_v);

	I64 result = static_cast<U64>(ah*bh)<<32;
	result += (static_cast<U64>(ah)*bl) 
                   + (static_cast<U64>(bh)*al);
	result += ( static_cast<U64>(al)*bl)>>32;
		
	result -= ((ah >> 31) * (other.m_v<<32)  ) ;
	result -= ((bh >> 31) * (m_v<<32) );

	return CFixed32L::FromRaw(result);	

}

But im a bit stumped how to make a divide, as the top part of the fraction would need to be shifted into a 128bit space, anybody tried this or have any ideas? [edit] all that code can prob be made MUCH faster but this is a proof of concept at the mo, in the long run i think its all going to be too slow to pursue but im interested if its been done cos ive been playing with it for a while now [edit 2] High performance is the name of the game btw :) so fast suggestions please ;-)

JohnBolton

Mar 30, 2005 at 1:19 PM

If your processor does not have a 128-bit divide, then you will have to do the divide in software and that will be slow no matter what you do. However, if speed is more important than precision, you might tweak the precision of the values enough to use a 64-bit divide (which I assume your processor has).

Structure

Mar 30, 2005 at 1:37 PM

@JohnBolton

Thanks for the reply,

yes I can see this is going to be a bit slow, what I was hoping for was some crafty way of separating the divide into smaller chunks so it will fit inside the 64bit data type, a bit like with the multiply above, but I cant think of a way of doing this…

Tjoppen

Mar 30, 2005 at 2:52 PM

Just brainstorming;

__int64 c = a / b,    //integer part        d = a % b;    //fraction * bresult = (c << 32) + (((d << 31) / b) << 1);

Something like that might work..
For instance:

3.0/2.0 = 1.5

a = 3<<32, b = 2<<32
c = 1
d = 1<<32

result = (1 << 32) (=1.0 fixed pt) + (((d << 31) (=1<<63) / (2<<32) (=1<<30)) << 1) (=0.1(base 2) fixed pt)
result == (1<<32) + (1<<31) = 1.1(base 2) = 1.5(base 10)

There are however precision issues. And it probably doesn't work with b > 2.0

Try experimenting with the expression
A/B = (a+b)/(c+d) = a/(c+d) + b/(c+d)

Where a and c are the integer part and b and d are the fractional part of A and B
Possibly one could take the individual bits(ones only needed. 0/X = 0) of A, divide them by B and add them back together to form the result. I dunno, experiment.

JohnBolton

Mar 30, 2005 at 3:52 PM

Quote:
Original post by JohnBolton
If your processor does not have a 128-bit divide, then you will have to do the divide in software and that will be slow no matter what you do. However, if speed is more important than precision, you might tweak the precision of the values enough to use a 64-bit divide (which I assume your processor has).

Further elaboration: with a 32.32 number, you need a 96 bit dividend for a full precision operation. If all you have is a 64-bit divide, then you need to drop 32 bits. You can drop the high-order 0 bits without loss of precision, but the rest would have to be dropped from the low end and that will result in loss of precision.

Now, it just occurred to me that a/b is a * 1.0/b and 1.0/b can be calculated with the loss of only 2 bits of precision by calculating 0x4000000000000000 / b and shifting the result appropriately (right 30 bits). Right? Oh, except if the divisor can only be 32 bits, I guess.

Eelco

Mar 30, 2005 at 4:11 PM

doesnt multiplying also require 128 bit for the result? you could argue you might never go that high, but why would you be using 64 bit for input if you dont expect to use it?

Structure

Mar 31, 2005 at 12:04 PM

Thanks for all the input,

On further investigation I found this

http://webster.cs.ucr.edu/AoA/Windows/HTML/AdvancedArithmetica2.html#1007619

Which basically states you cant do it quickly in full precision, which is a bummer, but to break in to song “that’s life, that’s what all the people say…”

[Edited by - Structure on March 31, 2005 2:04:21 PM]

Nice Coder

Apr 1, 2005 at 1:16 AM

Multiplying two numbers together is pretty easy.

Heres my code for it. (this is for integers. But it should be easy enough to convert to fixed point).

//For x * y = zMinmax(*x, *y); //This makes x the smaller number, and y the bigger number.while(x) {if (x & 1) {z += y}x =>> 1;y <<= 1;}

This is ok. 32 iterations max. (o(log2 n) time).

But there are much better algorithms avalible.

From,
Nice coder

Nice Coder

Apr 1, 2005 at 1:23 AM

For division:

//For x / y = z//Stage 1Int Max = 1;Int T = y;While (T < x) {T <<= 1;Max <<= 1;};//Stage 2Int Min = max >> 1;int g;int k;while(1) {g = (min + max) >> 1;if (g = max) {break;}k = g * x;if (k > x) {max = g;}else {min = g;}}

This is pretty fast.
Not that optimised tho.

From,
Nice coder

Structure

Apr 1, 2005 at 6:03 AM

@Nice coder

Thanks for your examples (btw I always enjoy reading your posts, gets the brain working so to speak)

I think the problem is that your multiply will shift the significant part of the fixed point number out of a 64 bit data type so would require extra work to store the higher part, also working on things bit by bit is not going to give me enough speed in time critical physics code

Ive not worked out the divide yet but I assume it will have the same problem, though I guess you could modify the function so it dosnt need shifting first, but again I think this approach would be too slow for my needs

Topic Locked

This topic has been locked by a moderator. New replies are not allowed.