def test_math_sqrt(self):
     def f(x):
         try:
             return math.sqrt(x)
         except ValueError:
             return -INFINITY
     
     res = self.interp_operations(f, [0.0])
     assert res == 0.0
     self.check_operations_history(call_pure=1)
     #
     res = self.interp_operations(f, [25.0])
     assert res == 5.0
     self.check_operations_history(call_pure=1)
     #
     res = self.interp_operations(f, [-0.0])
     assert str(res) == '-0.0'
     self.check_operations_history(call_pure=1)
     #
     res = self.interp_operations(f, [1000000.0])
     assert res == 1000.0
     self.check_operations_history(call_pure=1)
     #
     res = self.interp_operations(f, [-1.0])
     assert res == -INFINITY
     self.check_operations_history(call_pure=0)
     #
     res = self.interp_operations(f, [INFINITY])
     assert isinf(res) and not isnan(res) and res > 0.0
     self.check_operations_history(call_pure=0)
     #
     res = self.interp_operations(f, [NAN])
     assert isnan(res) and not isinf(res)
     self.check_operations_history(call_pure=0)

def c_cosh(x, y):
    if not isfinite(x) or not isfinite(y):
        if isinf(x) and isfinite(y) and y != 0.:
            if x > 0:
                real = copysign(INF, math.cos(y))
                imag = copysign(INF, math.sin(y))
            else:
                real = copysign(INF, math.cos(y))
                imag = -copysign(INF, math.sin(y))
            r = (real, imag)
        else:
            r = cosh_special_values[special_type(x)][special_type(y)]

        # need to raise ValueError if y is +/- infinity and x is not
        # a NaN
        if isinf(y) and not isnan(x):
            raise ValueError("math domain error")
        return r

    if fabs(x) > CM_LOG_LARGE_DOUBLE:
        # deal correctly with cases where cosh(x) overflows but
        # cosh(z) does not.
        x_minus_one = x - copysign(1., x)
        real = math.cos(y) * math.cosh(x_minus_one) * math.e
        imag = math.sin(y) * math.sinh(x_minus_one) * math.e
    else:
        real = math.cos(y) * math.cosh(x)
        imag = math.sin(y) * math.sinh(x)
    if isinf(real) or isinf(imag):
        raise OverflowError("math range error")
    return real, imag

def c_rect(r, phi):
    if not isfinite(r) or not isfinite(phi):
        # if r is +/-infinity and phi is finite but nonzero then
        # result is (+-INF +-INF i), but we need to compute cos(phi)
        # and sin(phi) to figure out the signs.
        if isinf(r) and isfinite(phi) and phi != 0.:
            if r > 0:
                real = copysign(INF, math.cos(phi))
                imag = copysign(INF, math.sin(phi))
            else:
                real = -copysign(INF, math.cos(phi))
                imag = -copysign(INF, math.sin(phi))
            z = (real, imag)
        else:
            z = rect_special_values[special_type(r)][special_type(phi)]

        # need to raise ValueError if r is a nonzero number and phi
        # is infinite
        if r != 0. and not isnan(r) and isinf(phi):
            raise ValueError("math domain error")
        return z

    real = r * math.cos(phi)
    imag = r * math.sin(phi)
    return real, imag

def c_sinh(x, y):
    # special treatment for sinh(+/-inf + iy) if y is finite and nonzero
    if not isfinite(x) or not isfinite(y):
        if isinf(x) and isfinite(y) and y != 0.:
            if x > 0:
                real = copysign(INF, math.cos(y))
                imag = copysign(INF, math.sin(y))
            else:
                real = -copysign(INF, math.cos(y))
                imag = copysign(INF, math.sin(y))
            r = (real, imag)
        else:
            r = sinh_special_values[special_type(x)][special_type(y)]

        # need to raise ValueError if y is +/- infinity and x is not
        # a NaN
        if isinf(y) and not isnan(x):
            raise ValueError("math domain error")
        return r

    if fabs(x) > CM_LOG_LARGE_DOUBLE:
        x_minus_one = x - copysign(1., x)
        real = math.cos(y) * math.sinh(x_minus_one) * math.e
        imag = math.sin(y) * math.cosh(x_minus_one) * math.e
    else:
        real = math.cos(y) * math.sinh(x)
        imag = math.sin(y) * math.cosh(x)
    if isinf(real) or isinf(imag):
        raise OverflowError("math range error")
    return real, imag

