Python sympy.series函数代码示例

def test_sec():
    x = symbols('x', real=True)
    z = symbols('z')

    assert sec.nargs == FiniteSet(1)

    assert sec(0) == 1
    assert sec(pi) == -1
    assert sec(pi/2) == zoo
    assert sec(-pi/2) == zoo
    assert sec(pi/6) == 2*sqrt(3)/3
    assert sec(pi/3) == 2
    assert sec(5*pi/2) == zoo
    assert sec(9*pi/7) == -sec(2*pi/7)
    assert sec(I) == 1/cosh(1)
    assert sec(x*I) == 1/cosh(x)
    assert sec(-x) == sec(x)

    assert sec(x).rewrite(exp) == 1/(exp(I*x)/2 + exp(-I*x)/2)
    assert sec(x).rewrite(sin) == sec(x)
    assert sec(x).rewrite(cos) == 1/cos(x)
    assert sec(x).rewrite(tan) == (tan(x/2)**2 + 1)/(-tan(x/2)**2 + 1)
    assert sec(x).rewrite(pow) == sec(x)
    assert sec(x).rewrite(sqrt) == sec(x)

    assert sec(z).conjugate() == sec(conjugate(z))

    assert (sec(z).as_real_imag() ==
    (cos(re(z))*cosh(im(z))/(sin(re(z))**2*sinh(im(z))**2 +
                             cos(re(z))**2*cosh(im(z))**2),
     sin(re(z))*sinh(im(z))/(sin(re(z))**2*sinh(im(z))**2 +
                             cos(re(z))**2*cosh(im(z))**2)))

    assert sec(x).expand(trig=True) == 1/cos(x)
    assert sec(2*x).expand(trig=True) == 1/(2*cos(x)**2 - 1)

    assert sec(x).is_real == True
    assert sec(z).is_real == None

    assert sec(x).as_leading_term() == sec(x)

    assert sec(0).is_bounded == True
    assert sec(x).is_bounded == None
    assert sec(pi/2).is_bounded == False

    assert series(sec(x), x, x0=0, n=6) == 1 + x**2/2 + 5*x**4/24 + O(x**6)

    # https://github.com/sympy/sympy/issues/7166
    assert series(sqrt(sec(x))) == 1 + x**2/4 + 7*x**4/96 + O(x**6)

    assert sec(x).diff(x) == tan(x)*sec(x)

    # Taylor Term checks
    assert sec(z).taylor_term(4, z) == 5*z**4/24
    assert sec(z).taylor_term(6, z) == 61*z**6/720
    assert sec(z).taylor_term(5, z) == 0

def test_csc():
    x = symbols('x', real=True)
    z = symbols('z')

    # https://github.com/sympy/sympy/issues/6707
    cosecant = csc('x')
    alternate = 1/sin('x')
    assert cosecant.equals(alternate) == True
    assert alternate.equals(cosecant) == True

    assert csc.nargs == FiniteSet(1)

    assert csc(0) == zoo
    assert csc(pi) == zoo

    assert csc(pi/2) == 1
    assert csc(-pi/2) == -1
    assert csc(pi/6) == 2
    assert csc(pi/3) == 2*sqrt(3)/3
    assert csc(5*pi/2) == 1
    assert csc(9*pi/7) == -csc(2*pi/7)
    assert csc(I) == -I/sinh(1)
    assert csc(x*I) == -I/sinh(x)
    assert csc(-x) == -csc(x)

    assert csc(x).rewrite(exp) == 2*I/(exp(I*x) - exp(-I*x))
    assert csc(x).rewrite(sin) == 1/sin(x)
    assert csc(x).rewrite(cos) == csc(x)
    assert csc(x).rewrite(tan) == (tan(x/2)**2 + 1)/(2*tan(x/2))

    assert csc(z).conjugate() == csc(conjugate(z))

    assert (csc(z).as_real_imag() ==
            (sin(re(z))*cosh(im(z))/(sin(re(z))**2*cosh(im(z))**2 +
                                     cos(re(z))**2*sinh(im(z))**2),
             -cos(re(z))*sinh(im(z))/(sin(re(z))**2*cosh(im(z))**2 +
                          cos(re(z))**2*sinh(im(z))**2)))

    assert csc(x).expand(trig=True) == 1/sin(x)
    assert csc(2*x).expand(trig=True) == 1/(2*sin(x)*cos(x))

    assert csc(x).is_real == True
    assert csc(z).is_real == None

    assert csc(x).as_leading_term() == csc(x)

    assert csc(0).is_bounded == False
    assert csc(x).is_bounded == None
    assert csc(pi/2).is_bounded == True

    assert series(csc(x), x, x0=pi/2, n=6) == \
        1 + (x - pi/2)**2/2 + 5*(x - pi/2)**4/24 + O((x - pi/2)**6, (x, pi/2))
    assert series(csc(x), x, x0=0, n=6) == \
            1/x + x/6 + 7*x**3/360 + 31*x**5/15120 + O(x**6)

    assert csc(x).diff(x) == -cot(x)*csc(x)

