def test_issue_841():
    a, b, c, d = symbols('a:d', positive=True, bounded=True)
    assert integrate(
        exp(-x**2 + I*c*x), x) == sqrt(pi)*erf(x - I*c/2)*exp(-c**S(2)/4)/2
    assert integrate(exp(a*x**2 + b*x + c), x) == \
          I*sqrt(pi)*erf(-I*x*sqrt(
              a) - I*b/(2*sqrt(a)))*exp(c)*exp(-b**2/(4*a))/(2*sqrt(a))

def test_series():
    x = symbols('x')
    R, Dx = DifferentialOperators(ZZ.old_poly_ring(x), 'Dx')
    p = HolonomicFunction(Dx**2 + 2*x*Dx, x, 0, [0, 1]).series(n=10)
    q = x - x**3/3 + x**5/10 - x**7/42 + x**9/216 + O(x**10)
    assert p == q
    p = HolonomicFunction(Dx - 1, x).composition(x**2, 0, [1])  # e^(x**2)
    q = HolonomicFunction(Dx**2 + 1, x, 0, [1, 0])  # cos(x)
    r = (p * q).series(n=10)  # expansion of cos(x) * exp(x**2)
    s = 1 + x**2/2 + x**4/24 - 31*x**6/720 - 179*x**8/8064 + O(x**10)
    assert r == s
    t = HolonomicFunction((1 + x)*Dx**2 + Dx, x, 0, [0, 1])  # log(1 + x)
    r = (p * t + q).series(n=10)
    s = 1 + x - x**2 + 4*x**3/3 - 17*x**4/24 + 31*x**5/30 - 481*x**6/720 +\
     71*x**7/105 - 20159*x**8/40320 + 379*x**9/840 + O(x**10)
    assert r == s
    p = HolonomicFunction((6+6*x-3*x**2) - (10*x-3*x**2-3*x**3)*Dx + \
        (4-6*x**3+2*x**4)*Dx**2, x, 0, [0, 1]).series(n=7)
    q = x + x**3/6 - 3*x**4/16 + x**5/20 - 23*x**6/960 + O(x**7)
    assert p == q
    p = HolonomicFunction((6+6*x-3*x**2) - (10*x-3*x**2-3*x**3)*Dx + \
        (4-6*x**3+2*x**4)*Dx**2, x, 0, [1, 0]).series(n=7)
    q = 1 - 3*x**2/4 - x**3/4 - 5*x**4/32 - 3*x**5/40 - 17*x**6/384 + O(x**7)
    assert p == q
    p = expr_to_holonomic(erf(x) + x).series(n=10)
    C_3 = symbols('C_3')
    q = (erf(x) + x).series(n=10)
    assert p.subs(C_3, -2/(3*sqrt(pi))) == q
    assert expr_to_holonomic(sqrt(x**3 + x)).series(n=10) == sqrt(x**3 + x).series(n=10)
    assert expr_to_holonomic((2*x - 3*x**2)**(S(1)/3)).series() == ((2*x - 3*x**2)**(S(1)/3)).series()
    assert  expr_to_holonomic(sqrt(x**2-x)).series() == (sqrt(x**2-x)).series()
    assert expr_to_holonomic(cos(x)**2/x**2, y0={-2: [1, 0, -1]}).series(n=10) == (cos(x)**2/x**2).series(n=10)
    assert expr_to_holonomic(cos(x)**2/x**2, x0=1).series(n=10) == (cos(x)**2/x**2).series(n=10, x0=1)
    assert expr_to_holonomic(cos(x-1)**2/(x-1)**2, x0=1, y0={-2: [1, 0, -1]}).series(n=10) \
        == (cos(x-1)**2/(x-1)**2).series(x0=1, n=10)

def test__erfs():
    assert _erfs(z).diff(z) == -2/sqrt(S.Pi)+2*z*_erfs(z)

    assert _erfs(1/z).series(z) == z/sqrt(pi) - z**3/(2*sqrt(pi)) + 3*z**5/(4*sqrt(pi)) + O(z**6)

    assert expand(erf(z).rewrite('tractable').diff(z).rewrite('intractable')) == erf(z).diff(z)
    assert _erfs(z).rewrite("intractable") == (-erf(z) + 1)*exp(z**2)