def round(space, number, w_ndigits=0):
    """round(number[, ndigits]) -> floating point number

Round a number to a given precision in decimal digits (default 0 digits).
This always returns a floating point number.  Precision may be negative."""
    # Algorithm copied directly from CPython

    # interpret 2nd argument as a Py_ssize_t; clip on overflow
    ndigits = space.getindex_w(w_ndigits, None)

    # nans, infinities and zeros round to themselves
    if number == 0 or isinf(number) or isnan(number):
        return space.wrap(number)

    # Deal with extreme values for ndigits. For ndigits > NDIGITS_MAX, x
    # always rounds to itself.  For ndigits < NDIGITS_MIN, x always
    # rounds to +-0.0.
    if ndigits > NDIGITS_MAX:
        return space.wrap(number)
    elif ndigits < NDIGITS_MIN:
        # return 0.0, but with sign of x
        return space.wrap(0.0 * number)

    # finite x, and ndigits is not unreasonably large
    z = round_double(number, ndigits)
    if isinf(z):
        raise OperationError(space.w_OverflowError, space.wrap("rounded value too large to represent"))
    return space.wrap(z)

def c_exp(x, y):
    if not isfinite(x) or not isfinite(y):
        if isinf(x) and isfinite(y) and y != 0.:
            if x > 0:
                real = copysign(INF, math.cos(y))
                imag = copysign(INF, math.sin(y))
            else:
                real = copysign(0., math.cos(y))
                imag = copysign(0., math.sin(y))
            r = (real, imag)
        else:
            r = exp_special_values[special_type(x)][special_type(y)]

        # need to raise ValueError if y is +/- infinity and x is not
        # a NaN and not -infinity
        if isinf(y) and (isfinite(x) or (isinf(x) and x > 0)):
            raise ValueError("math domain error")
        return r

    if x > CM_LOG_LARGE_DOUBLE:
        l = math.exp(x-1.)
        real = l * math.cos(y) * math.e
        imag = l * math.sin(y) * math.e
    else:
        l = math.exp(x)
        real = l * math.cos(y)
        imag = l * math.sin(y)
    if isinf(real) or isinf(imag):
        raise OverflowError("math range error")
    return real, imag

def test_nan_and_special_values():
    from pypy.rlib.rfloat import isnan, isinf, isfinite, copysign
    inf = 1e300 * 1e300
    assert isinf(inf)
    nan = inf/inf
    assert isnan(nan)

    for value, checker in [
            (inf,   lambda x: isinf(x) and x > 0.0),
            (-inf,  lambda x: isinf(x) and x < 0.0),
            (nan,   isnan),
            (42.0,  isfinite),
            (0.0,   lambda x: not x and copysign(1., x) == 1.),
            (-0.0,  lambda x: not x and copysign(1., x) == -1.),
            ]:
        def f():
            return value
        f1 = compile(f, [])
        res = f1()
        assert checker(res)

        l = [value]
        def g(x):
            return l[x]
        g2 = compile(g, [int])
        res = g2(0)
        assert checker(res)

        l2 = [(-value, -value), (value, value)]
        def h(x):
            return l2[x][1]
        h3 = compile(h, [int])
        res = h3(1)
        assert checker(res)

