Fast pow2 for pitch to frequency mapping

Ok, so I posted a thread about this in KVR, and obvious I'm not the only one to have thought about it (which was to be expected ofcourse), but I'm still posting this here for future reference, since I don't have anything else to post about today.

In synthesizers, being able to convert from pitch to frequency is necessary. In equal tempered scale this involves the relation f=b2^p/12, where f is the frequency, b is the tuning base, and p is the pitch in semitones relative to the base tuning. This is trivially solved in C++ by f*pow(2.0,p/12). Trouble is, pow() takes an awful lot of time, so it's not really realistic if you want to do it on per-sample, per-voice basis.

One would expect to find some nice code that solves this problem reasonably accurate for the purposes of synthesizers, yet still reasonably fast. I've tried doing this a few times, but always failed. By reasonably accurate I mean something that gives resolution of about 1/100 semitones over the range of hearing (around 10 octaves). Most solutions I've seen are either

1. wildly inaccurate or only accurate near 0
2. only marginally faster (often slower than) than what gcc gives for pow

So all hope is lost? Not really, if one is willing to pay for a few lookup tables.

Now, looking at the facts, the accuracy requirements are measured in fractions of semitones. This suggests rounding before exponentiation should not cause problems, if the exponential itself can be calculated accurately. So how much accuracy does one need?

First of all, one only needs to consider one side, since 2^-x = 1/2^x. Let's say we say we want to go up and down 10 octaves, we need at least 4 bits for the integer part. As for the fractional part, 1 octave = 12 semitones, and 1 semitones = 100 cents. Hence there are 1200 cents to an octave. For 12 bits one gets 4096 tunings between octave. That's a total of 16 bits, and more than enough.

Full lookup table for 16 bits would ofcourse take 65536 entries, which for 32-bit floats takes a 256kB. Not huge these days, but exponentials have the nice property that a^(x+y)=a^xa^y. And this ofcourse means one can do a separate table for the upper 8 bits and the lower 8 bits. That's 256 entries per table, both tables together taking 2kB of memory.

So the algorithm for pow2(x) with two tables of 8 bits each then goes:

1. if x<0 recurse for 1/2^-x
   
2. multiply x by 2¹²
3. round/floor/whatever x into an integer (in C++ cast works fine)
4. do table lookups with 8 lowest, and 8 next bits, and multiply them together
5. return result

Ofcourse one can use more than two tables for more accuracy or lower memory requirements. How to populate the lookup tables is left as an exercise for the reader.