8 ms·

16-bit math look-up tables – the unexpected power of scaled-integer math

114 points by mirceasoaica 10y ago

kosma 10y ago
I used this exact approach just yesterday when I needed to calculate some sine values on an 8-bit MCU. Turns out including a precalculated table takes way less program memory than linking softfloat and trig libraries. 8051 is not dead, folks!
Animats 10y ago
Cute. However, unless you're using a CPU from the 6502 era, it's probably not worth the trouble for multiplication and division. Today's low-end CPUs have good multiply hardware, and for a few dollars more, you get a decent FPU. Trig functions, though, may be worth precomputing. The standard libraries for trig functions often grind their way out to far more precision than you need for graphics or control, and that takes time. Interpolating from a modest sized table can be faster.
- kstenerud 10y ago
  Except for the cache hit.
  - Animats 10y ago
    If you're interpolating and only need 16 bits of precision, you need only a small table. 32 or 64 entries should be enough for sine and cosine.
    - raphlinus 10y ago
      An example of the linear-interpolated lookup table approach using fixed-point math is the sine generation[0] in my music synthesizer. It generates a table of 1024 entries, then lerps between those. This gives you precision beyond the least perceptible audio error, and is pretty fast. When I first wrote that code, I was targeting 32 bit ARM, with possibly very weak floating point units. These days, on the types of chips you'd see in phones (or computers in the $9-$40 range like rpi 3), the floating point calculation is _so_ fast (especially with simd) that the cost of the memory lookup starts dominating. So, the NEON implementation computes 4 sin values in a gulp, using a floating point polynomial approximation about the same precision as above[1]. [0] https://github.com/google/music-synthesizer-for-android/blob/master/app/src/main/jni/sin.cc https://github.com/google/music-synthesizer-for-android/blob... [1] https://github.com/google/music-synthesizer-for-android/blob/master/app/src/main/jni/fm_op_kernel.cc#L36 https://github.com/google/music-synthesizer-for-android/blob...
    - jacobolus 10y ago
      It’s probably a generally better idea to use a polynomial approximation if you have a floating point unit; a degree 8 polynomial for sin(x) on the range [0, π/2] gets you to just about the limits of single precision floating point. If you only need about 4 digits of precision, you can use a degree 5 polynomial.
      - adrianratnapala 10y ago
        Couldn't you do better by only using the range [0, π/4] and noting that sin(x) = cos(x - π/2)?
        jacobolus 10y ago
        Sure, then you can get away with an even smaller degree polynomial.
      - sehugg 10y ago
        FWIW, Applesoft for the 6502 uses a 5th order poly.
- weinzierl 10y ago
  > However, unless you're using a CPU from the 6502 era, it's probably not worth the trouble for multiplication and division. When we talk PC, fixed point math was popular a few generations longer than the 6502 era. The 6502 had no multiply and division instructions at all, and up to the 80386 there was only integer multiply and division and that was slow as molasses. Before the 80486 fixed point wasn't a matter of speed, it was a matter of survival (at least for graphics and such).
  - MegaDeKay 10y ago
    I remember my old 8086 PC turbo'd to 8 MHz would crawl along as it rendered the wireframe space shuttle demo image that came with Autocad. I somehow got my hands on an 8087 math coprocessor, plugged that bad boy in, and the difference was phenomenal. Good times.
  - louthy 10y ago
    indeed, I was using fixed point maths on PlayStation 1 games in the mid to late 90s. It was often responsible for the gaps you'd see between polygons on many PS1 games.
    - marktangotango 10y ago
      I was using it for z-80, game boy color homebrew. Getting down in the bits like this is.... bracing!
      - louthy 10y ago
        I remember doing a Gameboy Color game after working on a PS1 title. And even though I came from an 8 bit 6502 background (BBC Micro), going back to it was hell. 'Bracing' is an understatement! I couldn't imagine much worse these days.
    - rasz 10y ago
      Gaps were mostly a sign of sloppy developers and rushed schedules. Compare Tomb Raider 1 on both PC/PSX with Destruction Derby (reflections) and Gran Turismo 1/2 (polyphony).
