Python sympy.sqrt函数代码示例

def root_mul_rule(integral):
    integrand, symbol = integral
    a = sympy.Wild("a", exclude=[symbol])
    b = sympy.Wild("b", exclude=[symbol])
    c = sympy.Wild("c")
    match = integrand.match(sympy.sqrt(a * symbol + b) * c)

    if not match:
        return

    a, b, c = match[a], match[b], match[c]
    d = sympy.Wild("d", exclude=[symbol])
    e = sympy.Wild("e", exclude=[symbol])
    f = sympy.Wild("f")
    recursion_test = c.match(sympy.sqrt(d * symbol + e) * f)
    if recursion_test:
        return

    u = sympy.Dummy("u")
    u_func = sympy.sqrt(a * symbol + b)
    integrand = integrand.subs(u_func, u)
    integrand = integrand.subs(symbol, (u ** 2 - b) / a)
    integrand = integrand * 2 * u / a
    next_step = integral_steps(integrand, u)
    if next_step:
        return URule(u, u_func, None, next_step, integrand, symbol)

def test_pretty_functions():
    f = Function('f')

    # Simple
    assert pretty( (2*x + exp(x)) ) in [' x      \ne  + 2*x', '       x\n2*x + e ']
    assert pretty(abs(x)) == '|x|'
    assert pretty(abs(x/(x**2+1))) in [
            '|  x   |\n|------|\n|     2|\n|1 + x |',
            '|  x   |\n|------|\n| 2    |\n|x  + 1|']
    assert pretty(conjugate(x)) == '_\nx'
    assert pretty(conjugate(f(x+1))) in [
            '________\nf(1 + x)',
            '________\nf(x + 1)']

    # Univariate/Multivariate functions
    assert pretty(f(x)) == 'f(x)'
    assert pretty(f(x, y)) == 'f(x, y)'
    assert pretty(f(x/(y+1), y)) in [
            ' /  x     \\\nf|-----, y|\n \\1 + y   /',
            ' /  x     \\\nf|-----, y|\n \\y + 1   /',
            ]

    # Nesting of square roots
    assert pretty( sqrt((sqrt(x+1))+1) ) in [
            '   _______________\n  /       _______ \n\\/  1 + \\/ 1 + x  ',
            '   _______________\n  /   _______     \n\\/  \\/ x + 1  + 1 ']
    # Function powers
    assert pretty( sin(x)**2 ) == '   2   \nsin (x)'

    # Conjugates
    a,b = map(Symbol, 'ab')

def test_director_circle():
    x, y, a, b = symbols('x y a b')
    e = Ellipse((x, y), a, b)
    # the general result
    assert e.director_circle() == Circle((x, y), sqrt(a**2 + b**2))
    # a special case where Ellipse is a Circle
    assert Circle((3, 4), 8).director_circle() == Circle((3, 4), 8*sqrt(2))

def test_as_ordered_terms():
    f, g = symbols("f,g", cls=Function)

    assert x.as_ordered_terms() == [x]
    assert (sin(x) ** 2 * cos(x) + sin(x) * cos(x) ** 2 + 1).as_ordered_terms() == [
        sin(x) ** 2 * cos(x),
        sin(x) * cos(x) ** 2,
        1,
    ]

    args = [f(1), f(2), f(3), f(1, 2, 3), g(1), g(2), g(3), g(1, 2, 3)]
    expr = Add(*args)

    assert expr.as_ordered_terms() == args

    assert (1 + 4 * sqrt(3) * pi * x).as_ordered_terms() == [4 * pi * x * sqrt(3), 1]

    assert (2 + 3 * I).as_ordered_terms() == [2, 3 * I]
    assert (-2 + 3 * I).as_ordered_terms() == [-2, 3 * I]
    assert (2 - 3 * I).as_ordered_terms() == [2, -3 * I]
    assert (-2 - 3 * I).as_ordered_terms() == [-2, -3 * I]

    assert (4 + 3 * I).as_ordered_terms() == [4, 3 * I]
    assert (-4 + 3 * I).as_ordered_terms() == [-4, 3 * I]
    assert (4 - 3 * I).as_ordered_terms() == [4, -3 * I]
    assert (-4 - 3 * I).as_ordered_terms() == [-4, -3 * I]

    f = x ** 2 * y ** 2 + x * y ** 4 + y + 2

    assert f.as_ordered_terms(order="lex") == [x ** 2 * y ** 2, x * y ** 4, y, 2]
    assert f.as_ordered_terms(order="grlex") == [x * y ** 4, x ** 2 * y ** 2, y, 2]
    assert f.as_ordered_terms(order="rev-lex") == [2, y, x * y ** 4, x ** 2 * y ** 2]
    assert f.as_ordered_terms(order="rev-grlex") == [2, y, x ** 2 * y ** 2, x * y ** 4]