def exercise_6_43():

    z = sp.Symbol('z')

    b_of_z = 1/(1 - z)

    sp.pprint('\n' + 'b_of_z:')
    sp.pprint(b_of_z)

    e_b_prime_of_z = z**2/(1 - z)**2

    sp.pprint('\n' +'e_b_prime_of_z:')
    sp.pprint(e_b_prime_of_z)

    e_b_prime_of_n = sp.series(e_b_prime_of_z, z, 0, 10)

    sp.pprint('\n' +'e_b_prime_of_n:')
    sp.pprint(e_b_prime_of_n)

    e_b_of_z = sp.integrate(e_b_prime_of_z, z)
    constant = e_b_of_z.subs(z, 0)
    e_b_of_z = e_b_of_z - constant

    sp.pprint('\n' +'e_b_of_z:')
    sp.pprint(e_b_of_z)

    e_b_of_n = sp.series(e_b_of_z, z, 0, 10)

    sp.pprint('\n' +'e_b_of_n:')
    sp.pprint(e_b_of_n)

    e_b_prime_of_z = sp.diff(e_b_of_z, z)

    sp.pprint('\n' +'e_b_prime_of_z:')
    sp.pprint(e_b_prime_of_z)

    f_of_z = sp.together((1 - z)**2*e_b_prime_of_z)

    sp.pprint('\n' +'f_of_z:')
    sp.pprint(f_of_z)

    F_of_z = sp.integrate(f_of_z, z)

    sp.pprint('\n' +'F_of_z:')
    sp.pprint(F_of_z)

    c_b_of_z = (e_b_of_z.subs(z, 0) + F_of_z)/(1 - z)**2

    sp.pprint('\n' +'c_b_of_z:')
    sp.pprint(c_b_of_z)

    c_b_of_n = sp.series(c_b_of_z, z, 0, 10)

    sp.pprint('\n' +'c_b_of_n:')
    sp.pprint(c_b_of_n)

def test_sec():
    x = symbols('x', real=True)
    z = symbols('z')

    assert sec.nargs == 1

    assert sec(0) == 1
    assert sec(pi) == -1
    assert sec(pi/2) == oo
    assert sec(-pi/2) == oo
    assert sec(pi/6) == 2*sqrt(3)/3
    assert sec(pi/3) == 2
    assert sec(5*pi/2) == oo
    assert sec(9*pi/7) == -sec(2*pi/7)
    assert sec(I) == 1/cosh(1)
    assert sec(x*I) == 1/cosh(x)
    assert sec(-x) == sec(x)

    assert sec(x).rewrite(exp) == 1/(exp(I*x)/2 + exp(-I*x)/2)
    assert sec(x).rewrite(sin) == sec(x)
    assert sec(x).rewrite(cos) == 1/cos(x)
    assert sec(x).rewrite(tan) == (tan(x/2)**2 + 1)/(-tan(x/2)**2 + 1)
    assert sec(x).rewrite(pow) == sec(x)
    assert sec(x).rewrite(sqrt) == sec(x)

    assert sec(z).conjugate() == sec(conjugate(z))

    assert (sec(z).as_real_imag() ==
    (cos(re(z))*cosh(im(z))/(sin(re(z))**2*sinh(im(z))**2 +
                             cos(re(z))**2*cosh(im(z))**2),
     sin(re(z))*sinh(im(z))/(sin(re(z))**2*sinh(im(z))**2 +
                             cos(re(z))**2*cosh(im(z))**2)))

    assert sec(x).expand(trig=True) == 1/cos(x)
    assert sec(2*x).expand(trig=True) == 1/(2*cos(x)**2 - 1)

    assert sec(x).is_real == True
    assert sec(z).is_real == None

    assert sec(x).as_leading_term() == sec(x)

    assert sec(0).is_bounded == True
    assert sec(x).is_bounded == None
    assert sec(pi/2).is_bounded == False

    assert series(sec(x), x, x0=0, n=6) == 1 + x**2/2 + 5*x**4/24 + O(x**6)

    # https://code.google.com/p/sympy/issues/detail?id=4067
    assert series(sqrt(sec(x))) == 1 + x**2/4 + 7*x**4/96 + O(x**6)

    # https://code.google.com/p/sympy/issues/detail?id=4068
    assert (series(sqrt(sec(x)), x, x0=pi*3/2, n=4) ==
            1/sqrt(x) +x**(S(3)/2)/12 + x**(S(7)/2)/160 + O(x**4))

    assert sec(x).diff(x) == tan(x)*sec(x)