- jacobolus 10y ago
  It’s usually possible to re-frame problems to not require trig functions at all. For instance, you can represent rotations as unit magnitude complex numbers, compose them using complex multiplication, and trivially get whatever trig functions you want out. If you need to compress them for I/O, take the stereographic projection (requires 1 division per point for both forward and inverse transform) and then optionally reduce the precision of the result.
  - doubleunplussed 10y ago
    Construction of a unit complex number though, given an angle, requires trigonometry. Precomputing this and re-using it is identical to precomputing the sine and cosine of that angle and reusing them instead - the complex number itself doesn't simplify anything here other than storing both the sine and cosine in one variable.
    - jacobolus 10y ago
      The only time you need to start with an angle is if a human or other external system is feeding it to you, and the only time you need to convert to an angle is when you need to present the data to a human or other external system. The point here is that you can usually get rid of evaluating transcendental functions in the middle of the number crunching part of your code, which is where you would care most about saving operations.
      - doubleunplussed 10y ago
        I don't think that's right. A rotation implies an angle. If it's not coming from a human, then you are either hard-coding an angle to rotate by, in which case you need to take sine and cosine of that angle to generate a complex number, or you're getting the rotation that is implied from some coordinates, hard-coded or otherwise. In which case, you need to normalise a complex number - still requiring a square root. I don't think there's a way to do rotations without transcendental functions (or square roots) that isn't equivalent to merely precomputing the transcendental functions or square roots. Maybe square roots are faster, in which case you have a point that avoiding transcendental functions is the way to go.
shpx 10y ago
Related to 16 bit lookup tables, unum is a format for representing numbers as sets of number ranges as an alternative to floats. http://www.johngustafson.net/presentations/Unums2.0.pptx http://www.johngustafson.net/presentations/Unums2.0.pptx http://www.johngustafson.net/presentations/Unums2.0.pdf http://www.johngustafson.net/presentations/Unums2.0.pdf https://arxiv.org/pdf/1701.00722.pdf https://arxiv.org/pdf/1701.00722.pdf https://en.wikipedia.org/wiki/Interval_arithmetic https://en.wikipedia.org/wiki/Interval_arithmetic
- Mindless2112 10y ago
  Unums are unlikely to gain much usage. Posits[1][2], also by Gustafson, are a more reasonable alternative to IEEE-754 floating point (but will still have a difficult time displacing IEEE-754, if they can at all). [1] http://web.stanford.edu/class/ee380/Abstracts/170201-slides.pdf http://web.stanford.edu/class/ee380/Abstracts/170201-slides.... [2] https://www.youtube.com/watch?v=aP0Y1uAA-2Y https://www.youtube.com/watch?v=aP0Y1uAA-2Y
  - adrianratnapala 10y ago
    Posits seem impressive. My only concern is the lack of NaNs seems like a bug rather than a feature. It's true that some programmers do the silliest things when faced with NaNs. But the fact is they are useful. You often want to do calculations over big matrices, where some elements simply don't have a mathematically defined answer (usually because of div0s but also because input data might have holes). It would be a royal pain if the whole calculation stopped every time you had to add two infinities.
    - dnautics 10y ago
      >is the lack of NaNs seems like a bug rather than a feature. I think this is a reasonable concern. I'll propose to John that we make there be an optional "mode" where the infinity token is treated as "NaN". In reality, this mode just amounts to "ignore NaN traps", because the way that it's done in my hardware models, it requires almost no extra hardware.
weinzierl 10y ago
Back in the day of 80286's and 80386's there was a popular Fractal generating software that did all calculations in fixed point arithmetic because PCs usually had no FPU back then. Fractint[1], the site has a certificate error, but otherwise works fine, even the software still gets updated from time to time. [1] https://fractint.org/ https://fractint.org/
- quakeguy 10y ago
  Used it to generate uniform random plasma fields for game textures. Great program!
  - weinzierl 10y ago
    With the right color map they made a nice sky with clouds. Color cycling them through a smooth full color map was also nice.
emersonrsantos 10y ago
Āryabhaṭa's sine table, the first known mathematical LUT, was made around 500 DC: https://en.wikipedia.org/wiki/Āryabhaṭa's_sine_table https://en.wikipedia.org/wiki/Āryabhaṭa's_sine_table
- jacobolus 10y ago
  > first known mathematical LUT Babylonian multiplication and reciprocal tables are about 2500 years older. For that matter, Hipparchus and Ptolemy’s chord tables are also centuries older, https://en.wikipedia.org/wiki/Ptolemy's_table_of_chords https://en.wikipedia.org/wiki/Ptolemy's_table_of_chords
