What is the numerical stability of std::pow() comp

2019-06-21 11:27发布

问题:

What sort of stability issues arise or are resolved by using std::pow()?

  • Will it be more stable (or faster, or at all different) in general to implement a simple function to perform log(n) iterated multiplies if the exponent is known to be an integer?

  • How does std::sqrt(x) compare, stability-wise, to something of the form std::pow(x, k/2)? Would it make sense to choose the method preferred for the above to raise to an integer power, then multiply in a square root, or should I assume that std::pow() is fast and accurate to machine precision for this? If k = 1, is there a difference from std::sqrt()?

  • How would std::pow(x, k/2) or the method above compare, stability-wise, to an integer exponentiation of std::sqrt(x)?

And as a bonus, what are the speed differences likely to be?

回答1:

Will it be more stable (or faster, or at all different) in general to implement a simple function to perform log(n) iterated multiplies if the exponent is known to be an integer?

The result of exponentiation by squaring for integer exponents is in general less accurate than pow, but both are stable in the sense that close inputs produce close results. You can expect exponentiation by squaring to introduce 0.5 ULP of relative error by multiplication (for instance, 1 ULP of error for computing x3 as x * x * x).

When the second argument n is statically known to be 2, then by all means implement xn as x * x. In that case it is faster and more accurate than any possible alternative.

How does std::sqrt(x) compare, stability-wise, to something of the form std::pow(x, k/2)

First, the accuracy of sqrt cannot be beat for an IEEE 754 implementation, because sqrt is one of the basic operations that this standard mandates to be as accurate as possible.

But you are not asking about sqrt, you are asking (I think) about <computation of xn> * sqrt(x) as opposed to pow(x, n + 0.5). Again, in general, for a quality implementation of pow, you can expect pow(x, n + 0.5) to be more accurate than the alternatives. Although sqrt(x) would be computed to 0.5 ULP, the multiplication introduces its own approximation of up to 0.5 ULP, and all in all, it is better to obtain the result you are interested in in a single call to a well-implemented function. A quality implementation of pow will give you 1 ULP of accuracy for its result, and the best implementations will “guarantee” 0.5 ULP.

And as a bonus, what are the speed differences likely to be?

If you know in advance that the exponent is going to be a small integer or multiple of 0.5, then you have information that the implementer of pow did not have, so you can beat them by at least the cost of the test to determine that the second argument is a small integer. Plus, the implementer of a quality implementation is aiming for a more accurate result than simple exponentiation by squaring provides. On the other hand, the implementer of pow can use extremely sophisticated techniques to minimize the average execution time despite the better accuracy: see for instance CRlibm's implementation. I put the verb “guarantee” above inside quotes when talking about the best implementations of pow because pow is one function for which CRlibm's 0.5 ULP accuracy guarantee is only “with astronomical probability”.