def test_2124():
    assert series(1, x) == 1
    assert S(0).lseries(x).next() == 0
    assert cos(x).series() == cos(x).series(x)
    raises(ValueError, lambda: cos(x + y).series())
    raises(ValueError, lambda: x.series(dir=""))

    assert (cos(x).series(x, 1).removeO().subs(x, x - 1) -
            cos(x + 1).series(x).removeO().subs(x, x - 1)).expand() == 0
    e = cos(x).series(x, 1, n=None)
    assert [e.next() for i in range(2)] == [cos(1), -((x - 1)*sin(1))]
    e = cos(x).series(x, 1, n=None, dir='-')
    assert [e.next() for i in range(2)] == [cos(1), (1 - x)*sin(1)]
    # the following test is exact so no need for x -> x - 1 replacement
    assert abs(x).series(x, 1, dir='-') == x
    assert exp(x).series(x, 1, dir='-', n=3).removeO().subs(x, x - 1) == \
           E + E*(x - 1) + E*(x - 1)**2/2

    D = Derivative
    assert D(x**2 + x**3*y**2, x, 2, y, 1).series(x).doit() == 12*x*y
    assert D(cos(x), x).lseries().next() == D(1, x)
    assert D(
        exp(x), x).series(n=3) == D(1, x) + D(x, x) + D(x**2/2, x) + O(x**3)

    assert Integral(x, (x, 1, 3), (y, 1, x)).series(x) == -4 + 4*x

    assert (1 + x + O(x**2)).getn() == 2
    assert (1 + x).getn() is None

    assert ((1/sin(x))**oo).series() == oo
    logx = Symbol('logx')
    assert ((sin(x))**y).nseries(x, n=1, logx=logx) \
           == exp(y*logx) + O(x*exp(y*logx), x)

    raises(NotImplementedError, lambda: series(Function("f")(x)))

    assert sin(1/x).series(x, oo, n=5) == 1/x - 1/(6*x**3)
    assert abs(x).series(x, oo, n=5, dir='+') == x
    assert abs(x).series(x, -oo, n=5, dir='-') == -x
    assert abs(-x).series(x, oo, n=5, dir='+') == x
    assert abs(-x).series(x, -oo, n=5, dir='-') == -x

    assert exp(x*log(x)).series(n=3) == \
           1 + x*log(x) + x**2*log(x)**2/2 + O(x**3*log(x)**3)
    # XXX is this right? If not, fix "ngot > n" handling in expr.
    p = Symbol('p', positive=True)
    assert exp(sqrt(p)**3*log(p)).series(n=3) == \
        1 + p**S('3/2')*log(p) + O(p**3*log(p)**3)

    assert exp(sin(x)*log(x)).series(n=2) == 1 + x*log(x) + O(x**2*log(x)**2)

def test_airybiprime():
    z = Symbol('z', real=False)
    t = Symbol('t', negative=True)
    p = Symbol('p', positive=True)

    assert isinstance(airybiprime(z), airybiprime)

    assert airybiprime(0) == 3**(S(1)/6)/gamma(S(1)/3)
    assert airybiprime(oo) == oo
    assert airybiprime(-oo) == 0

    assert diff(airybiprime(z), z) == z*airybi(z)

    assert series(airybiprime(z), z, 0, 3) == (
        3**(S(1)/6)/gamma(S(1)/3) + 3**(S(5)/6)*z**2/(6*gamma(S(2)/3)) + O(z**3))

    assert airybiprime(z).rewrite(hyper) == (
        3**(S(5)/6)*z**2*hyper((), (S(5)/3,), z**S(3)/9)/(6*gamma(S(2)/3)) +
        3**(S(1)/6)*hyper((), (S(1)/3,), z**S(3)/9)/gamma(S(1)/3))

    assert isinstance(airybiprime(z).rewrite(besselj), airybiprime)
    assert airyai(t).rewrite(besselj) == (
        sqrt(-t)*(besselj(-S(1)/3, 2*(-t)**(S(3)/2)/3) +
                  besselj(S(1)/3, 2*(-t)**(S(3)/2)/3))/3)
    assert airybiprime(z).rewrite(besseli) == (
        sqrt(3)*(z**2*besseli(S(2)/3, 2*z**(S(3)/2)/3)/(z**(S(3)/2))**(S(2)/3) +
                 (z**(S(3)/2))**(S(2)/3)*besseli(-S(2)/3, 2*z**(S(3)/2)/3))/3)
    assert airybiprime(p).rewrite(besseli) == (
        sqrt(3)*p*(besseli(-S(2)/3, 2*p**(S(3)/2)/3) + besseli(S(2)/3, 2*p**(S(3)/2)/3))/3)

    assert expand_func(airybiprime(2*(3*z**5)**(S(1)/3))) == (
        sqrt(3)*(z**(S(5)/3)/(z**5)**(S(1)/3) - 1)*airyaiprime(2*3**(S(1)/3)*z**(S(5)/3))/2 +
        (z**(S(5)/3)/(z**5)**(S(1)/3) + 1)*airybiprime(2*3**(S(1)/3)*z**(S(5)/3))/2)

