def eval(cls, arg):
        from sympy import asin
        arg = sympify(arg)

        if arg.is_Number:
            if arg is S.NaN:
                return S.NaN
            elif arg is S.Infinity:
                return S.Infinity
            elif arg is S.NegativeInfinity:
                return S.NegativeInfinity
            elif arg is S.Zero:
                return S.Zero
            elif arg is S.One:
                return log(sqrt(2) + 1)
            elif arg is S.NegativeOne:
                return log(sqrt(2) - 1)
            elif arg.is_negative:
                return -cls(-arg)
        else:
            if arg is S.ComplexInfinity:
                return S.ComplexInfinity

            i_coeff = arg.as_coefficient(S.ImaginaryUnit)

            if i_coeff is not None:
                return S.ImaginaryUnit * asin(i_coeff)
            else:
                if _coeff_isneg(arg):
                    return -cls(-arg)

def get_math_macros():
    """ Returns a dictionary with math-related macros from math.h/cmath

    Note that these macros are not strictly required by the C/C++-standard.
    For MSVC they are enabled by defining "_USE_MATH_DEFINES" (preferably
    via a compilation flag).

    Returns
    =======

    Dictionary mapping sympy expressions to strings (macro names)

    """
    from sympy.codegen.cfunctions import log2, Sqrt
    from sympy.functions.elementary.exponential import log
    from sympy.functions.elementary.miscellaneous import sqrt

    return {
        S.Exp1: 'M_E',
        log2(S.Exp1): 'M_LOG2E',
        1/log(2): 'M_LOG2E',
        log(2): 'M_LN2',
        log(10): 'M_LN10',
        S.Pi: 'M_PI',
        S.Pi/2: 'M_PI_2',
        S.Pi/4: 'M_PI_4',
        1/S.Pi: 'M_1_PI',
        2/S.Pi: 'M_2_PI',
        2/sqrt(S.Pi): 'M_2_SQRTPI',
        2/Sqrt(S.Pi): 'M_2_SQRTPI',
        sqrt(2): 'M_SQRT2',
        Sqrt(2): 'M_SQRT2',
        1/sqrt(2): 'M_SQRT1_2',
        1/Sqrt(2): 'M_SQRT1_2'
    }

def swinnerton_dyer_poly(n, x=None, **args):
    """Generates n-th Swinnerton-Dyer polynomial in `x`.  """
    if n <= 0:
        raise ValueError(
            "can't generate Swinnerton-Dyer polynomial of order %s" % n)

    if x is not None:
        x, cls = sympify(x), Poly
    else:
        x, cls = Dummy('x'), PurePoly

    p, elts = 2, [[x, -sqrt(2)],
                  [x, sqrt(2)]]

    for i in xrange(2, n + 1):
        p, _elts = nextprime(p), []

        neg_sqrt = -sqrt(p)
        pos_sqrt = +sqrt(p)

        for elt in elts:
            _elts.append(elt + [neg_sqrt])
            _elts.append(elt + [pos_sqrt])

        elts = _elts

    poly = []

    for elt in elts:
        poly.append(Add(*elt))

    if not args.get('polys', False):
        return Mul(*poly).expand()
    else:
        return PurePoly(Mul(*poly), x)

def test_round():
    from sympy.abc import x

    assert round(Float('0.1249999'), 2) == 0.12
    d20 = 12345678901234567890
    ans = round(S(d20), 2)
    assert ans.is_Float and ans == d20
    ans = round(S(d20), -2)
    assert ans.is_Float and ans == 12345678901234567900
    assert round(S('1/7'), 4) == 0.1429
    assert round(S('.[12345]'), 4) == 0.1235
    assert round(S('.1349'), 2) == 0.13
    n = S(12345)
    ans = round(n)
    assert ans.is_Float
    assert ans == n
    ans = round(n, 1)
    assert ans.is_Float
    assert ans == n
    ans = round(n, 4)
    assert ans.is_Float
    assert ans == n
    assert round(n, -1) == 12350

    r = round(n, -4)
    assert r == 10000
    # in fact, it should equal many values since __eq__
    # compares at equal precision
    assert all(r == i for i in range(9984, 10049))

    assert round(n, -5) == 0

    assert round(pi + sqrt(2), 2) == 4.56
    assert round(10*(pi + sqrt(2)), -1) == 50
    raises(TypeError, 'round(x + 2, 2)')
    assert round(S(2.3), 1) == 2.3
    e = round(12.345, 2)
    assert e == _pyround(12.345, 2)
    assert type(e) is float

    assert round(Float(.3, 3) + 2*pi) == 7
    assert round(Float(.3, 3) + 2*pi*100) == 629
    assert round(Float(.03, 3) + 2*pi/100, 5) == 0.09283
    assert round(Float(.03, 3) + 2*pi/100, 4) == 0.0928

    assert round(S.Zero) == 0

    a = (Add(1, Float('1.'+'9'*27, ''), evaluate=0))
    assert round(a, 10) == Float('3.0000000000','')
    assert round(a, 25) == Float('3.0000000000000000000000000','')
    assert round(a, 26) == Float('3.00000000000000000000000000','')
    assert round(a, 27) == Float('2.999999999999999999999999999','')
    assert round(a, 30) == Float('2.999999999999999999999999999','')

    raises(TypeError, 'round(x)')
    raises(TypeError, 'round(1 + 3*I)')

    # exact magnitude of 10
    assert str(round(S(1))) == '1.'
    assert str(round(S(100))) == '100.'

