Live data from Hacker News

Faster asin() was hiding in plain sight

16bpp.net

131–140 of 140 posts

Re: Faster asin() was hiding in plain sight

#131

I'm pretty sure it's not faster, but it was fun to write: float asin(float x) { float x2 = 1.0f-fabs(x); u32 i = bitcast(x2); i = 0x5f3759df - (i>>1); float inv = bitcast(i); return copysign(pi/2-pi/2*(x2*inv),x); } Courtesy of evil floating point bithacks.

The bad thing about this method is that it's slower than native CPU instructions. The good thing is that the result is very precise for at least 2 values of x, namely 1.0 and -1.0 JK

Yeah, it's an interesting approximation. It follows the shape of the curve pretty well except around 0, much better than low degree polynomials. Polynomial approximations have much better relative error bounds though.

Re: Faster asin() was hiding in plain sight

#132
post #44

To be accurate, this is originally from Hastings 1955, Princeton "APPROXIMATIONS FOR DIGITAL COMPUTERS BY CECIL HASTINGS", page 159-163, there are actually multiple versions of the approximation with different constants used. So the original work was done with the goal of being performant for computers of the 1950's. Then the famous Abramowitz and Stegun guys put that in formula 4.4.45 with permission, then the nvidi…

I ran this down, because I have a particular interest in vectorizable function approximations. Particular those that exploit bit-banging to handle range normalization. (Anyone have a good reference for that?) Regrettably, this is NOT from Hastings 1955. Hastings provides Taylor series and Chebyshev polynomial approximations. The OP's solution is a Pade approximation, which are not covered at all in Hastings.

>Particular those that exploit bit-banging to handle range normalization

https://userpages.cs.umbc.edu/phatak/645/supl/Ng-ArgReductio...

That's tiny but weird.

Re: Faster asin() was hiding in plain sight

#133
post #58

> After all of the above work and that talk in mind, I decided to ask an LLM. Impressive that an LLM managed to produce the answer from a 7 year old stack overflow answer all on its own! [1] This would have been the first search result for “fast asin” before this article was published. [1]: https://stackoverflow.com/a/26030435

[deleted]

Re: Faster asin() was hiding in plain sight

#134
post #101

Does anyone knows the resources for the algos used in the HW implementations of math functions? I mean the algos inside the CPUs and GPUs. How they make a tradeoff between transistor number, power consumption, cycles, which algos allow this.

Here's one way to do it. https://en.wikipedia.org/wiki/CORDIC

Thanks, but this seems to be optimized for the smallest number of gates, so it applies for simple microcontrollers and FPGA, and with limited precision. I was interested in actual state of the art used in modern CPUs and GPUs.

Re: Faster asin() was hiding in plain sight

#136

Isn't the faster approach SIMD [edit: or GPU]? A 1.05x to 1.90x speedup is great. A 16x speedup is better! They could be orthogonal improvements, but if I were prioritizing, I'd go for SIMD first. I searched for asin on Intel's intrinsics guide. They have a AVX-512 instrinsic `_mm512_asin_ps` but it says "sequence" rather than single-instruction. Presumably the actual sequence they use is in some header file somewher…

The issue is that the algorithm is only half the story. The implementation (e.g. bytecode) is the other.

I've been trying to find ways to make the original graphics renderer of the CGA version of Elite faster as there have been dozens of little optimizations found over the decades since it was written.

I was buoyed by a video of Super Mario 64/Zelda optimizations where it was pointed out that sometimes an approx calculation of a trig value can be quicker than a table lookup depending on the architecture.

Based on that I had conversations with LLMs over what fast trig algorithms there are, but for 8088 you are cooked most of the time on implementing them at speed.

Re: Faster asin() was hiding in plain sight

#137
The glibc implementation already has tests for several ranges and hacks for them including Taylor series:

https://github.com/lattera/glibc/blob/master/sysdeps/ieee754...

The smallest range is |x| So right there, if we neglect this detail in our own wrapper, we may be able to get a speedup, at the cost of sending very small values to the Tayor series.

The next smallest range tested is |x| The authors are cetainly not missing any brilliant trick hiding in plain sight; they are doing a more assiduous job.

Re: Faster asin() was hiding in plain sight

#138
In that Padé approximant I think you can save a couple multiplications.

As written it does this:

  n = 1 - 367/714 * x**2
  d = 1 - 81/119 * x**2 + 183/4760 * x**4
  return x * (n/d)
That's got 7 multiplies (I'm counting divide as a multiply) and 3 additions. (I'm assuming the optimizer only computes x^2 once and computes x^2 by squaring x^2, and that all the constants are calculated at compile time).

Replace n/d with 1/(d/n) and then replace d/n with q + r/n where q is the quotient of polynomial d divided by polynomial n and r is the remainder.

This is the result:

  n = 1 - 367/714 * x**2
  q = 1587627/1346890 - 549/7340 * x**2
  r = -240737/1346890
  return x / (q + r/n)
That's got 5 multiplies and 3 additions.

Re: Faster asin() was hiding in plain sight

#139
post #23

While I'm glad to see the OP got a good minimax solution at the end, it seems like the article missed clarifying one of the key points: error waveforms over a specified interval are critical, and if you don't see the characteristic minimax-like wiggle, you're wasting easy opportunity for improvement. Taylor series in general are a poor choice, and Pade approximants of Taylor series are equally poor. If you're going t…

Thanks so much for that!

I've been struggling to curve fit an aerodynamics equation relating Mach number to rocket nozzle exit/entrance area ratio for quite some time. It's a 5th or 6th degree polynomial whose inverse doesn't have a closed-form solution:

https://en.wikipedia.org/wiki/Abel%E2%80%93Ruffini_theorem

But I was able to use a Chebyshev fit that is within a few percent accurate at 3rd degree, and is effectively identical at 4th degree or higher. And 4th degree (quartic) polynomials do have a closed-form solution for their inverse. That lets me step up to higher abstractions for mass flow rate, power, etc without having to resort to tables. At most, I might need to use piecewise-smooth sections, which are far easier to work with since they can just be dropped into spreadsheets, used for derivatives/integrals, etc.

Anyway, I also discovered (ok AI mentioned) that the Chebyshev approximation is based on the discrete cosine transform (DCT):

https://en.wikipedia.org/wiki/Discrete_Chebyshev_transform#R...

https://en.wikipedia.org/wiki/Discrete_cosine_transform#Appl...

That's why it's particularly good at curve-fitting in just 2 or 3 terms. Which is why they use the DCT for image compression in JPG etc:

https://www.mathworks.com/help/images/discrete-cosine-transf...

The secret sauce is that Chebyshev approximation spreads the error as ripples across the function, rather than at its edges like with Taylor series approximation. That helps it fit more intricate curves and arbitrary data points, as well as mesh better with neighboring approximations.

Re: Faster asin() was hiding in plain sight

#140

Earlier quoted context omitted.

Yeah, the only big problem with approx. sqrt is that it's not consistent across systems, for example Intel and AMD implement RSQRT differently... Fine for graphics, but if you need consistency, that messes things up.

Wait, what? Do you have a resource I could read up on about that? That is moderately concerning if your math isn't portable across chips.

Take a look at the "rsqrt_rcp" section of reference [6] in the accuracy report by Gladman et al referenced above. I did that work 10 years ago because some people at CERN had reported getting different results from certain programs depending on whether the exact same executables were run on Intel or AMD cpus. The result of the investigation was that the differing results were due to different implementations of the rsqrt instruction on the different cpus.
Post reply on HN