def _gamma(x):
    if rfloat.isnan(x) or (rfloat.isinf(x) and x > 0.):
        return x
    if rfloat.isinf(x):
        raise ValueError("math domain error")
    if x == 0.:
        raise ValueError("math domain error")
    if x == math.floor(x):
        if x < 0.:
            raise ValueError("math domain error")
        if x < len(_gamma_integrals):
            return _gamma_integrals[int(x) - 1]
    absx = abs(x)
    if absx < 1e-20:
        r = 1. / x
        if rfloat.isinf(r):
            raise OverflowError("math range error")
        return r
    if absx > 200.:
        if x < 0.:
            return 0. / -_sinpi(x)
        else:
            raise OverflowError("math range error")
    y = absx + _lanczos_g_minus_half
    if absx > _lanczos_g_minus_half:
        q = y - absx
        z = q - _lanczos_g_minus_half
    else:
        q = y - _lanczos_g_minus_half
        z = q - absx
    z = z * _lanczos_g / y
    if x < 0.:
        r = -math.pi / _sinpi(absx) / absx * math.exp(y) / _lanczos_sum(absx)
        r -= z * r
        if absx < 140.:
            r /= math.pow(y, absx - .5)
        else:
            sqrtpow = math.pow(y, absx / 2. - .25)
            r /= sqrtpow
            r /= sqrtpow
    else:
        r = _lanczos_sum(absx) / math.exp(y)
        r += z * r
        if absx < 140.:
            r *= math.pow(y, absx - .5)
        else:
            sqrtpow = math.pow(y, absx / 2. - .25)
            r *= sqrtpow
            r *= sqrtpow
    if rfloat.isinf(r):
        raise OverflowError("math range error")
    return r

 def format_float(self, w_value, char):
     space = self.space
     x = space.float_w(maybe_float(space, w_value))
     if isnan(x):
         if char in 'EFG':
             r = 'NAN'
         else:
             r = 'nan'
     elif isinf(x):
         if x < 0:
             if char in 'EFG':
                 r = '-INF'
             else:
                 r = '-inf'
         else:
             if char in 'EFG':
                 r = 'INF'
             else:
                 r = 'inf'
     else:
         prec = self.prec
         if prec < 0:
             prec = 6
         if char in 'fF' and x/1e25 > 1e25:
             char = chr(ord(char) + 1)     # 'f' => 'g'
         flags = 0
         if self.f_alt:
             flags |= DTSF_ALT
         r = formatd(x, char, prec, flags)
     self.std_wp_number(r)

 def push_primitive_constant(self, TYPE, value):
     ilasm = self.ilasm
     if TYPE is ootype.Void:
         pass
     elif TYPE is ootype.Bool:
         ilasm.opcode("ldc.i4", str(int(value)))
     elif TYPE is ootype.Char or TYPE is ootype.UniChar:
         ilasm.opcode("ldc.i4", ord(value))
     elif TYPE is ootype.Float:
         if isinf(value):
             if value < 0.0:
                 ilasm.opcode("ldc.r8", "(00 00 00 00 00 00 f0 ff)")
             else:
                 ilasm.opcode("ldc.r8", "(00 00 00 00 00 00 f0 7f)")
         elif isnan(value):
             ilasm.opcode("ldc.r8", "(00 00 00 00 00 00 f8 ff)")
         else:
             ilasm.opcode("ldc.r8", repr(value))
     elif isinstance(value, CDefinedIntSymbolic):
         ilasm.opcode("ldc.i4", DEFINED_INT_SYMBOLICS[value.expr])
     elif TYPE in (ootype.Signed, ootype.Unsigned, rffi.SHORT):
         ilasm.opcode("ldc.i4", str(value))
     elif TYPE in (ootype.SignedLongLong, ootype.UnsignedLongLong):
         ilasm.opcode("ldc.i8", str(value))
     elif TYPE in (ootype.String, ootype.Unicode):
         if value._str is None:
             ilasm.opcode("ldnull")
         else:
             ilasm.opcode("ldstr", string_literal(value._str))
     else:
         assert False, "Unexpected constant type"

def ll_math_fmod(x, y):
    if isinf(y):
        if isinf(x):
            raise ValueError("math domain error")
        return x  # fmod(x, +/-Inf) returns x for finite x (or if x is a NaN).

    _error_reset()
    r = math_fmod(x, y)
    errno = rposix.get_errno()
    if isnan(r):
        if isnan(x) or isnan(y):
            errno = 0
        else:
            errno = EDOM
    if errno:
        _likely_raise(errno, r)
    return r

def isfinitejs(ctx, args, this):
    if len(args) < 1:
        return newbool(True)
    n = args[0].ToNumber(ctx)
    if  isinf(n) or isnan(n):
        return newbool(False)
    else:
        return newbool(True)

def c_abs(x, y):
    if not isfinite(x) or not isfinite(y):
        # C99 rules: if either the real or the imaginary part is an
        # infinity, return infinity, even if the other part is a NaN.
        if isinf(x):
            return INF
        if isinf(y):
            return INF

        # either the real or imaginary part is a NaN,
        # and neither is infinite. Result should be NaN.
        return NAN

    result = math.hypot(x, y)
    if not isfinite(result):
        raise OverflowError("math range error")
    return result

def format_float(x, code, precision):
    # like float2string, except that the ".0" is not necessary
    if isinf(x):
        if x > 0.0:
            return "inf"
        else:
            return "-inf"
    elif isnan(x):
        return "nan"
    else:
        return formatd(x, code, precision)

def float2string(x, code, precision):
    # we special-case explicitly inf and nan here
    if isfinite(x):
        s = formatd(x, code, precision, DTSF_ADD_DOT_0)
    elif isinf(x):
        if x > 0.0:
            s = "inf"
        else:
            s = "-inf"
    else:  # isnan(x):
        s = "nan"
    return s

