Smooth approximations of 2^x

[andy-cytomic]

That said, I would argue you probably don't need anywhere NEAR that good accuracy for typical tuning purposes. I used to have a table based approximation with about 0.02 cents and never noticed any issues (even before I eventually added lerps to the thing).

The potentially sensitive areas would be detuned oscillator beating and converting audio-rate modulation signals, where subtle spectral differences can become audible. I'm not sure what's the required precision range for that off the top of my head.

Exactly. I'm using mine in applications where you can set the beat detune of an oscillator in Hz or BPM, so two oscillators beat exactly in time. I also use this for circuit simulation, where for example in an MS20 mk1 filter the signal is passing through the exp of a transistor for tuning and non-linearities.

I've found that computing a 7th order poly works great in f32, but I actually use f64 for the main part of the filter. The trick is to keep precision at f64 until you've subtracted the "dc" part and gotten your input to the poly in the range -1/2 to 1/2 (or 0 to 1), and then convert to f32, computed the poly and power parts and promote back to f64. This means you can do twice as many parallel exp/pow, and still have enough precision for what I need.

I'll double check the 7th order vs 5th order against the standard pow2 math results and let you know the differences at f32 precision.

[andy-cytomic]

My approximations are in the range 0 to 1, but the errors are identical either way. Here is the plot of the error when using double precision to calculate the x value, then f32 to compute the polynomial, as compared to the full double precision math function for pow(2.0, x).

Code: Select all

7th order poly(0,1):
{1.000000000000000, 0.693147180559945, 0.2402265069591007, 0.0555044941070, 0.009615262656, 0.001341316600, 0.000143623130, 0.000021615988}

7th order poly(-1/2,1/2):
{0.999999999615382, 0.69314718125408, 0.24022653230653, 0.0555041019307, 0.0096178619722, 0.0013333699936, 0.0001550537306, 0.000015284812}
 
error min max:
{-1.63515e-07, 1.51524e-07}

pow2 7th order error.png

Code: Select all

5th order poly(0,1):
{1.000000000000000, 0.69315168779795, 0.2401596780318, 0.055817593635, 0.008992164746, 0.001878875789}

5th order poly(-1/2,1/2):
{1.000000088673463, 0.69314693211407, 0.24022037362574, 0.0555072548370, 0.0096798351988, 0.0013285658116}

error min max:
{-2.23442e-07, 2.43317e-07}

pow2 5th order error.png

And here are the same things but done at all f32 and compered to the f32 pow2f

Code: Select all

7th order coeff as above
error min max:
{-1.74702e-07, 1.79029e-07}

pow2 7th order error f32.png

Code: Select all

5th order coeff as above
error min max:
{-2.12815e-07, 2.58052e-07}

pow2 5th order error f32.png

[mystran]

That said, I would argue you probably don't need anywhere NEAR that good accuracy for typical tuning purposes. I used to have a table based approximation with about 0.02 cents and never noticed any issues (even before I eventually added lerps to the thing).

The potentially sensitive areas would be detuned oscillator beating and converting audio-rate modulation signals, where subtle spectral differences can become audible. I'm not sure what's the required precision range for that off the top of my head.

Exactly. I'm using mine in applications where you can set the beat detune of an oscillator in Hz or BPM, so two oscillators beat exactly in time. I also use this for circuit simulation, where for example in an MS20 mk1 filter the signal is passing through the exp of a transistor for tuning and non-linearities.

For beating in Hz, wouldn't it make sense to add the linear offset separately after exponentiation?

As for beating in BPM... any two clocks will drift with regards to each other unless you compute the results in a bit-exact way, so if you want to match host beats, then the only option that truly works is to recompute (with or without gradual adjustment, phase locked loop, whatever) from the host time-stamps on per-block basis.

Even if you did your own computation in exact arithmetics, the host is going to have some rounding errors of it's own and you are going to drift, no matter what. Even if you somehow managed to match it in one host, it's going to drift in the next one.

[2DaT]

Yeah, for 1/octave single-precision control signal, with zero-point somewhere around the middle of the hearing range, the precision 4-5 octaves up or down (which is about the audible range) gets you about 0.0005 cents worth of accuracy before you even exponentiate. That said, I would argue you probably don't need anywhere NEAR that good accuracy for typical tuning purposes. I used to have a table based approximation with about 0.02 cents and never noticed any issues (even before I eventually added lerps to the thing).

That said, the approximation I've been using lately is 2DaT's exp_mp (this one I think: viewtopic.php?p=7161124#p7161124) which is about 3(?) ulps for single precision.

Well these are the coefficients for the linked approximation for 2^x in the range -1/2 to 1/2 (shown orange):
{1.000000119000000, 0.6931469440000000, 0.2402212024000000, 0.05550713092000000, 0.009675540961000000, 0.001327647245000000}

The above isn't even C0. If you want something with a lower error try this one which is at least C0 continuous, and has a more symmetric error (shown blue):
{1.000000088673463, 0.69314693211407, 0.24022037362574, 0.0555072548370, 0.0096798351988, 0.0013285658116}

With single precision coefficients, poly will drift a bit. C0 doesn't really matter at this point, because there are rounding errors in polynomial evaluation. You may get some precision by going 7th order, but error is still limited by eval errors (1-1.5ulp).

[KBonzo]

Sure memory is cheap, but cache misses are expensive. There are intrinsics to do all the "tricks", so there is no undefined behaviour, and you can do them on a f32x4 or higher width SIMD as well, so you get at least 4 lots of 2^x for the price of 1 using the same small number of coefficients for all of them, so it stays in the cache easily.

The best way to know for sure which is quicker is to test both methods. In my testing using polynomial approximations, or even rational approximations ends up being quicker than table based approaches, but I still use tables for pre-computing complicated shapes.

I take your point about SIMD but cache sizes continue to grow and use prefetching so it becomes less of an issue. It's better to be without the extra headache unless the speedup is significant. Like you say the proof is to benchmark it. Mytel above wrote "I was worried that the previous approximation I was using wasn't smooth". You need more certainty that that.

[Z1202]

the proof is to benchmark it.

Hmmm, I've yet to see a convincing strategy to benchmark the potential cache misses. It is so much context-dependent. Especially if you're running a loop with a large sample-by-sample body. Thus a table which works nice in a small effect might not work so well in a context of a large feedback loop. Plus, if the cache is shared between processors I don't see a reasonable way to predict this at all. Plus potential cache aliasing issues.

[Z1202]

Wrong button, deleted

[KBonzo]

Just curious about how the minimax algorithm is used to fit a polynomial? I usually use Octave but it doesn't have the matlab fminimax function. In what way would this be better than matlab's polyfit?

[juha_p]

Just curious about how the minimax algorithm is used to fit a polynomial? I usually use Octave but it doesn't have the matlab fminimax function. In what way would this be better than matlab's polyfit?

Isn't polyfitinf() (Optim -package) same as minimax algorithm?

[KBonzo]

Isn't polyfitinf() (Optim -package) same as minimax algorithm?

Looks like it is. : "That polynomial, also known as minimax polynomial, is obtained by the exchange algorithm"

Polyfit uses least squares.

Thanks for that.