GP2X Fixed Point Math

[sparr]

I am starting to make some progress on a decent fixed point library for C++. I know others exist, like the C fixed type in Allegro, but I plan to implement more and better mathematical functions than most of the simple alternatives. I want my type to be as close as possible to a drop-in replacement for float/double in an average game.

The negatives of fixed point math is that you have to decide on your range and precision ahead of time, which limits your options slightly when using a fixed point variable. For trig you may want a 0.31 type, giving a range of -1 to +1 with a precision of about .0000000005. For inverse trig, or other random handling of small numbers, 8.24 is a great type, -127 to +127 is the perfect range for binary angles, and 24 bits of decimal is one bit MORE than you get with a "0.23" float. 16.16 is your general purpose type, and the one I have implemented for benchmarking right now, which provides a range of 65536 and a precision of about .000015. 24.8 gives you very little decimal precision, about .004, but is handy for the few cases where you need 'big' numbes but also a couple of digits past the decimal place.

Eventually I would like a version of the library that defines multiple distinct types of different precisions and provides for correct promotion between them. That is a little ways off, definitely past having most of the functionality of libm provided.

Of course the positive of fixed point math is speed. Since we don't have a FPU we need to avoid floating point numbers like the plague. I wrote a little mandelbrot set generator for benchmarking my fixed type and the results are as follows for a 60x40 pixel image with 10000 iterations per pixel and an abort range of 4:

[root@gp2x sd]$./mandel-double
44008190
[root@gp2x sd]$./mandel-float
35858424
[root@gp2x sd]$./mandel-fixed
9254623

Translated... 44 seconds to draw the image using doubles, 36 seconds using floats, and 9 seconds using fixeds. With simple monochrome pixels the actual output of the fixed program differs from the float output on 8 pixels. That is a bit worrisome, until you consider that the float also differs from the double on 2 pixels. Chalk it up to rounding error. The fixed and float outputs are identical if I use 8.24 fixeds (compared to approx 2.21 floats), since the math in the mandelbrot algorithm doesn't require a large range, but I thought it would be fair to use 16.16 since that will be the norm.

I look forward to getting the library into a state where it can be used while porting existing games. Eventually I plan to revisit my Enigma port and see about moving the physics over to fixeds to resolve the FPS issue on physics-heavy levels like the meditation levels.

 

[dockthepod]

Are you planning on sharing this library when you get it finished? Sounds quite handy. I've been looking around for a nice C++ fixed point library.

 

[Puck2099]

Are you planning on sharing this library when you get it finished? Sounds quite handy. I've been looking around for a nice C++ fixed point library.

Have you tried "math-sll"? I used it on my Fenix port and it increased speed a lot in some cases

 

[rabidcow]

I've got a fixed point template somewhere that lets you pick arbitrary (within limits) fixed point precision. It tracks digits across operations as well (at compile time), so eg. a*b/c doesn't shift down in between.

Eg, you define Fixed<4,12> a; and it takes 2 bytes (16 bits). The expression a*a has type Fixed<8,24> and if you convert back to Fixed<4,12> it shifts down to match.

I've got sin, cos, log2, exp2, gamma, I forget what else. Have to find the code... If not, it's at least neat approach IMHO.

 

[sparr]

rabidcow, we call that kind of implementation "floating point"

Seriously though, it's a neat idea. You can still be faster than IEEE floats with a floating-fixed-point math type, but only barely so.

As to my library, yes it will be public once it is ready.

For simple precision testing I am iterating i=i*1.25f-sqrt(i)*.234375f fourty times for fixed, float, and double, then comparing the results. With a 16.16 fixed and no rounding (during operations requiring a >>, like * and /) float is about 100x as precise as fixed for numbers near 1 (since float functions like 1.22, 6 more bits of precision is about two decimal digits). With 16.16 and rounding, float is 33x as precise as fixed for numbers near 1. When I move up to 8.24 fixeds is when they really shine. As long as you don't need more range than 256 you get MORE precision than floats. Without rounding 8.24 fixed is about 2x as precise as float, and with rounding it is about 8x as precise. I can't easily compare fixed to double because double is my 'standard', it is assumed to be perfectly accurate.

 

[BradN]

If you're testing on an x86 machine, you could use the native 80 bit format and compare to that.

 

[rabidcow]

rabidcow at we call that kind of implementation "floating point"

Seriously though, it's a neat idea. You can still be faster than IEEE floats with a floating-fixed-point math type, but only barely so.

Hehe It's only really floating at compile time though. The number of fractional bits is encoded in the data type, at run time it looks like fixed point. Theoretically it should be marginally faster than a lib based on an actual fixed number of fractional bits, because it doesn't insist on always adjusting results. I haven't done any speed testing on it though...

Fixed<4,12>(1) is a single unsigned short with value 0x1000. It's all magic template sugar.

 

[sparr]

the shift amounts for operations have to be calculated at run time though. thats where your big hit will come in.

 

[BradN]

