Smooth approximations of 2^x
-
- KVRAF
- Topic Starter
- 1607 posts since 12 Apr, 2002
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.
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.
You do not have the required permissions to view the files attached to this post.
- KVRist
- 186 posts since 28 Jan, 2013 from Oakland
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.
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
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;
// ...
}
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
- KVRian
- 872 posts since 6 Aug, 2005 from England
Some might find this useful:
An integer double log2
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.
An integer double log2
Code: Select all
int32 ilog2(double x)
{
return (int32)((((reinterpret_cast<unsigned long long&>(x)) >> 52) & 0x7ff)-1023);
}
I use it to find the rough octave number of a given frequency.
Dave Hoskins. http://www.quikquak.com
-
- KVRAF
- Topic Starter
- 1607 posts since 12 Apr, 2002
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:
I got f^(n)(1) = 2^n f^(n)(2), probably that's the reason
- KVRAF
- 7892 posts since 12 Feb, 2006 from Helsinki, Finland
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.quikquak wrote: ↑Mon Feb 18, 2019 10:38 pm Some might find this useful:
An integer double log2I'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.Code: Select all
int32 ilog2(double x) { return (int32)((((reinterpret_cast<unsigned long long&>(x)) >> 52) & 0x7ff)-1023); }
I use it to find the rough octave number of a given frequency.
I never realised you could reinterpret_cast references like this, but I guess it kinda makes sense and apparently it works.
- KVRian
- 872 posts since 6 Aug, 2005 from England
-
- KVRAF
- Topic Starter
- 1607 posts since 12 Apr, 2002
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.
-
- KVRAF
- 2256 posts since 29 May, 2012
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.
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.
~stratum~
-
- KVRian
- 834 posts since 21 Feb, 2006 from FI
- KVRian
- 872 posts since 6 Aug, 2005 from England
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.
*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.
Dave Hoskins. http://www.quikquak.com
- KVRAF
- 7892 posts since 12 Feb, 2006 from Helsinki, Finland
With modern compilers it's more like what you get from that image when you replace all the empty squares with mines.
- KVRist
- 347 posts since 20 Apr, 2005 from Moscow, Russian Federation
quikquak
I don't think in context of this thread anyone assumes that the proposed code would work in a non-IEEE-754 environment.
I don't think in context of this thread anyone assumes that the proposed code would work in a non-IEEE-754 environment.
-
- KVRian
- 834 posts since 21 Feb, 2006 from FI
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:Z1202 wrote: ↑Fri Feb 08, 2019 1:39 pm ...
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.
...
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)?