Mark Dickinson <[email protected]> added the comment:
> Does the technique you had in mind involve testing 1 ulp up or down to see
> whether its square is closer to the input?
Kinda sorta. Below is some code: it's essentially just pure integer operations,
with a last minute conversion to float (implicit in the division in the case of
the second branch). And it would need to be better tested, documented, and
double-checked to be viable.
def isqrt_rto(n):
"""
Square root of n, rounded to the nearest integer using round-to-odd.
"""
a = math.isqrt(n)
return a | (a*a != n)
def isqrt_frac_rto(n, m):
"""
Square root of n/m, rounded to the nearest integer using round-to-odd.
"""
quotient, remainder = divmod(isqrt_rto(4*n*m), 2*m)
return quotient | bool(remainder)
def sqrt_frac(n, m):
"""
Square root of n/m as a float, correctly rounded.
"""
quantum = (n.bit_length() - m.bit_length() - 1) // 2 - 54
if quantum >= 0:
return float(isqrt_frac_rto(n, m << 2 * quantum) << quantum)
else:
return isqrt_frac_rto(n << -2 * quantum, m) / (1 << -quantum)
----------
_______________________________________
Python tracker <[email protected]>
<https://bugs.python.org/issue45876>
_______________________________________
_______________________________________________
Python-bugs-list mailing list
Unsubscribe:
https://mail.python.org/mailman/options/python-bugs-list/archive%40mail-archive.com