def test_to_expr():
    x = symbols('x')
    R, Dx = DifferentialOperators(ZZ.old_poly_ring(x), 'Dx')
    p = HolonomicFunction(Dx - 1, x, 0, 1).to_expr()
    q = exp(x)
    assert p == q
    p = HolonomicFunction(Dx**2 + 1, x, 0, [1, 0]).to_expr()
    q = cos(x)
    assert p == q
    p = HolonomicFunction(Dx**2 - 1, x, 0, [1, 0]).to_expr()
    q = cosh(x)
    assert p == q
    p = HolonomicFunction(2 + (4*x - 1)*Dx + \
        (x**2 - x)*Dx**2, x, 0, [1, 2]).to_expr()
    q = 1/(x**2 - 2*x + 1)
    assert p == q
    p = expr_to_holonomic(sin(x)**2/x).integrate((x, 0, x)).to_expr()
    q = (sin(x)**2/x).integrate((x, 0, x))
    assert p == q
    C_1, C_2, C_3 = symbols('C_1, C_2, C_3')
    p = expr_to_holonomic(log(1+x**2)).to_expr()
    q = C_2*log(x**2 + 1)
    assert p == q
    p = expr_to_holonomic(log(1+x**2)).diff().to_expr()
    q = C_1*x/(x**2 + 1)
    assert p == q
    p = expr_to_holonomic(erf(x) + x).to_expr()
    q = 3*C_3*x - 3*sqrt(pi)*C_3*erf(x)/2 + x + 2*x/sqrt(pi)
    assert p == q
    p = expr_to_holonomic(sqrt(x), x0=1).to_expr()
    assert p == sqrt(x)
    p = expr_to_holonomic(sqrt(1 + x**2)).to_expr()
    assert p == sqrt(1+x**2)
    p = expr_to_holonomic((2*x**2 + 1)**(S(2)/3)).to_expr()
    assert p == (2*x**2 + 1)**(S(2)/3)

def test_issue_1791():
    z = Symbol("z", positive=True)
    assert integrate(exp(-log(x) ** 2), x) == pi ** (S(1) / 2) * erf(-S(1) / 2 + log(x)) * exp(S(1) / 4) / 2
    assert integrate(exp(log(x) ** 2), x) == -I * pi ** (S(1) / 2) * erf(I * log(x) + I / 2) * exp(-S(1) / 4) / 2
    assert integrate(exp(-z * log(x) ** 2), x) == pi ** (S(1) / 2) * erf(
        z ** (S(1) / 2) * log(x) - 1 / (2 * z ** (S(1) / 2))
    ) * exp(S(1) / (4 * z)) / (2 * z ** (S(1) / 2))