raymondh 10y ago
For those who are interested, CORDIC (for COordinate Rotation DIgital Computer) is another technique from that era. It "is a simple and efficient algorithm to calculate hyperbolic and trigonometric functions, typically converging with one digit (or bit) per iteration." See https://en.wikibooks.org/wiki/Digital_Circuits/CORDIC https://en.wikibooks.org/wiki/Digital_Circuits/CORDIC For instructional purposes, here is a simple Python implementation that uses only adds and shifts in the inner loop, followed by a single scaled-multiply to finish it up: https://code.activestate.com/recipes/576792-polar-to-rectangular-conversions-using-cordic https://code.activestate.com/recipes/576792-polar-to-rectang...
- bsder 10y ago
  IIRC, CORDIC has some issues as you scale the number of bits you want--you can't just run CORDIC 64 times for 64-bits of accuracy, for example. I think the point of diminishing returns is 14-16 bits, but it's been a long time and I'd be happy to be proven wrong.
- planteen 10y ago
  I used CORDIC on a FPGA a few years ago to compute atan2.
marktangotango 10y ago
>> Jack Crenshaw addresses it in chapter 5 of his fine book "Math Toolkit for Real-Time Programming." Ah shit, can't believe I missed this, all these years. Jack Crenshaw ie Let's Build a Compiler. Love that guy.
marktangotango 10y ago
>> Jack Crenshaw addresses it in chapter 5 of his fine book "Math Toolkit for Real-Time Programming." Ah shit, can't believe I missed this, all these years. Jack Crenshaw ie Let's Build a Compiler. Love that guy.
otempomores 10y ago
Is there a stdclib which allows for easy replacement of expensive funtions like fsin or fcos with lutables? Also do those not poison the cache?
rasz 10y ago
>And by the way, real-world I/O in control situations is never floating-point, whether timers, counters, A/D and D/A converters, servos, etc.. well, actually ;) Adlib/Sound Blaster Yamaha OPL2 and OPL3 both used external Floating Point DACs (YM3014B,YAC512)

114 points by mirceasoaica 10y ago

kosma 10y ago

I used this exact approach just yesterday when I needed to calculate some sine values on an 8-bit MCU. Turns out including a precalculated table takes way less program memory than linking softfloat and trig libraries. 8051 is not dead, folks!

Animats 10y ago

Cute. However, unless you're using a CPU from the 6502 era, it's probably not worth the trouble for multiplication and division. Today's low-end CPUs have good multiply hardware, and for a few dollars more, you get a decent FPU. Trig functions, though, may be worth precomputing. The standard libraries for trig functions often grind their way out to far more precision than you need for graphics or control, and that takes time. Interpolating from a modest sized table can be faster.

kstenerud 10y ago
Except for the cache hit.
- Animats 10y ago
  If you're interpolating and only need 16 bits of precision, you need only a small table. 32 or 64 entries should be enough for sine and cosine.
  - raphlinus 10y ago
    An example of the linear-interpolated lookup table approach using fixed-point math is the sine generation[0] in my music synthesizer. It generates a table of 1024 entries, then lerps between those. This gives you precision beyond the least perceptible audio error, and is pretty fast. When I first wrote that code, I was targeting 32 bit ARM, with possibly very weak floating point units. These days, on the types of chips you'd see in phones (or computers in the $9-$40 range like rpi 3), the floating point calculation is _so_ fast (especially with simd) that the cost of the memory lookup starts dominating. So, the NEON implementation computes 4 sin values in a gulp, using a floating point polynomial approximation about the same precision as above[1]. [0] https://github.com/google/music-synthesizer-for-android/blob/master/app/src/main/jni/sin.cc https://github.com/google/music-synthesizer-for-android/blob... [1] https://github.com/google/music-synthesizer-for-android/blob/master/app/src/main/jni/fm_op_kernel.cc#L36 https://github.com/google/music-synthesizer-for-android/blob...
  - jacobolus 10y ago
    It’s probably a generally better idea to use a polynomial approximation if you have a floating point unit; a degree 8 polynomial for sin(x) on the range [0, π/2] gets you to just about the limits of single precision floating point. If you only need about 4 digits of precision, you can use a degree 5 polynomial.
    - adrianratnapala 10y ago
      Couldn't you do better by only using the range [0, π/4] and noting that sin(x) = cos(x - π/2)?
      - jacobolus 10y ago
        Sure, then you can get away with an even smaller degree polynomial.
    - sehugg 10y ago
      FWIW, Applesoft for the 6502 uses a 5th order poly.