def test_issue_6990():
    x = Symbol('x')
    a = Symbol('a')
    b = Symbol('b')
    assert (sqrt(a + b*x + x**2)).series(x, 0, 3).removeO() == \
        b*x/(2*sqrt(a)) + x**2*(1/(2*sqrt(a)) - \
        b**2/(8*a**(S(3)/2))) + sqrt(a)

    def intersection(self, o):
        if isinstance(o, Circle):
            if o.center == self.center:
                if o.radius == self.radius:
                    return o
                return []
            dx,dy = o.center - self.center
            d = sqrt( simplify(dy**2 + dx**2) )
            R = o.radius + self.radius
            if d > R or d < abs(self.radius - o.radius):
                return []

            a = simplify((self.radius**2 - o.radius**2 + d**2) / (2*d))

            x2 = self.center[0] + (dx * a/d)
            y2 = self.center[1] + (dy * a/d)

            h = sqrt( simplify(self.radius**2 - a**2) )
            rx = -dy * (h/d)
            ry =  dx * (h/d)

            xi_1 = simplify(x2 + rx)
            xi_2 = simplify(x2 - rx)
            yi_1 = simplify(y2 + ry)
            yi_2 = simplify(y2 - ry)

            ret = [Point(xi_1, yi_1)]
            if xi_1 != xi_2 or yi_1 != yi_2:
                ret.append(Point(xi_2, yi_2))
            return ret

        return Ellipse.intersection(self, o)

    def eval(cls, arg):
        arg = sympify(arg)

        if arg.is_Number:
            if arg is S.NaN:
                return S.NaN
            elif arg is S.Infinity:
                return S.One
            elif arg is S.NegativeInfinity:
                return S.NegativeOne
            elif arg is S.Zero:
                return S.Zero
            elif arg.is_negative:
                return -cls(-arg)
        else:
            i_coeff = arg.as_coefficient(S.ImaginaryUnit)

            if i_coeff is not None:
                return S.ImaginaryUnit * C.tan(i_coeff)
            else:
                coeff, terms = arg.as_coeff_terms()

                if coeff.is_negative:
                    return -cls(-arg)

            if isinstance(arg, asinh):
                x = arg.args[0]
                return x/sqrt(1+x**2)

            if isinstance(arg, acosh):
                x = arg.args[0]
                return sqrt(x-1) * sqrt(x+1) / x

            if isinstance(arg, atanh):
                return arg.args[0]

def test_issue_3891():
    x = Symbol("x")
    a = Symbol("a")
    b = Symbol("b")
    assert (sqrt(a + b * x + x ** 2)).series(x, 0, 3).removeO() == b * x / (2 * sqrt(a)) + x ** 2 * (
        1 / (2 * sqrt(a)) - b ** 2 / (8 * a ** (S(3) / 2))
    ) + sqrt(a)

    def eval(cls, arg):
        arg = sympify(arg)

        if arg.is_Number:
            if arg is S.NaN:
                return S.NaN
            elif arg is S.Infinity:
                return S.Infinity
            elif arg is S.NegativeInfinity:
                return S.NegativeInfinity
            elif arg is S.Zero:
                return S.Zero
            elif arg.is_negative:
                return -cls(-arg)
        else:
            i_coeff = arg.as_coefficient(S.ImaginaryUnit)

            if i_coeff is not None:
                return S.ImaginaryUnit * C.sin(i_coeff)
            else:
                coeff, terms = arg.as_coeff_terms()

                if coeff.is_negative:
                    return -cls(-arg)

            if arg.func == asinh:
                return arg.args[0]

            if arg.func == acosh:
                x = arg.args[0]
                return sqrt(x - 1) * sqrt(x + 1)

            if arg.func == atanh:
                x = arg.args[0]
                return x / sqrt(1 - x ** 2)

