Smooth approximations of 2^x

[Z1202]

Following the neighbor thread viewtopic.php?f=33&t=519728
I thought I post some approximations of 2^x which I did a couple of years ago. The key feature of these approximations are that they are smooth at the seams, which makes them particularly suitable for quite a few music DSP needs even at a relatively low precision, which means saving some CPU.

I post the results in the form that I have them already (PDF). My apologies for the small inconvenience for having to download and unpack the file.

[mtytel]

Seems to give pretty good results. I was worried that the previous approximation I was using wasn't smooth so I'm trying this out for now.

I also created a similar method for log2.
It uses the fact that you can compute the log2 of the exponent of a floating point number (the floored log) separately from the mantissa and then add the results.

Code: Select all

float log2(float f) {
  union {
    int as_int;
    float as_float;
  } v;
  v.as_float = f;
  int floored_log = (v.as_int >> 23) - 127;
  
  v.as_int = (v.as_int & 0x7fffff) | (0x7f << 23);
  float mantissa = v.as_float;
  
  // ...
}

Then use a similar polynomial approximation for the mantissa value. This time you keep the derivatives of the function evaluated at 2 half of the derivatives of the function evaluated at 1.
f'(1) = 2f'(2)

Here are the coefficients for the N=5 case:
Coefficient0 = -1081 / 541
Coefficient1 = 1495 / 541
Coefficient2 = -510 / 541
Coefficient3 = 110 / 541
Coefficient4 = -15 / 541
Coefficient5 = 1 / 541

[quikquak]

Some might find this useful:
An integer double log2

Code: Select all

        int32 ilog2(double x)
	{
		return (int32)((((reinterpret_cast<unsigned long long&>(x)) >> 52) & 0x7ff)-1023);
	}

I'm not sure if it'll work on Mac, you might need to do the same union trick, but it's fine on VS-2017.
I use it to find the rough octave number of a given frequency.

[Z1202]

I also created a similar method for log2.

Hmmm, I ended up with somewhat different coefficients which seem to give precision better by a magnitude order (ca. 0.0005 vs. 0.004 in your case).
x^5/31-x^4/3+(10*x^3)/7-(10*x^2)/3+5*x-1819/651
Interestingly, the precision is much worse than for a comparable order polynomial for 2^x.

Edit:

This time you keep the derivatives of the function evaluated at 2 half of the derivatives of the function evaluated at 1.
f'(1) = 2f'(2)

I got f^(n)(1) = 2^n f^(n)(2), probably that's the reason

[mystran]

Some might find this useful:
An integer double log2

Code: Select all

        int32 ilog2(double x)
	{
		return (int32)((((reinterpret_cast<unsigned long long&>(x)) >> 52) & 0x7ff)-1023);
	}

I'm not sure if it'll work on Mac, you might need to do the same union trick, but it's fine on VS-2017.
I use it to find the rough octave number of a given frequency.

I checked this with clang (on macOS using -Ofast) and it works fine (and no warnings with -Wall), although I replaced "unsigned long long" with "uint64_t" and "int32" with "int32_t" to use standard stdint.h types of predictable sizes.

I never realised you could reinterpret_cast references like this, but I guess it kinda makes sense and apparently it works.

[quikquak]

Thanks for checking, mystran.

[Z1202]

Some might find this useful:
...

I checked this with clang (on macOS using -Ofast) and it works fine (and no warnings with -Wall)...

I'd highly recommend using the "may_alias" attribute whenever doing this kind of reinterpret_cast tricks. C++ standard defines the results of such reinterpret casting as undefined behavior, and some compilers take this as a permission to produce completely random results in this case. This comes from my real experience with exactly this kind of code.

[Z1202]

While we're on the subject... My take on the C++ standard

[2DaT]

IMO, type punning through unions is better. Though union type punning is illegal in c++, AFAIK, major compilers extend the c++ standard to explicitly allow it.

[stratum]

About may_alias:

The fact that this behavior is undefined is pretty annoying. (https://blog.qt.io/blog/2011/06/10/type ... -aliasing/ )

I have always assumed that it is defined exactly as you would expect: pointer dereferenced from an incompatible type, i.e. reading the same memory region, the result being dependent upon the endianness of the processor.

[juha_p]

https://www.youtube.com/watch?v=sCjZuvtJd-k

[quikquak]

Interesting, so digging around inside the float using a union is also undefined.
*edit* I've worked with GPUs quite a bit and it amazes me that people just presume all cards use IEEE standard storage or even accuracy.

[mystran]

While we're on the subject... My take on the C++ standard

With modern compilers it's more like what you get from that image when you replace all the empty squares with mines.

[Max M.]

quikquak
I don't think in context of this thread anyone assumes that the proposed code would work in a non-IEEE-754 environment.

[juha_p]

...
I thought I post some approximations of 2^x which I did a couple of years ago. The key feature of these approximations are that they are smooth at the seams, which makes them particularly suitable for quite a few music DSP needs even at a relatively low precision, which means saving some CPU.
...

In comparison I calculated two 2^x approximations using Maclaurin series (coefficients are actually for sqrt(2^x) ... at different points) as shown in plot:
https://www.desmos.com/calculator/q49errotql

... but, how do I see/know if it's "smooth" or "smoother" than the one found in your linked PDF (equal degree approximation plotted as comparison)?