weinzierl 10y ago
> However, unless you're using a CPU from the 6502 era, it's probably not worth the trouble for multiplication and division. When we talk PC, fixed point math was popular a few generations longer than the 6502 era. The 6502 had no multiply and division instructions at all, and up to the 80386 there was only integer multiply and division and that was slow as molasses. Before the 80486 fixed point wasn't a matter of speed, it was a matter of survival (at least for graphics and such).
- MegaDeKay 10y ago
  I remember my old 8086 PC turbo'd to 8 MHz would crawl along as it rendered the wireframe space shuttle demo image that came with Autocad. I somehow got my hands on an 8087 math coprocessor, plugged that bad boy in, and the difference was phenomenal. Good times.
- louthy 10y ago
  indeed, I was using fixed point maths on PlayStation 1 games in the mid to late 90s. It was often responsible for the gaps you'd see between polygons on many PS1 games.
  - marktangotango 10y ago
    I was using it for z-80, game boy color homebrew. Getting down in the bits like this is.... bracing!
    - louthy 10y ago
      I remember doing a Gameboy Color game after working on a PS1 title. And even though I came from an 8 bit 6502 background (BBC Micro), going back to it was hell. 'Bracing' is an understatement! I couldn't imagine much worse these days.
  - rasz 10y ago
    Gaps were mostly a sign of sloppy developers and rushed schedules. Compare Tomb Raider 1 on both PC/PSX with Destruction Derby (reflections) and Gran Turismo 1/2 (polyphony).
jacobolus 10y ago
It’s usually possible to re-frame problems to not require trig functions at all. For instance, you can represent rotations as unit magnitude complex numbers, compose them using complex multiplication, and trivially get whatever trig functions you want out. If you need to compress them for I/O, take the stereographic projection (requires 1 division per point for both forward and inverse transform) and then optionally reduce the precision of the result.
- doubleunplussed 10y ago
  Construction of a unit complex number though, given an angle, requires trigonometry. Precomputing this and re-using it is identical to precomputing the sine and cosine of that angle and reusing them instead - the complex number itself doesn't simplify anything here other than storing both the sine and cosine in one variable.
  - jacobolus 10y ago
    The only time you need to start with an angle is if a human or other external system is feeding it to you, and the only time you need to convert to an angle is when you need to present the data to a human or other external system. The point here is that you can usually get rid of evaluating transcendental functions in the middle of the number crunching part of your code, which is where you would care most about saving operations.
    - doubleunplussed 10y ago
      I don't think that's right. A rotation implies an angle. If it's not coming from a human, then you are either hard-coding an angle to rotate by, in which case you need to take sine and cosine of that angle to generate a complex number, or you're getting the rotation that is implied from some coordinates, hard-coded or otherwise. In which case, you need to normalise a complex number - still requiring a square root. I don't think there's a way to do rotations without transcendental functions (or square roots) that isn't equivalent to merely precomputing the transcendental functions or square roots. Maybe square roots are faster, in which case you have a point that avoiding transcendental functions is the way to go.

shpx 10y ago

Related to 16 bit lookup tables, unum is a format for representing numbers as sets of number ranges as an alternative to floats. http://www.johngustafson.net/presentations/Unums2.0.pptx http://www.johngustafson.net/presentations/Unums2.0.pptx http://www.johngustafson.net/presentations/Unums2.0.pdf http://www.johngustafson.net/presentations/Unums2.0.pdf https://arxiv.org/pdf/1701.00722.pdf https://arxiv.org/pdf/1701.00722.pdf https://en.wikipedia.org/wiki/Interval_arithmetic https://en.wikipedia.org/wiki/Interval_arithmetic

