Message272782
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) |
|
Date |
User |
Action |
Args |
2016-08-15 17:11:47 | mark.dickinson | set | recipients:
+ mark.dickinson, tim.peters, rhettinger, steven.daprano, martin.panter |
2016-08-15 17:11:47 | mark.dickinson | set | messageid: <1471281107.9.0.887121560789.issue27761@psf.upfronthosting.co.za> |
2016-08-15 17:11:47 | mark.dickinson | link | issue27761 messages |
2016-08-15 17:11:47 | mark.dickinson | create | |
|