# Correctly rounded integer division

from math import ldexp

def nbits(n, correction = {
        '0': 4, '1': 3, '2': 2, '3': 2,
        '4': 1, '5': 1, '6': 1, '7': 1,
        '8': 0, '9': 0, 'a': 0, 'b': 0,
        'c': 0, 'd': 0, 'e': 0, 'f': 0}):
    """Number of bits in binary representation of the positive integer n,
    or 0 if n == 0.
    """
    if n < 0:
        raise ValueError("The argument to nbits should be nonnegative.")
    hex_n = "%x" % n
    return 4*len(hex_n) - correction[hex_n[0]]

def truediv(a, b):
    """Find the closest float to a/b, for integers a and b."""
    sign = -1 if (a > 0) ^ (b > 0) else 1
    a = abs(a)
    b = abs(b)
    # NBITS_WANTED is an upper bound for number of bits in the
    # mantissa of a float, plus 1; see _PyLong_AsScaledDouble in
    # longobject.c
    NBITS_WANTED = 57 
    shift = nbits(b)-nbits(a)+NBITS_WANTED
    if shift >= 0:
        q, r = divmod(a << shift, b)
    else:
        q, r = divmod(a, b << -shift)
    return ldexp(float(sign*(2*q+bool(r))), -(shift+1))