def test_to_expr():
    x = symbols('x')
    R, Dx = DifferentialOperators(ZZ.old_poly_ring(x), 'Dx')
    p = HolonomicFunction(Dx - 1, x, 0, [1]).to_expr()
    q = exp(x)
    assert p == q
    p = HolonomicFunction(Dx**2 + 1, x, 0, [1, 0]).to_expr()
    q = cos(x)
    assert p == q
    p = HolonomicFunction(Dx**2 - 1, x, 0, [1, 0]).to_expr()
    q = cosh(x)
    assert p == q
    p = HolonomicFunction(2 + (4*x - 1)*Dx + \
        (x**2 - x)*Dx**2, x, 0, [1, 2]).to_expr().expand()
    q = 1/(x**2 - 2*x + 1)
    assert p == q
    p = expr_to_holonomic(sin(x)**2/x).integrate((x, 0, x)).to_expr()
    q = (sin(x)**2/x).integrate((x, 0, x))
    assert p == q
    C_0, C_1, C_2, C_3 = symbols('C_0, C_1, C_2, C_3')
    p = expr_to_holonomic(log(1+x**2)).to_expr()
    q = C_2*log(x**2 + 1)
    assert p == q
    p = expr_to_holonomic(log(1+x**2)).diff().to_expr()
    q = C_0*x/(x**2 + 1)
    assert p == q
    p = expr_to_holonomic(erf(x) + x).to_expr()
    q = 3*C_3*x - 3*sqrt(pi)*C_3*erf(x)/2 + x + 2*x/sqrt(pi)
    assert p == q
    p = expr_to_holonomic(sqrt(x), x0=1).to_expr()
    assert p == sqrt(x)
    assert expr_to_holonomic(sqrt(x)).to_expr() == sqrt(x)
    p = expr_to_holonomic(sqrt(1 + x**2)).to_expr()
    assert p == sqrt(1+x**2)
    p = expr_to_holonomic((2*x**2 + 1)**(S(2)/3)).to_expr()
    assert p == (2*x**2 + 1)**(S(2)/3)
    p = expr_to_holonomic(sqrt(-x**2+2*x)).to_expr()
    assert p == sqrt(x)*sqrt(-x + 2)
    p = expr_to_holonomic((-2*x**3+7*x)**(S(2)/3)).to_expr()
    q = x**(S(2)/3)*(-2*x**2 + 7)**(S(2)/3)
    assert p == q
    p = from_hyper(hyper((-2, -3), (S(1)/2, ), x))
    s = hyperexpand(hyper((-2, -3), (S(1)/2, ), x))
    D_0 = Symbol('D_0')
    C_0 = Symbol('C_0')
    assert (p.to_expr().subs({C_0:1, D_0:0}) - s).simplify() == 0
    p.y0 = {0: [1], S(1)/2: [0]}
    assert p.to_expr() == s
    assert expr_to_holonomic(x**5).to_expr() == x**5
    assert expr_to_holonomic(2*x**3-3*x**2).to_expr().expand() == \
        2*x**3-3*x**2
    a = symbols("a")
    p = (expr_to_holonomic(1.4*x)*expr_to_holonomic(a*x, x)).to_expr()
    q = 1.4*a*x**2
    assert p == q
    p = (expr_to_holonomic(1.4*x)+expr_to_holonomic(a*x, x)).to_expr()
    q = x*(a + 1.4)
    assert p == q
    p = (expr_to_holonomic(1.4*x)+expr_to_holonomic(x)).to_expr()
    assert p == 2.4*x

def test_manualintegrate_special():
    f, F = 4*exp(-x**2/3), 2*sqrt(3)*sqrt(pi)*erf(sqrt(3)*x/3)
    assert manualintegrate(f, x) == F and F.diff(x).equals(f)
    f, F = 3*exp(4*x**2), 3*sqrt(pi)*erfi(2*x)/4
    assert manualintegrate(f, x) == F and F.diff(x).equals(f)
    f, F = x**(S(1)/3)*exp(-x/8), -16*uppergamma(S(4)/3, x/8)
    assert manualintegrate(f, x) == F and F.diff(x).equals(f)
    f, F = exp(2*x)/x, Ei(2*x)
    assert manualintegrate(f, x) == F and F.diff(x).equals(f)
    f, F = exp(1 + 2*x - x**2), sqrt(pi)*exp(2)*erf(x - 1)/2
    assert manualintegrate(f, x) == F and F.diff(x).equals(f)
    f = sin(x**2 + 4*x + 1)
    F = (sqrt(2)*sqrt(pi)*(-sin(3)*fresnelc(sqrt(2)*(2*x + 4)/(2*sqrt(pi))) +
        cos(3)*fresnels(sqrt(2)*(2*x + 4)/(2*sqrt(pi))))/2)
    assert manualintegrate(f, x) == F and F.diff(x).equals(f)
    f, F = cos(4*x**2), sqrt(2)*sqrt(pi)*fresnelc(2*sqrt(2)*x/sqrt(pi))/4
    assert manualintegrate(f, x) == F and F.diff(x).equals(f)
    f, F = sin(3*x + 2)/x, sin(2)*Ci(3*x) + cos(2)*Si(3*x)
    assert manualintegrate(f, x) == F and F.diff(x).equals(f)
    f, F = sinh(3*x - 2)/x, -sinh(2)*Chi(3*x) + cosh(2)*Shi(3*x)
    assert manualintegrate(f, x) == F and F.diff(x).equals(f)
    f, F = 5*cos(2*x - 3)/x, 5*cos(3)*Ci(2*x) + 5*sin(3)*Si(2*x)
    assert manualintegrate(f, x) == F and F.diff(x).equals(f)
    f, F = cosh(x/2)/x, Chi(x/2)
    assert manualintegrate(f, x) == F and F.diff(x).equals(f)
    f, F = cos(x**2)/x, Ci(x**2)/2
    assert manualintegrate(f, x) == F and F.diff(x).equals(f)
    f, F = 1/log(2*x + 1), li(2*x + 1)/2
    assert manualintegrate(f, x) == F and F.diff(x).equals(f)
    f, F = polylog(2, 5*x)/x, polylog(3, 5*x)
    assert manualintegrate(f, x) == F and F.diff(x).equals(f)
    f, F = 5/sqrt(3 - 2*sin(x)**2), 5*sqrt(3)*elliptic_f(x, S(2)/3)/3
    assert manualintegrate(f, x) == F and F.diff(x).equals(f)
    f, F = sqrt(4 + 9*sin(x)**2), 2*elliptic_e(x, -S(9)/4)
    assert manualintegrate(f, x) == F and F.diff(x).equals(f)

