61–70 of 140 posts
[flagged]
And similarly, entire generations of programmers were never taught Horner's scheme. You can see it in the article, where they write stuff like A * x * x * x * x * x * x + B * x * x * x * x + C * x * x + D (10 muls, 3 muladds) instead of the faster tmp = x * x; ((A * tmp + B) * tmp + C) * tmp + D (1 mul, 3 muladds)
In general, I find that minimax approximation is an underappreciated tool, especially the quite simple Remez algorithm to generate an optimal polynomial approximation [0]. With some modifications, you can adapt it to optimize for either absolute or relative error within an interval, or even come up with rational-function approximations. (Though unfortunately, many presentations of the algorithm use overly-simple form…
They teach a lot of Taylor/Maclaurin series in Math classes (and trig functions are sometimes called "CORDIC" which is an old method too) but these are not used much in actual FPUs and libraries. Maybe we should update the curricula so people know better ways.
[flagged]
And similarly, entire generations of programmers were never taught Horner's scheme. You can see it in the article, where they write stuff like A * x * x * x * x * x * x + B * x * x * x * x + C * x * x + D (10 muls, 3 muladds) instead of the faster tmp = x * x; ((A * tmp + B) * tmp + C) * tmp + D (1 mul, 3 muladds)
But if you're writing in an interpreted language that doesn't have a good JIT, or for a platform with a custom compiler, it might be worth hand-tweaking expressions with an eye towards performance and precision.
Earlier quoted context omitted.
I experimented a bit with the code. Various tables with different datatypes. There is enough noise from the Monte Carlo to not make a difference if you use smaller data types than double or float. Even dropping interpolation worked fine, and got the speed to be on par with the best in the article, but not faster.
Does your benchmark use sequential or randomly ordered inputs? That would make a substantial difference with an LUT, I would think. But I'm guessing. Maybe 32K is so small it doesn't matter (if almost all of the LUT sits in the cache and is never displaced). > if you use smaller data types than double or float. Even dropping interpolation worked fine, That's kinda tautological isn't it? Of course reduced precision is…
About dropping the interpolation: Yes you are right of course. I was thinking about the speed. No noticable speed improvement by dropping interpolation. The asin calls are only a small fraction of everything.
[flagged]
And similarly, entire generations of programmers were never taught Horner's scheme. You can see it in the article, where they write stuff like A * x * x * x * x * x * x + B * x * x * x * x + C * x * x + D (10 muls, 3 muladds) instead of the faster tmp = x * x; ((A * tmp + B) * tmp + C) * tmp + D (1 mul, 3 muladds)
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.
And mostly it has been OK, except for some cases like games or simulations that want to get bitwise identical results across HW, which (if they're lucky) just don't use these operations or (if they're unlucky) use them and have to handle mismatches somehow. Compilers never generate these operations implicitly unless you're compiling with some sort of fast-math flag, so you mostly only get to them by explicitly using an intrinsic, and in theory you know what you're signing up for if you do that.
However, this did make them unusable for some scenarios where you would otherwise like to use them, so a bunch of graphics and scientific computing and math library developers said "please fully specify these operations next time" and now NEON/SVE and AVX512 have fully-specified reciprocal estimates,¹ which solves the problem unless you have to interoperate between x86 and ARM.
¹ e.g. Intel "specifies" theirs here: https://www.intel.com/content/www/us/en/developer/articles/c...
ARM's is a little more readable: https://developer.arm.com/documentation/ddi0596/2021-03/Shar...
Earlier quoted context omitted.
And similarly, entire generations of programmers were never taught Horner's scheme. You can see it in the article, where they write stuff like A * x * x * x * x * x * x + B * x * x * x * x + C * x * x + D (10 muls, 3 muladds) instead of the faster tmp = x * x; ((A * tmp + B) * tmp + C) * tmp + D (1 mul, 3 muladds)
Thinking about speed like this used to be necessary in C and C++ but these days you should feel free to write the most legible thing (Horner's form) and let the compiler find the optimal code for it (probably similar to Horner's form but broken up to have a shallower dependency chain). But if you're writing in an interpreted language that doesn't have a good JIT, or for a platform with a custom compiler, it might be…
Earlier quoted context omitted.
And similarly, entire generations of programmers were never taught Horner's scheme. You can see it in the article, where they write stuff like A * x * x * x * x * x * x + B * x * x * x * x + C * x * x + D (10 muls, 3 muladds) instead of the faster tmp = x * x; ((A * tmp + B) * tmp + C) * tmp + D (1 mul, 3 muladds)
The common subexpression elimination (CSE) pass in compilers takes care of that.
Earlier quoted context omitted.
And similarly, entire generations of programmers were never taught Horner's scheme. You can see it in the article, where they write stuff like A * x * x * x * x * x * x + B * x * x * x * x + C * x * x + D (10 muls, 3 muladds) instead of the faster tmp = x * x; ((A * tmp + B) * tmp + C) * tmp + D (1 mul, 3 muladds)
The problem with Horner’s scheme is that it creates a long chain of data dependencies, instead of making full use of all execution units. Usually you’d want more of a binary tree than a chain.
_Can_ you even make a reasonable high-ILP scheme for a polynomial, unless it's of extremely high degree?
Earlier quoted context omitted.
And similarly, entire generations of programmers were never taught Horner's scheme. You can see it in the article, where they write stuff like A * x * x * x * x * x * x + B * x * x * x * x + C * x * x + D (10 muls, 3 muladds) instead of the faster tmp = x * x; ((A * tmp + B) * tmp + C) * tmp + D (1 mul, 3 muladds)
Didn't know this technique had a name, but I would think a modern compiler could make this optimization on its own, no?