If it inlines every operation as it's used, it shouldn't have to calculate conversion aspects at runtime.

 

[andrew_j_w]

[root@gp2x sd]$./mandel-double
44008190
[root@gp2x sd]$./mandel-float
35858424
[root@gp2x sd]$./mandel-fixed
9254623

Are you comparing against soft-floats and not hardware floats using a version of GCC later than 4.0? I've recently been writing my own C based fixed point library and I'm not seeing anywhere near as big a speed up as you do.

Andrew

 

[sparr]

I am using oopo's toolchain. gcc 4.0 with soft floats. (note: static linking against libm and libgcc makes for a huge binary).

what kind of speedup are you seeing? from a purely mathematical standpoint, 4x seems about right. multiplying two fixeds is 3 instructions, multiplying two floats is 7-23 (by my ill-remembered count of __mulf* from libgcc.a) instructions depending on branches based on differences in scale.

 

[andrew_j_w]

I am using oopo's toolchain. gcc 4.0 with soft floats. (note: static linking against libm and libgcc makes for a huge binary).

what kind of speedup are you seeing? from a purely mathematical standpoint, 4x seems about right. multiplying two fixeds is 3 instructions, multiplying two floats is 7-23 (by my ill-remembered count of __mulf* from libgcc.a) instructions depending on branches based on differences in scale.

Basically I'm not seeing much of a speed up at all. The difference in speed between by 16.16 fixed point type and float for addition and multiplication is virtually nonexistant.

I'm wondering if it's actually my timing code that is the problem, I'm sure I'll work it out

Andrew

 

[lemon]

rabidcow, I find your idea of not shifting during concatenated operations very interesting. Anyway, I'm finding a problem implementing it. My Fixed<> class is a template where only the fractional part is specified (it uses always a 32bit integer to strore the values so a Fixed<16> would be 16.16)
For the conversion I should define then something like:
template<int T1> template<int T2> Fixed<T1>& operator=(const Fixed<T2>& f)

but... it must do something different (a richt or a left shift) depending on if T1 is greater or lower than T2. How do you specify this kind of behaviour?
I know there are more problems like the necessity to use a 64bit integer for multiplication results, but first I should find a solution to this problem.

 

[screetch]

template<int shift, bool IsSuperiorZero>
unsigned int shiftWhateverTheSignIs(unsigned int);

template<int shift>
unsigned int shiftWhateverTheSignIs<shift, true>(unsigned int val)
{
	return val << shift;
}

template<int shift>
unsigned int shiftWhateverTheSignIs<shift,false>(unsigned int);
{
	return val >> -shift;
}

and call shiftWhateverTheSignIs< T1-T2,(T1>=T2) >(val)
(something like this should work, i'm not a GCC compiler though

 

[rabidcow]

rabidcow, I find your idea of not shifting during concatenated operations very interesting. Anyway, I'm finding a problem implementing it. My Fixed<> class is a template where only the fractional part is specified (it uses always a 32bit integer to strore the values so a Fixed<16> would be 16.16)
For the conversion I should define then something like:
template<int T1> template<int T2> Fixed<T1>& operator=(const Fixed<T2>& f)

but... it must do something different (a richt or a left shift) depending on if T1 is greater or lower than T2. How do you specify this kind of behaviour?
I know there are more problems like the necessity to use a 64bit integer for multiplication results, but first I should find a solution to this problem.

Pretty much as screetch says, except I use a more generic If template:

  template <bool CONDITION, typename IFTRUE, typename IFFALSE>
  struct If						 { typedef IFTRUE ValueT; };

  template <typename IFTRUE, typename IFFALSE>
  struct If<false, IFTRUE, IFFALSE> { typedef IFFALSE ValueT; };

  struct ShiftRightF
  {
	template <typename ValueT>
	ValueT operator()(ValueT value, int shift) const
	{ return value >> shift; }
  };

  struct ShiftLeftNegF
  {
	template <typename ValueT>
	ValueT operator()(ValueT value, int shift) const
	{ return value << -shift; }
  };

	template <int IRHS, unsigned FRHS>
	static BaseT convertFrom(Fixed<IRHS, FRHS> rhs)
	{
	  typedef Fixed<IRHS, FRHS> FromT;
	  typedef typename FromT::BaseT FromBaseT;

	  /* use the larger base type */ 
	  typedef typename
		  If< (FromT::nTotalBits > nTotalBits),
			  FromBaseT,
			  BaseT
			>::ValueT
		ValueT;

	  ValueT value = rhs.data_;

	  /* GCC whines if this is an acual if */ 
	  typename
		  If< (FromT::nIntegerShift > nIntegerShift),
			  ShiftRightF,
			  ShiftLeftNegF
			>::ValueT
		shiftF;
	  value = shiftF(value, FromT::nIntegerShift - nIntegerShift);

	  /* need the cast to eliminate a warning if we used FromBaseT */ 
	  return static_cast<BaseT>(value);
	}

I've been meaning to put up the full code, but I've had hosting issues and I haven't had the time to complain about it properly...