def test_simplify_other():
    assert simplify(sin(x) ** 2 + cos(x) ** 2) == 1
    assert simplify(gamma(x + 1) / gamma(x)) == x
    assert simplify(sin(x) ** 2 + cos(x) ** 2 + factorial(x) / gamma(x)) == 1 + x
    assert simplify(Eq(sin(x) ** 2 + cos(x) ** 2, factorial(x) / gamma(x))) == Eq(x, 1)
    nc = symbols("nc", commutative=False)
    assert simplify(x + x * nc) == x * (1 + nc)
    # issue 6123
    # f = exp(-I*(k*sqrt(t) + x/(2*sqrt(t)))**2)
    # ans = integrate(f, (k, -oo, oo), conds='none')
    ans = I * (
        -pi
        * x
        * exp(-3 * I * pi / 4 + I * x ** 2 / (4 * t))
        * erf(x * exp(-3 * I * pi / 4) / (2 * sqrt(t)))
        / (2 * sqrt(t))
        + pi * x * exp(-3 * I * pi / 4 + I * x ** 2 / (4 * t)) / (2 * sqrt(t))
    ) * exp(-I * x ** 2 / (4 * t)) / (sqrt(pi) * x) - I * sqrt(pi) * (
        -erf(x * exp(I * pi / 4) / (2 * sqrt(t))) + 1
    ) * exp(
        I * pi / 4
    ) / (
        2 * sqrt(t)
    )
    assert simplify(ans) == -(-1) ** (S(3) / 4) * sqrt(pi) / sqrt(t)
    # issue 6370
    assert simplify(2 ** (2 + x) / 4) == 2 ** x

def test_branch_bug():
    assert hyperexpand(hyper((-S(1)/3, S(1)/2), (S(2)/3, S(3)/2), -z)) == \
        -z**S('1/3')*lowergamma(exp_polar(I*pi)/3, z)/5 \
        + sqrt(pi)*erf(sqrt(z))/(5*sqrt(z))
    assert hyperexpand(meijerg([S(7)/6, 1], [], [S(2)/3], [S(1)/6, 0], z)) == \
        2*z**S('2/3')*(2*sqrt(pi)*erf(sqrt(z))/sqrt(z) - 2*lowergamma(
                       S(2)/3, z)/z**S('2/3'))*gamma(S(2)/3)/gamma(S(5)/3)

def test_branch_bug():
    from sympy import powdenest, lowergamma
    # TODO combsimp cannot prove that the factor is unity
    assert powdenest(integrate(erf(x**3), x, meijerg=True).diff(x),
           polar=True) == 2*erf(x**3)*gamma(S(2)/3)/3/gamma(S(5)/3)
    assert integrate(erf(x**3), x, meijerg=True) == \
        2*x*erf(x**3)*gamma(S(2)/3)/(3*gamma(S(5)/3)) \
        - 2*gamma(S(2)/3)*lowergamma(S(2)/3, x**6)/(3*sqrt(pi)*gamma(S(5)/3))

def test_issue_1791():
    z = Symbol('z', positive=True)
    assert integrate(exp(-log(x)**2), x) == \
        sqrt(pi)*erf(-S(1)/2 + log(x))*exp(S(1)/4)/2
    assert integrate(exp(log(x)**2), x) == \
        -I*sqrt(pi)*exp(-S(1)/4)*erf(I*log(x) + I/2)/2
    assert integrate(exp(-z*log(x)**2), x) == sqrt(pi)*erf(sqrt(z)*log(x)
        - 1/(2*sqrt(z)))*exp(S(1)/(4*z))/(2*sqrt(z))