def test_airyai():
    z = Symbol('z', real=False)
    t = Symbol('t', negative=True)
    p = Symbol('p', positive=True)

    assert isinstance(airyai(z), airyai)

    assert airyai(0) == 3**(S(1)/3)/(3*gamma(S(2)/3))
    assert airyai(oo) == 0
    assert airyai(-oo) == 0

    assert diff(airyai(z), z) == airyaiprime(z)

    assert series(airyai(z), z, 0, 3) == (
        3**(S(5)/6)*gamma(S(1)/3)/(6*pi) - 3**(S(1)/6)*z*gamma(S(2)/3)/(2*pi) + O(z**3))

    assert airyai(z).rewrite(hyper) == (
        -3**(S(2)/3)*z*hyper((), (S(4)/3,), z**S(3)/9)/(3*gamma(S(1)/3)) +
         3**(S(1)/3)*hyper((), (S(2)/3,), z**S(3)/9)/(3*gamma(S(2)/3)))

    assert isinstance(airyai(z).rewrite(besselj), airyai)
    assert airyai(t).rewrite(besselj) == (
        sqrt(-t)*(besselj(-S(1)/3, 2*(-t)**(S(3)/2)/3) +
                  besselj(S(1)/3, 2*(-t)**(S(3)/2)/3))/3)
    assert airyai(z).rewrite(besseli) == (
        -z*besseli(S(1)/3, 2*z**(S(3)/2)/3)/(3*(z**(S(3)/2))**(S(1)/3)) +
         (z**(S(3)/2))**(S(1)/3)*besseli(-S(1)/3, 2*z**(S(3)/2)/3)/3)
    assert airyai(p).rewrite(besseli) == (
        sqrt(p)*(besseli(-S(1)/3, 2*p**(S(3)/2)/3) -
                 besseli(S(1)/3, 2*p**(S(3)/2)/3))/3)

    assert expand_func(airyai(2*(3*z**5)**(S(1)/3))) == (
        -sqrt(3)*(-1 + (z**5)**(S(1)/3)/z**(S(5)/3))*airybi(2*3**(S(1)/3)*z**(S(5)/3))/6 +
         (1 + (z**5)**(S(1)/3)/z**(S(5)/3))*airyai(2*3**(S(1)/3)*z**(S(5)/3))/2)

    def expectation(self, expr, var, evaluate=True, **kwargs):
        """ Expectation of expression over distribution """
        # TODO: support discrete sets with non integer stepsizes

        if evaluate:
            try:
                p = poly(expr, var)

                t = Dummy('t', real=True)

                mgf = self.moment_generating_function(t)
                deg = p.degree()
                taylor = poly(series(mgf, t, 0, deg + 1).removeO(), t)
                result = 0
                for k in range(deg+1):
                    result += p.coeff_monomial(var ** k) * taylor.coeff_monomial(t ** k) * factorial(k)

                return result

            except PolynomialError:
                return summation(expr * self.pdf(var),
                                 (var, self.set.inf, self.set.sup), **kwargs)

        else:
            return Sum(expr * self.pdf(var),
                         (var, self.set.inf, self.set.sup), **kwargs)

def t_expand(expr, max_deg):
    """
    Takes an expression f(x,y) and computes the Taylor expansion 
    in x and y up to bi-degree fixed by max_deg.
    For example: 
    x / (1 - y) 
    with value of max_deg = 3 will give 
    x + x * y
    up to bi-degree 2.
    """
    f = expr.subs([(x, t*x), (y, t*y)])
    return series(f, t, 0, max_deg).removeO().subs(t, 1)

