Python functions.binomial函数代码示例

def test_catalan():
    n = Symbol('n', integer=True)
    m = Symbol('n', integer=True, positive=True)

    catalans = [1, 1, 2, 5, 14, 42, 132, 429, 1430, 4862, 16796, 58786]
    for i, c in enumerate(catalans):
        assert catalan(i) == c
        assert catalan(n).rewrite(factorial).subs(n, i) == c
        assert catalan(n).rewrite(Product).subs(n, i).doit() == c

    assert catalan(x) == catalan(x)
    assert catalan(2*x).rewrite(binomial) == binomial(4*x, 2*x)/(2*x + 1)
    assert catalan(Rational(1, 2)).rewrite(gamma) == 8/(3*pi)
    assert catalan(Rational(1, 2)).rewrite(factorial).rewrite(gamma) ==\
        8 / (3 * pi)
    assert catalan(3*x).rewrite(gamma) == 4**(
        3*x)*gamma(3*x + Rational(1, 2))/(sqrt(pi)*gamma(3*x + 2))
    assert catalan(x).rewrite(hyper) == hyper((-x + 1, -x), (2,), 1)

    assert catalan(n).rewrite(factorial) == factorial(2*n) / (factorial(n + 1)
                                                              * factorial(n))
    assert isinstance(catalan(n).rewrite(Product), catalan)
    assert isinstance(catalan(m).rewrite(Product), Product)

    assert diff(catalan(x), x) == (polygamma(
        0, x + Rational(1, 2)) - polygamma(0, x + 2) + log(4))*catalan(x)

    assert catalan(x).evalf() == catalan(x)
    c = catalan(S.Half).evalf()
    assert str(c) == '0.848826363156775'
    c = catalan(I).evalf(3)
    assert str((re(c), im(c))) == '(0.398, -0.0209)'

def test_catalan():
    assert catalan(1) == 1
    assert catalan(2) == 2
    assert catalan(3) == 5
    assert catalan(4) == 14

    assert catalan(x) == catalan(x)
    assert catalan(2*x).rewrite(binomial) == binomial(4*x, 2*x)/(2*x + 1)
    assert catalan(Rational(1, 2)).rewrite(gamma) == 8/(3*pi)
    assert catalan(3*x).rewrite(gamma) == 4**(
        3*x)*gamma(3*x + Rational(1, 2))/(sqrt(pi)*gamma(3*x + 2))
    assert catalan(x).rewrite(hyper) == hyper((-x + 1, -x), (2,), 1)

    assert diff(catalan(x), x) == (polygamma(
        0, x + Rational(1, 2)) - polygamma(0, x + 2) + log(4))*catalan(x)

    c = catalan(0.5).evalf()
    assert str(c) == '0.848826363156775'

def rsolve_poly(coeffs, f, n, **hints):
    """
    Given linear recurrence operator `\operatorname{L}` of order
    `k` with polynomial coefficients and inhomogeneous equation
    `\operatorname{L} y = f`, where `f` is a polynomial, we seek for
    all polynomial solutions over field `K` of characteristic zero.

    The algorithm performs two basic steps:

        (1) Compute degree `N` of the general polynomial solution.
        (2) Find all polynomials of degree `N` or less
            of `\operatorname{L} y = f`.

    There are two methods for computing the polynomial solutions.
    If the degree bound is relatively small, i.e. it's smaller than
    or equal to the order of the recurrence, then naive method of
    undetermined coefficients is being used. This gives system
    of algebraic equations with `N+1` unknowns.

    In the other case, the algorithm performs transformation of the
    initial equation to an equivalent one, for which the system of
    algebraic equations has only `r` indeterminates. This method is
    quite sophisticated (in comparison with the naive one) and was
    invented together by Abramov, Bronstein and Petkovsek.

    It is possible to generalize the algorithm implemented here to
    the case of linear q-difference and differential equations.

    Lets say that we would like to compute `m`-th Bernoulli polynomial
    up to a constant. For this we can use `b(n+1) - b(n) = m n^{m-1}`
    recurrence, which has solution `b(n) = B_m + C`. For example:

    >>> from sympy import Symbol, rsolve_poly
    >>> n = Symbol('n', integer=True)

    >>> rsolve_poly([-1, 1], 4*n**3, n)
    C0 + n**4 - 2*n**3 + n**2

    References
    ==========

    .. [1] S. A. Abramov, M. Bronstein and M. Petkovsek, On polynomial
           solutions of linear operator equations, in: T. Levelt, ed.,
           Proc. ISSAC '95, ACM Press, New York, 1995, 290-296.

    .. [2] M. Petkovsek, Hypergeometric solutions of linear recurrences
           with polynomial coefficients, J. Symbolic Computation,
           14 (1992), 243-264.

    .. [3] M. Petkovsek, H. S. Wilf, D. Zeilberger, A = B, 1996.

    """
    f = sympify(f)

    if not f.is_polynomial(n):
        return None

    homogeneous = f.is_zero

    r = len(coeffs) - 1

    coeffs = [ Poly(coeff, n) for coeff in coeffs ]

    polys = [ Poly(0, n) ] * (r + 1)
    terms = [ (S.Zero, S.NegativeInfinity) ] *(r + 1)

    for i in xrange(0, r + 1):
        for j in xrange(i, r + 1):
            polys[i] += coeffs[j]*binomial(j, i)

        if not polys[i].is_zero:
            (exp,), coeff = polys[i].LT()
            terms[i] = (coeff, exp)

    d = b = terms[0][1]

    for i in xrange(1, r + 1):
        if terms[i][1] > d:
            d = terms[i][1]

        if terms[i][1] - i > b:
            b = terms[i][1] - i

    d, b = int(d), int(b)

    x = Dummy('x')

    degree_poly = S.Zero

    for i in xrange(0, r + 1):
        if terms[i][1] - i == b:
            degree_poly += terms[i][0]*FallingFactorial(x, i)

    nni_roots = roots(degree_poly, x, filter='Z',
        predicate=lambda r: r >= 0).keys()

    if nni_roots:
        N = [max(nni_roots)]
    else:
        N = []