Mindless2112 10y ago
Unums are unlikely to gain much usage. Posits[1][2], also by Gustafson, are a more reasonable alternative to IEEE-754 floating point (but will still have a difficult time displacing IEEE-754, if they can at all). [1] http://web.stanford.edu/class/ee380/Abstracts/170201-slides.pdf http://web.stanford.edu/class/ee380/Abstracts/170201-slides.... [2] https://www.youtube.com/watch?v=aP0Y1uAA-2Y https://www.youtube.com/watch?v=aP0Y1uAA-2Y
- adrianratnapala 10y ago
  Posits seem impressive. My only concern is the lack of NaNs seems like a bug rather than a feature. It's true that some programmers do the silliest things when faced with NaNs. But the fact is they are useful. You often want to do calculations over big matrices, where some elements simply don't have a mathematically defined answer (usually because of div0s but also because input data might have holes). It would be a royal pain if the whole calculation stopped every time you had to add two infinities.
  - dnautics 10y ago
    >is the lack of NaNs seems like a bug rather than a feature. I think this is a reasonable concern. I'll propose to John that we make there be an optional "mode" where the infinity token is treated as "NaN". In reality, this mode just amounts to "ignore NaN traps", because the way that it's done in my hardware models, it requires almost no extra hardware.

weinzierl 10y ago

Back in the day of 80286's and 80386's there was a popular Fractal generating software that did all calculations in fixed point arithmetic because PCs usually had no FPU back then. Fractint[1], the site has a certificate error, but otherwise works fine, even the software still gets updated from time to time. [1] https://fractint.org/ https://fractint.org/

quakeguy 10y ago
Used it to generate uniform random plasma fields for game textures. Great program!
- weinzierl 10y ago
  With the right color map they made a nice sky with clouds. Color cycling them through a smooth full color map was also nice.

emersonrsantos 10y ago

Āryabhaṭa's sine table, the first known mathematical LUT, was made around 500 DC: https://en.wikipedia.org/wiki/Āryabhaṭa's_sine_table https://en.wikipedia.org/wiki/Āryabhaṭa's_sine_table

jacobolus 10y ago
> first known mathematical LUT Babylonian multiplication and reciprocal tables are about 2500 years older. For that matter, Hipparchus and Ptolemy’s chord tables are also centuries older, https://en.wikipedia.org/wiki/Ptolemy's_table_of_chords https://en.wikipedia.org/wiki/Ptolemy's_table_of_chords

raymondh 10y ago

For those who are interested, CORDIC (for COordinate Rotation DIgital Computer) is another technique from that era. It "is a simple and efficient algorithm to calculate hyperbolic and trigonometric functions, typically converging with one digit (or bit) per iteration." See https://en.wikibooks.org/wiki/Digital_Circuits/CORDIC https://en.wikibooks.org/wiki/Digital_Circuits/CORDIC For instructional purposes, here is a simple Python implementation that uses only adds and shifts in the inner loop, followed by a single scaled-multiply to finish it up: https://code.activestate.com/recipes/576792-polar-to-rectangular-conversions-using-cordic https://code.activestate.com/recipes/576792-polar-to-rectang...

bsder 10y ago
IIRC, CORDIC has some issues as you scale the number of bits you want--you can't just run CORDIC 64 times for 64-bits of accuracy, for example. I think the point of diminishing returns is 14-16 bits, but it's been a long time and I'd be happy to be proven wrong.
planteen 10y ago
I used CORDIC on a FPGA a few years ago to compute atan2.

marktangotango 10y ago

>> Jack Crenshaw addresses it in chapter 5 of his fine book "Math Toolkit for Real-Time Programming." Ah shit, can't believe I missed this, all these years. Jack Crenshaw ie Let's Build a Compiler. Love that guy.

marktangotango 10y ago

otempomores 10y ago

Is there a stdclib which allows for easy replacement of expensive funtions like fsin or fcos with lutables? Also do those not poison the cache?

rasz 10y ago

>And by the way, real-world I/O in control situations is never floating-point, whether timers, counters, A/D and D/A converters, servos, etc.. well, actually ;) Adlib/Sound Blaster Yamaha OPL2 and OPL3 both used external Floating Point DACs (YM3014B,YAC512)