def test_root():
    from sympy.abc import x
    n = Symbol('n', integer=True)
    k = Symbol('k', integer=True)

    assert root(2, 2) == sqrt(2)
    assert root(2, 1) == 2
    assert root(2, 3) == 2**Rational(1, 3)
    assert root(2, 3) == cbrt(2)
    assert root(2, -5) == 2**Rational(4, 5)/2

    assert root(-2, 1) == -2

    assert root(-2, 2) == sqrt(2)*I
    assert root(-2, 1) == -2

    assert root(x, 2) == sqrt(x)
    assert root(x, 1) == x
    assert root(x, 3) == x**Rational(1, 3)
    assert root(x, 3) == cbrt(x)
    assert root(x, -5) == x**Rational(-1, 5)

    assert root(x, n) == x**(1/n)
    assert root(x, -n) == x**(-1/n)

    assert root(x, n, k) == x**(1/n)*(-1)**(2*k/n)

def test_power_rewrite_exp():
    assert (I**I).rewrite(exp) == exp(-pi/2)

    expr = (2 + 3*I)**(4 + 5*I)
    assert expr.rewrite(exp) == exp((4 + 5*I)*(log(sqrt(13)) + I*atan(S(3)/2)))
    assert expr.rewrite(exp).expand() == \
        169*exp(5*I*log(13)/2)*exp(4*I*atan(S(3)/2))*exp(-5*atan(S(3)/2))

    assert ((6 + 7*I)**5).rewrite(exp) == 7225*sqrt(85)*exp(5*I*atan(S(7)/6))

    expr = 5**(6 + 7*I)
    assert expr.rewrite(exp) == exp((6 + 7*I)*log(5))
    assert expr.rewrite(exp).expand() == 15625*exp(7*I*log(5))

    assert Pow(123, 789, evaluate=False).rewrite(exp) == 123**789
    assert (1**I).rewrite(exp) == 1**I
    assert (0**I).rewrite(exp) == 0**I

    expr = (-2)**(2 + 5*I)
    assert expr.rewrite(exp) == exp((2 + 5*I)*(log(2) + I*pi))
    assert expr.rewrite(exp).expand() == 4*exp(-5*pi)*exp(5*I*log(2))

    assert ((-2)**S(-5)).rewrite(exp) == (-2)**S(-5)

    x, y = symbols('x y')
    assert (x**y).rewrite(exp) == exp(y*log(x))
    assert (7**x).rewrite(exp) == exp(x*log(7), evaluate=False)
    assert ((2 + 3*I)**x).rewrite(exp) == exp(x*(log(sqrt(13)) + I*atan(S(3)/2)))
    assert (y**(5 + 6*I)).rewrite(exp) == exp(log(y)*(5 + 6*I))

    assert all((1/func(x)).rewrite(exp) == 1/(func(x).rewrite(exp)) for func in
                    (sin, cos, tan, sec, csc, sinh, cosh, tanh))

    def __new__(
        cls, center=None, hradius=None, vradius=None, eccentricity=None,
            **kwargs):
        hradius = sympify(hradius)
        vradius = sympify(vradius)

        eccentricity = sympify(eccentricity)

        if center is None:
            center = Point(0, 0)
        else:
            center = Point(center, dim=2)

        if len(center) != 2:
            raise ValueError('The center of "{0}" must be a two dimensional point'.format(cls))

        if len(list(filter(None, (hradius, vradius, eccentricity)))) != 2:
            raise ValueError('Exactly two arguments of "hradius", '
                '"vradius", and "eccentricity" must not be None."')

        if eccentricity is not None:
            if hradius is None:
                hradius = vradius / sqrt(1 - eccentricity**2)
            elif vradius is None:
                vradius = hradius * sqrt(1 - eccentricity**2)

        if hradius == vradius:
            return Circle(center, hradius, **kwargs)

        return GeometryEntity.__new__(cls, center, hradius, vradius, **kwargs)

    def __new__(cls, center=None, hradius=None, vradius=None, eccentricity=None,
                **kwargs):
        hradius = sympify(hradius)
        vradius = sympify(vradius)
        eccentricity = sympify(eccentricity)

        if len(filter(None, (hradius, vradius, eccentricity))) != 2:
            raise ValueError, 'Exactly two arguments between "hradius", '\
                '"vradius", and "eccentricity" must be not None."'

        if eccentricity is not None:
            if hradius is None:
                hradius = vradius / sqrt(1 - eccentricity**2)
            elif vradius is None:
                vradius = hradius * sqrt(1 - eccentricity**2)
        else:
            if hradius is None and vradius is None:
                raise ValueError("At least two arguments between hradius, "
                    "vradius and eccentricity must not be none.")

        if center is None:
            center = Point(0, 0)

        if not isinstance(center, Point):
            raise TypeError("center must be a Point")

        if hradius == vradius:
            return Circle(center, hradius, **kwargs)
        return GeometryEntity.__new__(cls, center, hradius, vradius, **kwargs)

 def _cholesky(self, hermitian=True):
     """Helper function of cholesky.
     Without the error checks.
     To be used privately.
     Implements the Cholesky-Banachiewicz algorithm.
     Returns L such that L*L.H == self if hermitian flag is True,
     or L*L.T == self if hermitian is False.
     """
     L = zeros(self.rows, self.rows)
     if hermitian:
         for i in range(self.rows):
             for j in range(i):
                 L[i, j] = (1 / L[j, j])*expand_mul(self[i, j] -
                     sum(L[i, k]*L[j, k].conjugate() for k in range(j)))
             Lii2 = expand_mul(self[i, i] -
                 sum(L[i, k]*L[i, k].conjugate() for k in range(i)))
             if Lii2.is_positive is False:
                 raise ValueError("Matrix must be positive-definite")
             L[i, i] = sqrt(Lii2)
     else:
         for i in range(self.rows):
             for j in range(i):
                 L[i, j] = (1 / L[j, j])*(self[i, j] -
                     sum(L[i, k]*L[j, k] for k in range(j)))
             L[i, i] = sqrt(self[i, i] -
                 sum(L[i, k]**2 for k in range(i)))
     return self._new(L)