def eval_trigsubstitution(theta, func, rewritten, substep, restriction, integrand, symbol):
    func = func.subs(sympy.sec(theta), 1/sympy.cos(theta))

    trig_function = list(func.find(TrigonometricFunction))
    assert len(trig_function) == 1
    trig_function = trig_function[0]
    relation = sympy.solve(symbol - func, trig_function)
    assert len(relation) == 1
    numer, denom = sympy.fraction(relation[0])

    if isinstance(trig_function, sympy.sin):
        opposite = numer
        hypotenuse = denom
        adjacent = sympy.sqrt(denom**2 - numer**2)
        inverse = sympy.asin(relation[0])
    elif isinstance(trig_function, sympy.cos):
        adjacent = numer
        hypotenuse = denom
        opposite = sympy.sqrt(denom**2 - numer**2)
        inverse = sympy.acos(relation[0])
    elif isinstance(trig_function, sympy.tan):
        opposite = numer
        adjacent = denom
        hypotenuse = sympy.sqrt(denom**2 + numer**2)
        inverse = sympy.atan(relation[0])

    substitution = [
        (sympy.sin(theta), opposite/hypotenuse),
        (sympy.cos(theta), adjacent/hypotenuse),
        (sympy.tan(theta), opposite/adjacent),
        (theta, inverse)
    ]
    return sympy.Piecewise(
        (_manualintegrate(substep).subs(substitution).trigsimp(), restriction)
    )

def test_issue_1364():
    assert solve(-a*x + 2*x*log(x), x) == [exp(a/2)]
    assert solve(a/x + exp(x/2), x) == [2*LambertW(-a/2)]
    assert solve(x**x) == []
    assert solve(x**x - 2) == [exp(LambertW(log(2)))]
    assert solve(((x - 3)*(x - 2))**((x - 3)*(x - 4))) == [2]
    assert solve((a/x + exp(x/2)).diff(x), x) == [4*LambertW(sqrt(2)*sqrt(a)/4)]