#.........这里部分代码省略.........

def test_euler_failing():
    # depends on dummy variables being implemented https://github.com/sympy/sympy/issues/5665
    assert euler(2*n).rewrite(Sum) == I*Sum(Sum((-1)**_j*2**(-_k)*I**(-_k)*(-2*_j + _k)**(2*n + 1)*binomial(_k, _j)/_k, (_j, 0, _k)), (_k, 1, 2*n + 1))

def trigintegrate(f, x, conds='piecewise'):
    """Integrate f = Mul(trig) over x

       >>> from sympy import Symbol, sin, cos, tan, sec, csc, cot
       >>> from sympy.integrals.trigonometry import trigintegrate
       >>> from sympy.abc import x

       >>> trigintegrate(sin(x)*cos(x), x)
       sin(x)**2/2

       >>> trigintegrate(sin(x)**2, x)
       x/2 - sin(x)*cos(x)/2

       >>> trigintegrate(tan(x)*sec(x), x)
       1/cos(x)

       >>> trigintegrate(sin(x)*tan(x), x)
       -log(sin(x) - 1)/2 + log(sin(x) + 1)/2 - sin(x)

       http://en.wikibooks.org/wiki/Calculus/Integration_techniques

    See Also
    ========

    sympy.integrals.integrals.Integral.doit
    sympy.integrals.integrals.Integral
    """
    from sympy.integrals.integrals import integrate
    pat, a, n, m = _pat_sincos(x)

    f = f.rewrite('sincos')
    M = f.match(pat)

    if M is None:
        return

    n, m = M[n], M[m]
    if n is S.Zero and m is S.Zero:
        return x
    zz = x if n is S.Zero else S.Zero

    a = M[a]

    if n.is_odd or m.is_odd:
        u = _u
        n_, m_ = n.is_odd, m.is_odd

        # take smallest n or m -- to choose simplest substitution
        if n_ and m_:
            n_ = n_ and (n < m)  # NB: careful here, one of the
            m_ = m_ and not (n < m)  # conditions *must* be true

        #  n      m       u=C        (n-1)/2    m
        # S(x) * C(x) dx  --> -(1-u^2)       * u  du
        if n_:
            ff = -(1 - u**2)**((n - 1)/2) * u**m
            uu = cos(a*x)

        #  n      m       u=S   n         (m-1)/2
        # S(x) * C(x) dx  -->  u  * (1-u^2)       du
        elif m_:
            ff = u**n * (1 - u**2)**((m - 1)/2)
            uu = sin(a*x)

        fi = integrate(ff, u)  # XXX cyclic deps
        fx = fi.subs(u, uu)
        if conds == 'piecewise':
            return Piecewise((zz, Eq(a, 0)), (fx / a, True))
        return fx / a

    # n & m are both even
    #
    #               2k      2m                         2l       2l
    # we transform S (x) * C (x) into terms with only S (x) or C (x)
    #
    # example:
    #  100     4       100        2    2    100          4         2
    # S (x) * C (x) = S (x) * (1-S (x))  = S (x) * (1 + S (x) - 2*S (x))
    #
    #                  104       102     100
    #               = S (x) - 2*S (x) + S (x)
    #       2k
    # then S   is integrated with recursive formula

    # take largest n or m -- to choose simplest substitution
    n_ = (abs(n) > abs(m))
    m_ = (abs(m) > abs(n))
    res = S.Zero

    if n_:
        #  2k         2 k             i             2i
        # C   = (1 - S )  = sum(i, (-) * B(k, i) * S  )
        if m > 0:
            for i in range(0, m//2 + 1):
                res += ((-1)**i * binomial(m//2, i) *
                        _sin_pow_integrate(n + 2*i, x))

        elif m == 0:
            res = _sin_pow_integrate(n, x)
        else:
#.........这里部分代码省略.........

#.........这里部分代码省略.........
            # S(x) * C(x) dx  -->  u  * (1-u^2)       du
            elif m_:
                ff = u**n * (1-u**2)**((m-1)/2)
                uu = sin(a*x)

            fi= sympy.integrals.integrate(ff, u)    # XXX cyclic deps
            fx= fi.subs(u, uu)
            return fx / a

        # n & m are even
        else:
            #               2k      2m                         2l       2l
            # we transform S (x) * C (x) into terms with only S (x) or C (x)
            #
            # example:
            #  100     4       100        2    2    100          4         2
            # S (x) * C (x) = S (x) * (1-S (x))  = S (x) * (1 + S (x) - 2*S (x))
            #
            #                  104       102     100
            #               = S (x) - 2*S (x) + S (x)
            #       2k
            # then S   is integrated with recursive formula

            # take largest n or m -- to choose simplest substitution
            n_ =  (abs(n) > abs(m))
            m_ =  (abs(m) > abs(n))
            res = S.Zero

            if n_:
                #  2k       2 k             i            2i
                # C   = (1-S )  = sum(i, (-) * B(k,i) * S  )
                if m > 0 :
                    for i in range(0,m/2+1):
                        res += (-1)**i * binomial(m/2,i) * _sin_pow_integrate(n+2*i, x)

                elif m == 0:
                    res=_sin_pow_integrate(n,x)
                else:
                    # m < 0 , |n| > |m|
                    #  /                                                           /
                    # |                                                           |
                    # |    m       n            -1        m+1     n-1     n - 1   |     m+2     n-2
                    # | cos (x) sin (x) dx =  ________ cos (x) sin (x) + _______  |  cos (x) sin (x) dx
                    # |                                                           |
                    # |                         m + 1                     m + 1   |
                    #/                                                           /
                    #
                    #
                    res=Rational(-1,m+1)*cos(x)**(m+1)*sin(x)**(n-1) + Rational(n-1,m+1)*trigintegrate(cos(x)**(m+2)*sin(x)**(n-2),x)


            elif m_:
                #  2k        2 k            i            2i
                # S   = (1 -C ) = sum(i, (-) * B(k,i) * C  )
                if n > 0:
                    #      /                            /
                    #     |                            |
                    #     |    m       n               |    -m         n
                    #     | cos (x)*sin (x) dx  or     | cos (x) * sin (x) dx
                    #     |                            |
                    #    /                            /
                    #
                    #    |m| > |n| ; m,n >0 ; m,n belong to Z - {0}
                    #       n                                        2
                    #    sin (x) term is expanded here interms of cos (x), and then integrated.
                    for i in range(0,n/2+1):

def trigintegrate(f, x):
    """Integrate f = Mul(trig) over x

       >>> from sympy import Symbol, sin, cos
       >>> from sympy.integrals.trigonometry import trigintegrate
       >>> from sympy.abc import x

       >>> trigintegrate(sin(x)*cos(x), x)
       sin(x)**2/2

       >>> trigintegrate(sin(x)**2, x)
       x/2 - cos(x)*sin(x)/2

       http://en.wikibooks.org/wiki/Calculus/Further_integration_techniques
    """

    pat, a,n,m = _pat_sincos(x)
    ##m - cos
    ##n - sin

    M = f.match(pat)

    if M is None:
        return

    n, m = M[n], M[m]   # should always be there
    if n is S.Zero and m is S.Zero:
        return x

    a = M[a]

    if n.is_integer and m.is_integer:

        if n.is_odd or m.is_odd:
            u = _u
            n_, m_ = n.is_odd, m.is_odd

            # take smallest n or m -- to choose simplest substitution
            if n_ and m_:
                n_ = n_ and     (n < m)  # NB: careful here, one of the
                m_ = m_ and not (n < m)  #     conditions *must* be true

            #  n      m       u=C        (n-1)/2    m
            # S(x) * C(x) dx  --> -(1-u^2)       * u  du
            if n_:
                ff = -(1-u**2)**((n-1)/2) * u**m
                uu = cos(a*x)

            #  n      m       u=S   n         (m-1)/2
            # S(x) * C(x) dx  -->  u  * (1-u^2)       du
            elif m_:
                ff = u**n * (1-u**2)**((m-1)/2)
                uu = sin(a*x)

            fi= sympy.integrals.integrate(ff, u)    # XXX cyclic deps
            fx= fi.subs(u, uu)
            return fx / a

        # n & m are even
        else:
            #               2k      2m                         2l       2l
            # we transform S (x) * C (x) into terms with only S (x) or C (x)
            #
            # example:
            #  100     4       100        2    2    100          4         2
            # S (x) * C (x) = S (x) * (1-S (x))  = S (x) * (1 + S (x) - 2*S (x))
            #
            #                  104       102     100
            #               = S (x) - 2*S (x) + S (x)
            #       2k
            # then S   is integrated with recursive formula

            # take largest n or m -- to choose simplest substitution
            n_ =  (abs(n) > abs(m))
            m_ =  (abs(m) > abs(n))
            res = S.Zero

            if n_:
                #  2k       2 k             i            2i
                # C   = (1-S )  = sum(i, (-) * B(k,i) * S  )
                if m > 0 :
                    for i in range(0,m/2+1):
                        res += (-1)**i * binomial(m/2,i) * sin_pow_integrate(n+2*i, x)

                elif m == 0:
                    res=sin_pow_integrate(n,x)
                else:
                    # m < 0 , |n| > |m|
                    #  /                                                           /
                    # |                                                           |
                    # |    m       n            -1        m+1     n-1     n - 1   |     m+2     n-2
                    # | cos (x) sin (x) dx =  ________ cos (x) sin (x) + _______  |  cos (x) sin (x) dx
                    # |                                                           |
                    # |                         m + 1                     m + 1   |
                    #/                                                           /
                    #
                    #
                    res=Rational(-1,m+1)*cos(x)**(m+1)*sin(x)**(n-1) + Rational(n-1,m+1)*trigintegrate(cos(x)**(m+2)*sin(x)**(n-2),x)


#.........这里部分代码省略.........

def trigintegrate(f, x):
    """Integrate f = Mul(trig) over x

       >>> from sympy import Symbol, sin, cos
       >>> from sympy.integrals.trigonometry import trigintegrate
       >>> x = Symbol('x')

       >>> trigintegrate(sin(x)*cos(x), x)
       1/2*sin(x)**2

       >>> trigintegrate(sin(x)**2, x)
       x/2 - cos(x)*sin(x)/2

       http://en.wikibooks.org/wiki/Calculus/Further_integration_techniques
    """

    pat, a,n,m = _pat_sincos(x)

    M = f.match(pat)
    if M is None:
        return

    n, m = M[n], M[m]   # should always be there
    if n is S.Zero and m is S.Zero:
        return x

    a = M[a]

    if n.is_integer and n.is_integer:

        if n.is_odd or m.is_odd:
            u = _u
            n_, m_ = n.is_odd, m.is_odd

            # take smallest n or m -- to choose simplest substitution
            if n_ and m_:
               n_ = n_ and     (n < m)  # NB: careful here, one of the
               m_ = m_ and not (n < m)  #     conditions *must* be true

            #  n      m       u=C        (n-1)/2    m
            # S(x) * C(x) dx  --> -(1-u^2)       * u  du
            if n_:
                ff = -(1-u**2)**((n-1)/2) * u**m
                uu = cos(a*x)

            #  n      m       u=S   n         (m-1)/2
            # S(x) * C(x) dx  -->  u  * (1-u^2)       du
            elif m_:
                ff = u**n * (1-u**2)**((m-1)/2)
                uu = sin(a*x)

            fi= sympy.integrals.integrate(ff, u)    # XXX cyclic deps
            fx= fi.subs(u, uu)
            return fx / a

        # n & m are even
        else:
            #               2k      2m                         2l       2l
            # we transform S (x) * C (x) into terms with only S (x) or C (x)
            #
            # example:
            #  100     4       100        2    2    100          4         2
            # S (x) * C (x) = S (x) * (1-S (x))  = S (x) * (1 + S (x) - 2*S (x))
            #
            #                  104       102     100
            #               = S (x) - 2*S (x) + S (x)
            #       2k
            # then S   is integrated with recursive formula

            # take largest n or m -- to choose simplest substitution
            n_ =     (n > m)    # NB: careful here, one of the
            m_ = not (n > m)    #     conditions *must* be true

            res = S.Zero

            if n_:
                #  2k       2 k             i            2i
                # C   = (1-S )  = sum(i, (-) * B(k,i) * S  )
                for i in range(0,m/2+1):
                    res += (-1)**i * binomial(m/2,i) * Sin_2k_integrate(n/2+i, x)

            elif m_:
                #  2k        2 k            i            2i
                # S   = (1 -C ) = sum(i, (-) * B(k,i) * C  )
                for i in range(0,n/2+1):
                    res += (-1)**i * binomial(n/2,i) * Cos_2k_integrate(m/2+i, x)

            return res.subs(x, a*x) / a