def swinnerton_dyer_poly(n, x=None, **args):
    """Generates n-th Swinnerton-Dyer polynomial in `x`.  """
    from numberfields import minimal_polynomial
    if n <= 0:
        raise ValueError(
            "can't generate Swinnerton-Dyer polynomial of order %s" % n)

    if x is not None:
        sympify(x)
    else:
        x = Dummy('x')

    if n > 3:
        p = 2
        a = [sqrt(2)]
        for i in xrange(2, n + 1):
            p = nextprime(p)
            a.append(sqrt(p))
        return minimal_polynomial(Add(*a), x, polys=args.get('polys', False))

    if n == 1:
        ex = x**2 - 2
    elif n == 2:
        ex = x**4 - 10*x**2 + 1
    elif n == 3:
        ex = x**8 - 40*x**6 + 352*x**4 - 960*x**2 + 576
    if not args.get('polys', False):
        return ex
    else:
        return PurePoly(ex, x)

 def eval(cls, n, m, z=None):
     if z is not None:
         n, z, m = n, m, z
         k = 2 * z / pi
         if n == S.Zero:
             return elliptic_f(z, m)
         elif n == S.One:
             return elliptic_f(z, m) + (sqrt(1 - m * sin(z) ** 2) * tan(z) - elliptic_e(z, m)) / (1 - m)
         elif k.is_integer:
             return k * elliptic_pi(n, m)
         elif m == S.Zero:
             return atanh(sqrt(n - 1) * tan(z)) / sqrt(n - 1)
         elif n == m:
             return elliptic_f(z, n) - elliptic_pi(1, z, n) + tan(z) / sqrt(1 - n * sin(z) ** 2)
         elif n in (S.Infinity, S.NegativeInfinity):
             return S.Zero
         elif m in (S.Infinity, S.NegativeInfinity):
             return S.Zero
         elif z.could_extract_minus_sign():
             return -elliptic_pi(n, -z, m)
     else:
         if n == S.Zero:
             return elliptic_k(m)
         elif n == S.One:
             return S.ComplexInfinity
         elif m == S.Zero:
             return pi / (2 * sqrt(1 - n))
         elif m == S.One:
             return -S.Infinity / sign(n - 1)
         elif n == m:
             return elliptic_e(n) / (1 - n)
         elif n in (S.Infinity, S.NegativeInfinity):
             return S.Zero
         elif m in (S.Infinity, S.NegativeInfinity):
             return S.Zero

    def intersection(self, o):
        if isinstance(o, Circle):
            dx,dy = o.center - self.center
            d = sqrt( simplify(dy**2 + dx**2) )
            a = simplify((self.radius**2 - o.radius**2 + d**2) / (2*d))

            x2 = self.center[0] + (dx * a/d)
            y2 = self.center[1] + (dy * a/d)

            h = sqrt( simplify(self.radius**2 - a**2) )
            rx = -dy * (h/d)
            ry =  dx * (h/d)

            xi_1 = simplify(x2 + rx)
            xi_2 = simplify(x2 - rx)
            yi_1 = simplify(y2 + ry)
            yi_2 = simplify(y2 - ry)

            ret = [Point(xi_1, yi_1)]
            if xi_1 != xi_2 or yi_1 != yi_2:
                ret.append(Point(xi_2, yi_2))
            return ret
        elif isinstance(o, Ellipse):
            a, b, r = o.hradius, o.vradius, self.radius
            x = a*sqrt(simplify((r**2 - b**2)/(a**2 - b**2)))
            y = b*sqrt(simplify((a**2 - r**2)/(a**2 - b**2)))
            return list(set([Point(x,y), Point(x,-y), Point(-x,y), Point(-x,-y)]))

        return Ellipse.intersection(self, o)

    def eval(cls, arg):
        arg = sympify(arg)

        if arg.is_Number:
            if arg is S.NaN:
                return S.NaN
            elif arg is S.Infinity:
                return S.Infinity
            elif arg is S.NegativeInfinity:
                return S.Infinity
            elif arg is S.Zero:
                return S.One
            elif arg.is_negative:
                return cls(-arg)
        else:
            i_coeff = arg.as_coefficient(S.ImaginaryUnit)

            if i_coeff is not None:
                return C.cos(i_coeff)
            else:
                if arg.as_coeff_mul()[0].is_negative:
                    return cls(-arg)

            if arg.func == asinh:
                return sqrt(1+arg.args[0]**2)

            if arg.func == acosh:
                return arg.args[0]

            if arg.func == atanh:
                return 1/sqrt(1-arg.args[0]**2)

            if arg.func == acoth:
                x = arg.args[0]
                return x/(sqrt(x-1) * sqrt(x+1))

    def eval(cls, arg):
        if arg.is_Number:
            if arg is S.NaN:
                return S.NaN
            elif arg is S.Infinity:
                return S.Infinity
            elif arg.is_Integer:
                if arg.is_positive:
                    return C.factorial(arg - 1)
                else:
                    return S.ComplexInfinity
            elif arg.is_Rational:
                if arg.q == 2:
                    n = abs(arg.p) // arg.q

                    if arg.is_positive:
                        k, coeff = n, S.One
                    else:
                        n = k = n + 1

                        if n & 1 == 0:
                            coeff = S.One
                        else:
                            coeff = S.NegativeOne

                    for i in range(3, 2 * k, 2):
                        coeff *= i

                    if arg.is_positive:
                        return coeff * sqrt(S.Pi) / 2 ** n
                    else:
                        return 2 ** n * sqrt(S.Pi) / coeff

