3

In C, I was wondering whether there's an 'optimal' way to compute half-integer powers. In a nutshell, the problem is computing x^(n/2) (assuming n is odd and decently small, and x is a some float). Is there a major difference in performance/accuracy between sqrt(pow(x, n)) and pow(x, 0.5 * n)? Or even reversed: pow(sqrt(x), n).

Is there some other implementation for handling this specific case of half-integers?

My first thought is that you would just use pow and compute the whole thing in one call, but I feel like with floating point roundoff and things I'm losing some of the precision of the question that comes from the fact that this is explicitly a half-integer. I thought then maybe there's better error performance if you use pow for raising to an integer power and let sqrt handle the (1/2) part.

I also noticed that GSL has functions for computing small integer powers; would combining those functions with sqrt be better than just using pow?

I'm pretty new to scientific programming with C so I'm not sure where I would even go to look for implementations of something like this, and Google hasn't really turned anything up.

phuclv
  • 37,963
  • 15
  • 156
  • 475
  • 1
    Usually the more direct approach is better: `pow()` in this case. Remember when doing math like this, `sqrt()` ends up as just a specialized case of `pow()`. Do be aware that because of the *joys* of floating point math, you will need to pay **very careful attention** to how you round your results. If you get values like 24.999999993851 you'll probably want to interpret that as 25. – tadman Nov 29 '22 at 15:53
  • 2
    For integers I don't recommend the floating point functions at all if you can avoid it. Unfortunately it's kind of hard to avoid with `sqrt` as it's not that easy to compute (but I'll bet it's relatively easy to find good integer algorithms for that too). – Some programmer dude Nov 29 '22 at 15:53
  • 3
    Does not seem to complicated to run some benchmark on performance and precision for the various methods. – klutt Nov 29 '22 at 15:58
  • @Someprogrammerdude sorry, I probably should have specified more clearly, but I'd like `x` in the problem to a be a general double; only the exponent is explicitly integer (divided by two) – Connor Hainje Nov 29 '22 at 15:59
  • Well, once you divide it by two there's a 50% chance of it no longer being an integer, assuming you don't do integer division or `n >> 1`. This is just *ye olde floating point math* at this point. – tadman Nov 29 '22 at 16:03
  • @klutt you're right, I'll try to build a benchmark. I just didn't know if there were any methods _like_ this that were already implemented in some math library I wasn't able to find, or if there was some obvious reasoning I was overlooking. – Connor Hainje Nov 29 '22 at 16:07
  • @tadman true, it's all floating point at this point, but I know that (a) there are more efficient and precise ways to compute small integer powers of floating numbers, and (b) sqrt(x) can sometimes be better than pow(x, 0.5), so I thought that probably there was a 'best' way to compute these half-integer powers. I'll try to run some tests and see if I find anything unexpected. – Connor Hainje Nov 29 '22 at 16:10
  • 5
    Before you go down this rabbit hole, do you have a *measurable performance problem* that *must* be resolved? If you don't, just use `pow()` and move on, live life. If you *are* having performance problems, then you need to dig deeper. `pow()` is probably not the core issue. – tadman Nov 29 '22 at 16:16
  • 1
    @tadman haha that is a good point; I'm absolutely just getting nerd-sniped over something that is probably not an issue at all. I haven't encountered any specific problems; I just got to a point where I was implementing this and thought "_surely_ someone has thought about this before..." – Connor Hainje Nov 29 '22 at 16:19
  • If you must have the most precise way of doing this then a bignum library might be the way to go. Keep in mind that it will be significantly slower than even the slowest native single or double precision math. – SafelyFast Nov 29 '22 at 16:20
  • 1
    Nothing wrong with asking. I think as a general rule here, consider that `sqrt(pow(..))` is two operations, both fairly expensive, while `pow(..., n/2)` is one expensive operation, one very cheap. Further, if you are having problems here because you're doing this operation a lot, consider SIMD or some other vectorized approach. – tadman Nov 29 '22 at 16:48
  • 1
    At least the glibc implementation of `pow` doesn't seem to do anything special for integer powers, so it might as well be called with a half-integer power rather than calling `sqrt` separately. – Ian Abbott Nov 29 '22 at 16:48
  • 3
    I did a small non-scientific test where I computed `pow(x, 3.5)` and `pow(sqrt(x), 7)` for various values of `x`, and compared the results (using glibc). Unsurprisingly, the first one proved to be more accurate in the majority of cases. (Unlike `pow`, `sqrt` is guaranteed to be as accurate as possible, But it still introduces some error, which is then magnified by `pow`.) I think accuracy is more important than speed for most calculations, although YMMV, but the first one is also likely to be faster. So I guess it's the one you should go with. You might want to do a more thorough test. – rici Nov 29 '22 at 17:20
  • Related: https://stackoverflow.com/questions/101439/the-most-efficient-way-to-implement-an-integer-based-power-function-powint-int – klutt Nov 29 '22 at 20:35

1 Answers1

5

Floating-point multiplication is a fairly cheap operation in typical modern processors, and multiplication of an integer by .5 introduces no rounding error when a binary-based floating-point format is used. (If the expression is written as n/2, where n is a floating-point type, I would expect a decent compiler to implement it as multiplication by .5. However, to be sure, it can be written as n*.5.)

pow is a complicated routine but its execution time1 is unlikely to be affected much by a difference between pow(x, n) and either pow(sqrt(x), n) or pow(x, n*.5). We can generally expect pow(x, n*.5) to be a good way of computing xn/2.

sqrt typically takes more execution time than floating-point multiplication and may introduce rounding error. We can expect each of sqrt(pow(x, n)) and pow(sqrt(x), n*.5) to at least as much time as pow(x, n*.5), and probably more, with no benefit in accuracy.

Thus, pow(x, n*.5) is preferred.

I also noticed that GSL has functions for computing small integer powers; would combining those functions with sqrt be better than just using pow?

Maybe. pow is an expensive routine, so customized routines for specific powers could out-perform it, even with sqrt added. This would be situation-dependent, and you would likely have to measure it to know.

Footnote

1 Execution time is not actually a single thing. Performing an operation not only consumes time for that operation but may affect other operations being performed in parallel in modern processors and may affect the start times of later operations, and its own start time may be affected by its relationship with prior operations.

Eric Postpischil
  • 195,579
  • 13
  • 168
  • 312