def test_issue_4890():
    z = Symbol('z', positive=True)
    assert integrate(exp(-log(x)**2), x) == \
        sqrt(pi)*exp(S(1)/4)*erf(log(x)-S(1)/2)/2
    assert integrate(exp(log(x)**2), x) == \
        sqrt(pi)*exp(-S(1)/4)*erfi(log(x)+S(1)/2)/2
    assert integrate(exp(-z*log(x)**2), x) == \
        sqrt(pi)*exp(1/(4*z))*erf(sqrt(z)*log(x) - 1/(2*sqrt(z)))/(2*sqrt(z))

def test_errorinverses():
    assert solveset_real(erf(x) - S.One/2, x) == \
        FiniteSet(erfinv(S.One/2))
    assert solveset_real(erfinv(x) - 2, x) == \
        FiniteSet(erf(2))
    assert solveset_real(erfc(x) - S.One, x) == \
        FiniteSet(erfcinv(S.One))
    assert solveset_real(erfcinv(x) - 2, x) == FiniteSet(erfc(2))

 def fdiff(self, argindex=2):
     if argindex == 2:
         x0, x1 = self.args
         return -2 * exp(-x1 ** 2) / (sqrt(pi) * (erf(x0) - erf(x1)))
     elif argindex == 1:
         x0, x1 = self.args
         return 2.0 * exp(-x0 ** 2) / (sqrt(pi) * (erf(x0) - erf(x1)))
     else:
         raise ArgumentIndexError(self, argindex)

def test_issue_841():
    a,b,c,d = symbols('a:d', positive=True, bounded=True)
    assert integrate(exp(-x**2 + I*c*x), x) == sqrt(pi)*erf(x - I*c/2)*exp(-c**S(2)/4)/2
    assert integrate(exp(a*x**2 + b*x + c), x) == I*sqrt(pi)*erf(-I*x*sqrt(a) - I*b/(2*sqrt(a)))*exp(c)*exp(-b**2/(4*a))/(2*sqrt(a))
    a,b,c,d = symbols('a:d', positive=True)
    i = integrate(exp(-a*x**2 + 2*d*x), (x, -oo, oo))
    ans = sqrt(pi)*exp(d**2/a)*(1 + erf(oo - d/sqrt(a)))/(2*sqrt(a))
    n, d = i.as_numer_denom()
    assert factor(n, expand=False)/d == ans

def test_erfinv():
    assert erfinv(0) == 0
    assert erfinv(1) == S.Infinity
    assert erfinv(nan) == S.NaN

    assert erfinv(erf(w)) == w
    assert erfinv(erf(-w)) == -w

    assert erfinv(x).diff() == sqrt(pi)*exp(erfinv(x)**2)/2

    assert erfinv(z).rewrite('erfcinv') == erfcinv(1-z)

def test_erf():
    assert erf(nan) == nan

    assert erf(oo) == 1
    assert erf(-oo) == -1

    assert erf(0) == 0

    assert erf(-2) == -erf(2)
    assert erf(-x*y) == -erf(x*y)
    assert erf(-x - y) == -erf(x + y)

def test_issue_841():
    from sympy import simplify
    a = Symbol('a', positive = True)
    b = Symbol('b')
    c = Symbol('c')
    d = Symbol('d', positive = True)
    assert integrate(exp(-x**2 + I*c*x), x) == pi**(S(1)/2)*erf(x - I*c/2)*exp(-c**S(2)/4)/2
    assert integrate(exp(a*x**2 + b*x + c), x) == \
          I*pi**(S(1)/2)*erf(-I*x*a**(S(1)/2) - I*b/(2*a**(S(1)/2)))*exp(c)*exp(-b**2/(4*a))/(2*a**(S(1)/2))
    assert simplify(integrate(exp(-a*x**2 + 2*d*x), (x, -oo, oo))) == pi**(S(1)/2)*(1 + erf(oo - d/a**(S(1)/2))) \
           *exp(d**2/a)/(2*a**(S(1)/2))