def c_tanh(x, y):
    # Formula:
    #
    #   tanh(x+iy) = (tanh(x)(1+tan(y)^2) + i tan(y)(1-tanh(x))^2) /
    #   (1+tan(y)^2 tanh(x)^2)
    #
    #   To avoid excessive roundoff error, 1-tanh(x)^2 is better computed
    #   as 1/cosh(x)^2.  When abs(x) is large, we approximate 1-tanh(x)^2
    #   by 4 exp(-2*x) instead, to avoid possible overflow in the
    #   computation of cosh(x).

    if not isfinite(x) or not isfinite(y):
        if isinf(x) and isfinite(y) and y != 0.:
            if x > 0:
                real = 1.0        # vv XXX why is the 2. there?
                imag = copysign(0., 2. * math.sin(y) * math.cos(y))
            else:
                real = -1.0
                imag = copysign(0., 2. * math.sin(y) * math.cos(y))
            r = (real, imag)
        else:
            r = tanh_special_values[special_type(x)][special_type(y)]

        # need to raise ValueError if y is +/-infinity and x is finite
        if isinf(y) and isfinite(x):
            raise ValueError("math domain error")
        return r

    if fabs(x) > CM_LOG_LARGE_DOUBLE:
        real = copysign(1., x)
        imag = 4. * math.sin(y) * math.cos(y) * math.exp(-2.*fabs(x))
    else:
        tx = math.tanh(x)
        ty = math.tan(y)
        cx = 1. / math.cosh(x)
        txty = tx * ty
        denom = 1. + txty * txty
        real = tx * (1. + ty*ty) / denom
        imag = ((ty / denom) * cx) * cx
    return real, imag

def _lgamma(x):
    if rfloat.isnan(x):
        return x
    if rfloat.isinf(x):
        return rfloat.INFINITY
    if x == math.floor(x) and x <= 2.:
        if x <= 0.:
            raise ValueError("math range error")
        return 0.
    absx = abs(x)
    if absx < 1e-20:
        return -math.log(absx)
    if x > 0.:
        r = (math.log(_lanczos_sum(x)) - _lanczos_g + (x - .5) *
             (math.log(x + _lanczos_g - .5) - 1))
    else:
        r = (math.log(math.pi) - math.log(abs(_sinpi(absx))) - math.log(absx) -
             (math.log(_lanczos_sum(absx)) - _lanczos_g +
              (absx - .5) * (math.log(absx + _lanczos_g - .5) - 1)))
    if rfloat.isinf(r):
        raise OverflowError("math domain error")
    return r

 def _push_double_constant(self, value):
     if isnan(value):
         jvm.DOUBLENAN.load(self)
     elif isinf(value):
         if value > 0: jvm.DOUBLEPOSINF.load(self)
         else: jvm.DOUBLENEGINF.load(self)
     elif value == 0.0:
         self.emit(jvm.DCONST_0)
     elif value == 1.0:
         self.emit(jvm.DCONST_1)
     else:
         # Big hack to avoid exponential notation:
         self.emit(jvm.LDC2, "%22.22f" % value)

def fsum(space, w_iterable):
    """Sum an iterable of floats, trying to keep precision."""
    w_iter = space.iter(w_iterable)
    inf_sum = special_sum = 0.0
    partials = []
    while True:
        try:
            w_value = space.next(w_iter)
        except OperationError, e:
            if not e.match(space, space.w_StopIteration):
                raise
            break
        v = _get_double(space, w_value)
        original = v
        added = 0
        for y in partials:
            if abs(v) < abs(y):
                v, y = y, v
            hi = v + y
            yr = hi - v
            lo = y - yr
            if lo != 0.0:
                partials[added] = lo
                added += 1
            v = hi
        del partials[added:]
        if v != 0.0:
            if rfloat.isinf(v) or rfloat.isnan(v):
                if (not rfloat.isinf(original) and
                    not rfloat.isnan(original)):
                    raise OperationError(space.w_OverflowError,
                                         space.wrap("intermediate overflow"))
                if rfloat.isinf(original):
                    inf_sum += original
                special_sum += original
                del partials[:]
            else:
                partials.append(v)

def c_phase(x, y):
    # Windows screws up atan2 for inf and nan, and alpha Tru64 5.1 doesn't
    # follow C99 for atan2(0., 0.).
    if isnan(x) or isnan(y):
        return NAN
    if isinf(y):
        if isinf(x):
            if copysign(1., x) == 1.:
                # atan2(+-inf, +inf) == +-pi/4
                return copysign(0.25 * math.pi, y)
            else:
                # atan2(+-inf, -inf) == +-pi*3/4
                return copysign(0.75 * math.pi, y)
        # atan2(+-inf, x) == +-pi/2 for finite x
        return copysign(0.5 * math.pi, y)
    if isinf(x) or y == 0.:
        if copysign(1., x) == 1.:
            # atan2(+-y, +inf) = atan2(+-0, +x) = +-0.
            return copysign(0., y)
        else:
            # atan2(+-y, -inf) = atan2(+-0., -x) = +-pi.
            return copysign(math.pi, y)
    return math.atan2(y, x)