On Sun, Aug 14, 2016 at 08:29:39AM +0000, Mark Dickinson wrote:
> Steven: can you explain why you think your code *should* be giving 
> exact results for exact powers? Do you have an error analysis that 
> says that should be the case?

No error analysis, only experimental results.

> One issue here is that libm pow functions vary hugely in quality, so 
> any algorithm that depends on ** or math.pow is going to be hard to 
> prove anything about.

pow() is just the starting point. An earlier version of the function 
started with an initial guess of x for the nth root of x, but it was 
very slow, and sometimes failed to converge in any reasonable time. By 
swapping to an initial guess calculated with pow(), I cut the number of 
iterations down to a much smaller number.

Because I'm not to worried about being super-fast, the nth root function 
goes through the following steps:

- use y = pow(x, 1/n) for an initial guess;
- if y**n == x, return y;
- improve the guess with the "Nth root algorithm";
- if that doesn't converge after a while, return y;
- finally, compare the epsilons abs(y**n - x) for the initial guess
  and improved version, returning whichever gives the smaller epsilon.

So I'm confident that nth_root() should never be worse than pow().

> I think the tests should simply be weakened: it's unreasonable to 
> expect perfect results in this case.

Okay. I'll weaken the tests.