def test_uselogcombine_1():
    assert solveset_real(log(x - 3) + log(x + 3), x) == \
        FiniteSet(sqrt(10))
    assert solveset_real(log(x + 1) - log(2*x - 1), x) == FiniteSet(2)
    assert solveset_real(log(x + 3) + log(1 + 3/x) - 3) == FiniteSet(
        -3 + sqrt(-12 + exp(3))*exp(S(3)/2)/2 + exp(3)/2,
        -sqrt(-12 + exp(3))*exp(S(3)/2)/2 - 3 + exp(3)/2)

    def eta_fil(self, x, V_app, apprx=(0, 0, 0, 0)):
        m_eff = self.m_r * const.electron_mass

        mpmath.mp.dps = 20
        x0 = Symbol('x0')  # eta_fil
        x1 = Symbol('x1')  # eta_ac
        x2 = Symbol('x2')  # eta_hop
        x3 = Symbol('x3')  # V_tunnel

        f0 = const.Boltzmann * self.T / (1 - self.alpha) / const.elementary_charge / self.z * \
             ln(self.A_fil/self.A_ac*(exp(- self.alpha * const.elementary_charge * self.z / const.Boltzmann / self.T * x0) - 1) + 1) - x1# eta_ac = f(eta_fil) x1 = f(x0)
        f1 = x*2*const.Boltzmann*self.T/self.a/self.z/const.elementary_charge*\
             asinh(self.j_0et/self.j_0hop*(exp(- self.alpha * const.elementary_charge * self.z / const.Boltzmann / self.T * x0) - 1)) - x2# eta_hop = f(eta_fil)
        f2 = x1 - x0 + x2 - x3

        f3 = -V_app + ((self.C * 3 * sqrt(2 * m_eff * ((4+x3/2)*const.elementary_charge)) / 2 / x * (const.elementary_charge / const.Planck)**2 * \
             exp(- 4 * const.pi * x / const.Planck * sqrt(2 * m_eff * ((4+x3/2)*const.elementary_charge))) * self.A_fil*x3)
                       + (self.j_0et*self.A_fil*(exp(-self.alpha*const.elementary_charge*self.z*x0/const.Boltzmann/self.T) - 1))) * (self.R_el + self.R_S + self.rho_fil*(self.L - x) / self.A_fil) \
             + x3

        eta_fil, eta_ac, eta_hop, V_tunnel = nsolve((f0, f1, f2, f3), [x0, x1, x2, x3], apprx)
        eta_fil = np.real(np.complex128(eta_fil))
        eta_ac = np.real(np.complex128(eta_ac))
        eta_hop = np.real(np.complex128(eta_hop))
        V_tunnel = np.real(np.complex128(V_tunnel))
        current = ((self.C * 3 * sqrt(2 * m_eff * ((4+V_tunnel)*const.elementary_charge)) / 2 / x * (const.elementary_charge / const.Planck)**2 * \
            exp(- 4 * const.pi * x / const.Planck * sqrt(2 * m_eff * ((4+V_tunnel)*const.elementary_charge))) * self.A_fil*V_tunnel)
                       + (self.j_0et*self.A_fil*(exp(-self.alpha*const.elementary_charge*self.z*eta_fil/const.Boltzmann/self.T) - 1)))
        print(eta_fil, eta_ac, eta_hop, V_tunnel)
        # print(eta_ac - eta_fil + eta_hop - V_tunnel)
        return eta_fil, eta_ac, eta_hop, V_tunnel, current

def test_cse_single2():
    # Simple substitution, test for being able to pass the expression directly
    e = Add(Pow(x+y,2), sqrt(x+y))
    substs, reduced = cse(e, optimizations=[])
    assert substs == [(x0, x+y)]
    assert reduced == [sqrt(x0) + x0**2]
    assert isinstance(cse(Matrix([[1]]))[1][0], Matrix)

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_roots_quartic():
    assert roots_quartic(Poly(x**4, x)) == [0, 0, 0, 0]
    assert roots_quartic(Poly(x**4 + x**3, x)) in [
        [-1,0,0,0],
        [0,-1,0,0],
        [0,0,-1,0],
        [0,0,0,-1]
    ]
    assert roots_quartic(Poly(x**4 - x**3, x)) in [
        [1,0,0,0],
        [0,1,0,0],
        [0,0,1,0],
        [0,0,0,1]
    ]

    lhs = roots_quartic(Poly(x**4 + x, x))
    rhs = [S.Half + I*sqrt(3)/2, S.Half - I*sqrt(3)/2, S.Zero, -S.One]

    assert sorted(lhs, key=hash) == sorted(rhs, key=hash)

    # test of all branches of roots quartic
    for i, (a, b, c, d) in enumerate([(1, 2, 3, 0),
                                      (3, -7, -9, 9),
                                      (1, 2, 3, 4),
                                      (1, 2, 3, 4),
                                      (-7, -3, 3, -6),
                                      (-3, 5, -6, -4)]):
        if i == 2:
            c = -a*(a**2/S(8) - b/S(2))
        elif i == 3:
            d = a*(a*(3*a**2/S(256) - b/S(16)) + c/S(4))
        eq = x**4 + a*x**3 + b*x**2 + c*x + d
        ans = roots_quartic(Poly(eq, x))
        assert all([eq.subs(x, ai).n(chop=True) == 0 for ai in ans])