def test_issue350():
    #test if powers are simplified correctly
    #see also issue 896
    x = Symbol('x')
    assert ((x**Rational(1,3))**Rational(2)) == x**Rational(2,3)
    assert ((x**Rational(3))**Rational(2,5)) == (x**Rational(3))**Rational(2,5)

    a = Symbol('a', real=True)
    b = Symbol('b', real=True)
    assert (a**2)**b == (abs(a)**b)**2
    assert sqrt(1/a) != 1/sqrt(a) # e.g. for a = -1
    assert (a**3)**Rational(1,3) != a
    assert (x**a)**b != x**(a*b) # e.g. x = -1, a=2, b=1/2
    assert (x**.5)**b == x**(.5*b)
    assert (x**.5)**.5 == x**.25
    assert (x**2.5)**.5 != x**1.25 # e.g. for x = 5*I

    k = Symbol('k',integer=True)
    m = Symbol('m',integer=True)
    assert (x**k)**m == x**(k*m)
    assert Number(5)**Rational(2,3)==Number(25)**Rational(1,3)

    assert (x**.5)**2 == x**1.0
    assert (x**2)**k == (x**k)**2 == x**(2*k)

    a = Symbol('a', positive=True)
    assert (a**3)**Rational(2,5) == a**Rational(6,5)
    assert (a**2)**b == (a**b)**2
    assert (a**Rational(2, 3))**x == (a**(2*x/3)) != (a**x)**Rational(2, 3)