def t_degree(expr, max_deg):
    """
    Takes an expression f(x,y) and computes the Taylor expansion 
    in x and y up to bi-degree fixed by max_deg. Then, picks the term 
    of highest bi-degree, and returns the value of the degree.
    For example: 
    x / (1 - y) 
    with value of max_deg = 3 will give 
    t * x + t**2 * x * y
    up to bi-degree 2. Will return 2.
    """
    f = expr.subs([(x, t*x), (y, t*y)])
    return degree(series(f, t, 0, max_deg).removeO(), t)

def t_min_degree(expr, max_deg, large_number=1000):
    """
    Takes an expression f(x,y) and computes the Taylor expansion 
    in x and y up to bi-degree fixed by max_deg. Then, picks the term 
    of lowest bi-degree.
    For example: 
    x / (1 - y) 
    with value of max_deg = 3 will give 
    t * x + t**2 * x * y
    up to bi-degree 2. Will return x.
    """
    f = expr.subs([(x, x / t), (y, y / t)])
    f_t_series = t**large_number * series(f, t, 0, max_deg).removeO() 
    leading_term = LT(f_t_series.expand(), t).subs(t, 1)
    return leading_term 

 def expectation(self, expr, var, evaluate=True, **kwargs):
     """ Expectation of expression over distribution """
     if evaluate:
         try:
             p = poly(expr, var)
             t = Dummy('t', real=True)
             mgf = self._moment_generating_function(t)
             if mgf is None:
                 return integrate(expr * self.pdf(var), (var, self.set), **kwargs)
             deg = p.degree()
             taylor = poly(series(mgf, t, 0, deg + 1).removeO(), t)
             result = 0
             for k in range(deg+1):
                 result += p.coeff_monomial(var ** k) * taylor.coeff_monomial(t ** k) * factorial(k)
             return result
         except PolynomialError:
             return integrate(expr * self.pdf(var), (var, self.set), **kwargs)
     else:
         return Integral(expr * self.pdf(var), (var, self.set), **kwargs)

    def _get_dict(self, other):
        """
        Returns a dictionary representing the (exponent, coefficient) pairs
        of the Puiseux series expansion of `other`.

        If `other` is a `PuiseuxSeries` then we take the minimum order
        of `other` and self. If `other is a SymPy expression then we
        compute the Taylor series in (x-alpha) of other up to `self.order`.
        """
        # if other is a Puiseux series (also centered at self.alpha)
        # then the order of the new puiseux series is the smaller of
        # the two
        if isinstance(other, PuiseuxSeries):
            if self.alpha != other.alpha:
                raise ValueError("Puiseux series must be centered at " + \
                                 "the same point.")         
            order = min(self.order, other.order)
            d = other.d

        # if other is a SymPy Expr then compute the series expansion
        # of other up to self's order and store the exponent, coeff
        # pairs in d
        elif isinstance(other, sympy.Expr):
            d = {}
            order = self.order
         
            s = sympy.series(other, self.var, x0=self.alpha, n=None)
            for term in s:
                # we have to shift because SymPy has trouble
                # collecting expressions such as c*(x-a) in terms of
                # (x-a). Example: 1-x should be -1*(x-1) but SymPy
                # cannot detect this.
                term = term.subs(x,x+self.alpha)
                coeff, exp = term.as_coeff_exponent(self.var)
                if exp < order:
                    d[exp] = coeff
                else:
                    break

        return order,d

def get_coefficients(function, p, x0, n):

	taylor_series = sp.series(expr=function, x=p, x0=x0, n=n).removeO()
	
	coefficients_dict_sympy_pow_keys = taylor_series.as_coefficients_dict()

	coefficients_dict_string_keys = {}
	for key in coefficients_dict_sympy_pow_keys.keys():
		coefficients_dict_string_keys[str(key)] = coefficients_dict_sympy_pow_keys[key]

	taylor_coefficients = []
	for i in range(n):
		if i == 0:
			key = '1'
		elif i == 1:
			key = 'p'
		else:
			key = "(p - 0.5)**{}".format(i)
		taylor_coefficients.append( coefficients_dict_string_keys[key] / (2**i) )


	coeffs = [ sp.N(c, precision) for c in taylor_coefficients ]

	return coeffs

def test_exp2():
    e1 = exp(cos(x)).series(x, 0)
    e2 = series(exp(cos(x)), x, 0)
    assert e1 == e2