def test_expr_sorting():
    f, g = symbols('f,g', cls=Function)

    exprs = [1/x**2, 1/x, sqrt(sqrt(x)), sqrt(x), x, sqrt(x)**3, x**2]
    assert sorted(exprs, key=default_sort_key) == exprs

    exprs = [x, 2*x, 2*x**2, 2*x**3, x**n, 2*x**n, sin(x), sin(x)**n, sin(x**2), cos(x), cos(x**2), tan(x)]
    assert sorted(exprs, key=default_sort_key) == exprs

    exprs = [x + 1, x**2 + x + 1, x**3 + x**2 + x + 1]
    assert sorted(exprs, key=default_sort_key) == exprs

    exprs = [S(4), x - 3*I/2, x + 3*I/2, x - 4*I + 1, x + 4*I + 1]
    assert sorted(exprs, key=default_sort_key) == exprs

    exprs = [f(1), f(2), f(3), f(1, 2, 3), g(1), g(2), g(3), g(1, 2, 3)]
    assert sorted(exprs, key=default_sort_key) == exprs

    exprs = [f(x), g(x), exp(x), sin(x), cos(x), factorial(x)]
    assert sorted(exprs, key=default_sort_key) == exprs

    exprs = [Tuple(x, y), Tuple(x, z), Tuple(x, y, z)]
    assert sorted(exprs, key=default_sort_key) == exprs

    exprs = [[3], [1, 2]]
    assert sorted(exprs, key=default_sort_key) == exprs

    exprs = [[1, 2], [2, 3]]
    assert sorted(exprs, key=default_sort_key) == exprs

    exprs = [[1, 2], [1, 2, 3]]
    assert sorted(exprs, key=default_sort_key) == exprs

    exprs = [{x: -y}, {x: y}]
    assert sorted(exprs, key=default_sort_key) == exprs

def test_rayleigh():
    sigma = Symbol("sigma", positive=True)

    X = Rayleigh('x', sigma)
    assert density(X)(x) ==  x*exp(-x**2/(2*sigma**2))/sigma**2
    assert E(X) == sqrt(2)*sqrt(pi)*sigma/2
    assert variance(X) == -pi*sigma**2/2 + 2*sigma**2

def test_lognormal():
    mean = Symbol('mu', real=True, finite=True)
    std = Symbol('sigma', positive=True, real=True, finite=True)
    X = LogNormal('x', mean, std)
    # The sympy integrator can't do this too well
    #assert E(X) == exp(mean+std**2/2)
    #assert variance(X) == (exp(std**2)-1) * exp(2*mean + std**2)

    # Right now, only density function and sampling works
    # Test sampling: Only e^mean in sample std of 0
    for i in range(3):
        X = LogNormal('x', i, 0)
        assert S(sample(X)) == N(exp(i))
    # The sympy integrator can't do this too well
    #assert E(X) ==

    mu = Symbol("mu", real=True)
    sigma = Symbol("sigma", positive=True)

    X = LogNormal('x', mu, sigma)
    assert density(X)(x) == (sqrt(2)*exp(-(-mu + log(x))**2
                                    /(2*sigma**2))/(2*x*sqrt(pi)*sigma))

    X = LogNormal('x', 0, 1)  # Mean 0, standard deviation 1
    assert density(X)(x) == sqrt(2)*exp(-log(x)**2/2)/(2*x*sqrt(pi))

def test_checking():
    assert set(solve(x*(x - y/x),x, check=False)) == set([sqrt(y), S(0), -sqrt(y)])
    assert set(solve(x*(x - y/x),x, check=True)) == set([sqrt(y), -sqrt(y)])
    # {x: 0, y: 4} sets denominator to 0 in the following so system should return None
    assert solve((1/(1/x + 2), 1/(y - 3) - 1)) is None
    # 0 sets denominator of 1/x to zero so [] is returned
    assert solve(1/(1/x + 2)) == []