Just for fun, here's a recipe for a correctly-rounded nth root operation for positive finite floats. I'm not suggesting using this in the business logic: it's likely way too slow (especially for large n), but it may have a use in the tests. The logic for the Newton iteration in floor_nroot follows that outlined in this Stack Overflow answer: http://stackoverflow.com/a/35276426
from math import frexp, ldexp
def floor_nroot(x, n):
""" For positive integers x, n, return the floor of the nth root of x. """
# Initial guess: here we use the smallest power of 2 that exceeds the nth
# root. But any value greater than or equal to the target result will do.
a = 1 << -(-x.bit_length() // n)
while True:
d = x // a**(n-1)
if a <= d:
return a
a = ((n-1) * a + d) // n
def nroot(x, n):
"""
Correctly-rounded nth root (n >= 2) of x, for a finite positive float x.
"""
if not (x > 0 and n >= 2):
raise ValueError("x should be positive; n should be at least 2")
m, e = frexp(x)
rootm = floor_nroot(int(m * 2**53) << (53*n + (e-1)%n - 52), n)
return ldexp(rootm + rootm % 2, (e-1)//n - 53) |