def test_laplace_transform():
    LT = laplace_transform
    a, b, c, = symbols('a b c', positive=True)
    t = symbols('t')
    w = Symbol("w")
    f = Function("f")

    # Test unevaluated form
    assert laplace_transform(f(t), t, w) == LaplaceTransform(f(t), t, w)
    assert inverse_laplace_transform(f(w), w, t, plane=0) == InverseLaplaceTransform(f(w), w, t, 0)

    # test a bug
    spos = symbols('s', positive=True)
    assert LT(exp(t), t, spos)[:2] == (1/(spos - 1), True)

    # basic tests from wikipedia

    assert LT((t-a)**b*exp(-c*(t-a))*Heaviside(t-a), t, s) \
           == ((s + c)**(-b - 1)*exp(-a*s)*gamma(b + 1), -c, True)
    assert LT(t**a, t, s) == (s**(-a - 1)*gamma(a + 1), 0, True)
    assert LT(Heaviside(t), t, s) == (1/s, 0, True)
    assert LT(Heaviside(t - a), t, s) == (exp(-a*s)/s, 0, True)
    assert LT(1 - exp(-a*t), t, s) == (a/(s*(a + s)), 0, True)

    assert LT((exp(2*t)-1)*exp(-b - t)*Heaviside(t)/2, t, s, noconds=True) \
           == exp(-b)/(s**2 - 1)

    assert LT(exp(t), t, s)[:2] == (1/(s-1), 1)
    assert LT(exp(2*t), t, s)[:2] == (1/(s-2), 2)
    assert LT(exp(a*t), t, s)[:2] == (1/(s-a), a)

    assert LT(log(t/a), t, s) == ((log(a) + log(s) + EulerGamma)/(-s), 0, True)

    assert LT(erf(t), t, s) == ((-erf(s/2) + 1)*exp(s**2/4)/s, 0, True)

    assert LT(sin(a*t), t, s) == (a/(a**2 + s**2), 0, True)
    assert LT(cos(a*t), t, s) == (s/(a**2 + s**2), 0, True)
    # TODO would be nice to have these come out better
    assert LT(exp(-a*t)*sin(b*t), t, s) == (1/b/(1 + (a + s)**2/b**2), -a, True)
    assert LT(exp(-a*t)*cos(b*t), t, s) == \
           (1/(s + a)/(1 + b**2/(a + s)**2), -a, True)
    # TODO sinh, cosh have delicate cancellation

    assert LT(besselj(0, t), t, s) == (1/sqrt(1 + s**2), 0, True)
    assert LT(besselj(1, t), t, s) == (1 - 1/sqrt(1 + 1/s**2), 0, True)
    # TODO general order works, but is a *mess*
    # TODO besseli also works, but is an even greater mess

    # test a bug in conditions processing
    # TODO the auxiliary condition should be recognised/simplified
    assert LT(exp(t)*cos(t), t, s)[:-1] in [
        ((s - 1)/(s**2 - 2*s + 2), -oo),
        ((s - 1)/((s - 1)**2 + 1), -oo),
        ]

def test_recursive():
    from sympy import symbols, exp_polar

    a, b, c = symbols("a b c", positive=True)
    assert simplify(integrate(exp(-(x - a) ** 2) * exp(-(x - b) ** 2), (x, 0, oo))) == sqrt(2 * pi) / 4 * (
        1 + erf(sqrt(2) / 2 * (a + b))
    ) * exp(-a ** 2 - b ** 2 + (a + b) ** 2 / 2)
    assert simplify(integrate(exp(-(x - a) ** 2) * exp(-(x - b) ** 2) * exp(c * x), (x, 0, oo))) == sqrt(2 * pi) / 4 * (
        1 + erf(sqrt(2) / 4 * (2 * a + 2 * b + c))
    ) * exp(-a ** 2 - b ** 2 + (2 * a + 2 * b + c) ** 2 / 8)
    assert simplify(integrate(exp(-(x - a - b - c) ** 2), (x, 0, oo))) == sqrt(pi) / 2 * (1 + erf(a + b + c))
    assert simplify(integrate(exp(-(x + a + b + c) ** 2), (x, 0, oo))) == sqrt(pi) / 2 * (1 - erf(a + b + c))