Python arith.binomial函数代码示例

 def cardinality(self):
     """
     EXAMPLES::
     
         sage: IntegerVectors(3,3, min_part=1).cardinality()
         1
         sage: IntegerVectors(5,3, min_part=1).cardinality()
         6
         sage: IntegerVectors(13, 4, min_part=2, max_part=4).cardinality()
         16
     """
     if not self.constraints:
         if self.n >= 0:
             return binomial(self.n+self.k-1,self.n)
         else:
             return 0
     else:
         if len(self.constraints) == 1 and 'max_part' in self.constraints and self.constraints['max_part'] != infinity:
             m = self.constraints['max_part']
             if m >= self.n:
                 return binomial(self.n+self.k-1,self.n)
             else: #do by inclusion / exclusion on the number
                   #i of parts greater than m
                 return sum( [(-1)**i * binomial(self.n+self.k-1-i*(m+1), self.k-1)*binomial(self.k,i) for i in range(0, self.n/(m+1)+1)])
         else:
             return len(self.list())

    def cardinality(self):
        r"""
        Return the number of Baxter permutations of size ``self._n``.

        For any positive integer `n`, the number of Baxter
        permutations of size `n` equals

        .. MATH::

            \sum_{k=1}^n \dfrac
            {\binom{n+1}{k-1} \binom{n+1}{k} \binom{n+1}{k+1}}
            {\binom{n+1}{1} \binom{n+1}{2}} .

        This is :oeis:`A001181`.

        EXAMPLES::

            sage: [BaxterPermutations(n).cardinality() for n in xrange(13)]
            [1, 1, 2, 6, 22, 92, 422, 2074, 10754, 58202, 326240, 1882960, 11140560]

            sage: BaxterPermutations(3r).cardinality()
            6
            sage: parent(_)
            Integer Ring
        """
        if self._n == 0:
            return 1
        from sage.rings.arith import binomial
        return sum((binomial(self._n + 1, k) *
                    binomial(self._n + 1, k + 1) *
                    binomial(self._n + 1, k + 2)) //
                   ((self._n + 1) * binomial(self._n + 1, 2))
                   for k in xrange(self._n))

    def parameters(self, t=None):
        """
        Returns `(t,v,k,lambda)`. Does not check if the input is a block
        design.

        INPUT:

        - ``t`` -- `t` such that the design is a `t`-design.

        EXAMPLES::

            sage: from sage.combinat.designs.block_design import BlockDesign
            sage: BD = BlockDesign(7,[[0,1,2],[0,3,4],[0,5,6],[1,3,5],[1,4,6],[2,3,6],[2,4,5]], name="FanoPlane")
            sage: BD.parameters(t=2)
            (2, 7, 3, 1)
            sage: BD.parameters(t=3)
            (3, 7, 3, 0)
        """
        if t is None:
            from sage.misc.superseded import deprecation
            deprecation(15664, "the 't' argument will become mandatory soon. 2"+
                        " is used when none is provided.")
            t = 2

        v = len(self.points())
        blks = self.blocks()
        k = len(blks[int(0)])
        b = len(blks)
        #A = self.incidence_matrix()
        #r = sum(A.rows()[0])
        lmbda = int(b/(binomial(v, t)/binomial(k, t)))
        return (t, v, k, lmbda)

    def __init__(self, R, elements):
        """
        Initialize ``self``.

        EXAMPLES::

            sage: R.<x,y,z> = QQ[]
            sage: K = KoszulComplex(R, [x,y])
            sage: TestSuite(K).run()
        """
        # Generate the differentials
        self._elements = elements
        n = len(elements)
        I = range(n)
        diff = {}
        zero = R.zero()
        for i in I:
            M = matrix(R, binomial(n,i), binomial(n,i+1), zero)
            j = 0
            for comb in itertools.combinations(I, i+1):
                for k,val in enumerate(comb):
                    r = rank(comb[:k] + comb[k+1:], n, False)
                    M[r,j] = (-1)**k * elements[val]
                j += 1
            M.set_immutable()
            diff[i+1] = M
        diff[0] = matrix(R, 0, 1, zero)
        diff[0].set_immutable()
        diff[n+1] = matrix(R, 1, 0, zero)
        diff[n+1].set_immutable()
        ChainComplex_class.__init__(self, ZZ, ZZ(-1), R, diff)

 def _basic_integral(self, a, j, twist=None):
     r"""
     Computes the integral
     
         .. MATH::
            
            \int_{a+p\ZZ_p}(z-\omega(a))^jd\mu_\chi.
     
     If ``twist`` is ``None``, `\\chi` is the trivial character. Otherwise, ``twist`` can be a primitive quadratic character of conductor prime to `p`.
     """
     #is this the negative of what we want?
     #if Phis is fixed for this p-adic L-function, we should make this method cached
     p = self._Phis.parent().prime()
     if twist is None:
         pass
     elif twist in ZZ:
         twist = kronecker_character(twist)
         if twist.is_trivial():
             twist = None
         else:
             D = twist.level()
             assert(D.gcd(p) == 1)
     else:
         if twist.is_trivial():
             twist = None
         else:
             assert((twist**2).is_trivial())
             twist = twist.primitive_character()
             D = twist.level()
             assert(D.gcd(p) == 1)
     
     onDa = self._on_Da(a, twist)#self._Phis(Da)
     aminusat = a - self._Phis.parent().base_ring().base_ring().teichmuller(a)
     #aminusat = a - self._coefficient_ring.base_ring().teichmuller(a)
     try:
         ap = self._ap
     except AttributeError:
         self._ap = self._Phis.Tq_eigenvalue(p) #catch exception if not eigensymbol
         ap = self._ap
     if not twist is None:
         ap *= twist(p)
     if j == 0:
         return (~ap) * onDa.moment(0)
     if a == 1:
         #aminusat is 0, so only the j=r term is non-zero
         return (~ap) * (p ** j) * onDa.moment(j)
     #print "j =", j, "a = ", a
     ans = onDa.moment(0) * (aminusat ** j)
     #ans = onDa.moment(0)
     #print "\tr =", 0, " ans =", ans
     for r in range(1, j+1):
         if r == j:
             ans += binomial(j, r) * (p ** r) * onDa.moment(r)
         else:
             ans += binomial(j, r) * (aminusat ** (j - r)) * (p ** r) * onDa.moment(r)
         #print "\tr =", r, " ans =", ans
     #print " "
     return (~ap) * ans

def from_rank(r, n, k):
    """
    Returns the combination of rank r in the subsets of range(n) of
    size k when listed in lexicographic order.

    The algorithm used is based on combinadics and James McCaffrey's
    MSDN article. See: http://en.wikipedia.org/wiki/Combinadic

    EXAMPLES::

        sage: import sage.combinat.choose_nk as choose_nk
        sage: choose_nk.from_rank(0,3,0)
        ()
        sage: choose_nk.from_rank(0,3,1)
        (0,)
        sage: choose_nk.from_rank(1,3,1)
        (1,)
        sage: choose_nk.from_rank(2,3,1)
        (2,)
        sage: choose_nk.from_rank(0,3,2)
        (0, 1)
        sage: choose_nk.from_rank(1,3,2)
        (0, 2)
        sage: choose_nk.from_rank(2,3,2)
        (1, 2)
        sage: choose_nk.from_rank(0,3,3)
        (0, 1, 2)
    """
    if k < 0:
        raise ValueError("k must be > 0")
    if k > n:
        raise ValueError("k must be <= n")

    a = n
    b = k
    x = binomial(n, k) - 1 - r  # x is the 'dual' of m
    comb = [None] * k

    for i in xrange(k):
        comb[i] = _comb_largest(a, b, x)
        x = x - binomial(comb[i], b)
        a = comb[i]
        b = b - 1

    for i in xrange(k):
        comb[i] = (n - 1) - comb[i]

    return tuple(comb)

    def rank(self, sub):
        """
        Returns the rank of sub as a subset of s.

        EXAMPLES::

            sage: Subsets(3).rank([])
            0
            sage: Subsets(3).rank([1,2])
            4
            sage: Subsets(3).rank([1,2,3])
            7
            sage: Subsets(3).rank([2,3,4]) == None
            True
        """
        subset = Set(sub)
        lset = __builtin__.list(self.s)
        lsubset = __builtin__.list(subset)

        try:
            index_list = sorted(map(lambda x: lset.index(x), lsubset))
        except ValueError:
            return None

        n = len(self.s)
        r = 0

        for i in range(len(index_list)):
            r += binomial(n,i)
        return r + choose_nk.rank(index_list,n)

    def rank(self, x):
        """
        Returns the position of a given element.

        INPUT:

        - ``x`` - a list with ``sum(x) == n`` and ``len(x) == k``

        TESTS::

            sage: IV = IntegerVectors(4,5) 
            sage: range(IV.cardinality()) == [IV.rank(x) for x in IV]
            True
        """

        if x not in self:
            raise ValueError("argument is not a member of IntegerVectors(%d,%d)" % (self.n, self.k))

        n = self.n
        k = self.k

        r = 0
        for i in range(k-1):
          k -= 1
          n -= x[i]
          r += binomial(k+n-1,k)

        return r

    def cardinality(self):
        r"""
        Return the number of words in the shuffle product
        of ``w1`` and ``w2``.

        This is understood as a multiset cardinality, not as a
        set cardinality; it does not count the distinct words only.

        It is given by `\binom{l_1+l_2}{l_1}`, where `l_1` is the
        length of ``w1`` and where `l_2` is the length of ``w2``.

        EXAMPLES::

            sage: from sage.combinat.words.shuffle_product import ShuffleProduct_w1w2
            sage: w, u = map(Words("abcd"), ["ab", "cd"])
            sage: S = ShuffleProduct_w1w2(w,u)
            sage: S.cardinality()
            6

            sage: w, u = map(Words("ab"), ["ab", "ab"])
            sage: S = ShuffleProduct_w1w2(w,u)
            sage: S.cardinality()
            6
        """
        return binomial(self._w1.length()+self._w2.length(), self._w1.length())

    def cardinality(self):
        """
        EXAMPLES::

            sage: Subsets(Set([1,2,3]), 2).cardinality()
            3
            sage: Subsets([1,2,3,3], 2).cardinality()
            3
            sage: Subsets([1,2,3], 1).cardinality()
            3
            sage: Subsets([1,2,3], 3).cardinality()
            1
            sage: Subsets([1,2,3], 0).cardinality()
            1
            sage: Subsets([1,2,3], 4).cardinality()
            0
            sage: Subsets(3,2).cardinality()
            3
            sage: Subsets(3,4).cardinality()
            0
        """
        if self.k not in range(len(self.s)+1):
            return 0
        else:
            return binomial(len(self.s),self.k)

    def __getitem__(self, key):
        """
        EXAMPLES::

            sage: a, b, c, q, z = var('a b c q z')
            sage: bhs = BasicHypergeometricSeries([a, b], [c], q, z)
            sage: for i in range(4): print bhs[i]
            1
            (b - 1)*(a - 1)*z/((q - 1)*(c - 1))
            (b - 1)*(a - 1)*(b*q - 1)*(a*q - 1)*z^2/((q - 1)*(c - 1)*(q^2 - 1)*(c*q - 1))
            (b - 1)*(a - 1)*(b*q - 1)*(a*q - 1)*(b*q^2 - 1)*(a*q^2 - 1)*z^3/((q - 1)*(c - 1)*(q^2 - 1)*(c*q - 1)*(q^3 - 1)*(c*q^2 - 1))
        
        """
        if key >= 0:
            j, k = len(self.list_a), len(self.list_b)

            nominator = qPochhammerSymbol(self.list_a, self.q, key).evaluate()
            if nominator == 0:
                return 0
            denominator = (
                qPochhammerSymbol(self.list_b, self.q, key).evaluate()
                * qPochhammerSymbol(self.q, self.q, key).evaluate()
            )

            return nominator / denominator * ((-1) ** key * self.q ** (binomial(key, 2))) ** (1 + k - j) * self.z ** key
        else:
            return 0

def log_gamma_binomial(p,gamma,z,n,M):
    r"""
    Returns the list of coefficients in the power series
    expansion (up to precision `M`) of `{\log_p(z)/\log_p(\gamma) \choose n}`

    INPUT:

        - ``p`` --  prime
        - ``gamma`` -- topological generator e.g., `1+p`
        - ``z`` -- variable
        - ``n`` -- nonnegative integer
        - ``M`` -- precision

    OUTPUT:

    The list of coefficients in the power series expansion of
    `{\log_p(z)/\log_p(\gamma) \choose n}`

    EXAMPLES:

        sage: R.<z> = QQ['z']
        sage: from sage.modular.pollack_stevens.padic_lseries import log_gamma_binomial
        sage: log_gamma_binomial(5,1+5,z,2,4)
        [0, -3/205, 651/84050, -223/42025]
        sage: log_gamma_binomial(5,1+5,z,3,4)
        [0, 2/205, -223/42025, 95228/25845375]
    """
    L = sum([ZZ(-1)**j / j*z**j for j in range (1,M)]) #log_p(1+z)
    loggam = L / (L(gamma - 1))                  #log_{gamma}(1+z)= log_p(1+z)/log_p(gamma)
    return z.parent()(binomial(loggam,n)).truncate(M).list()

def upper_bound(min_length, max_length, floor, ceiling, min_slope, max_slope):
    """
    Compute a coarse upper bound on the size of a vector satisfying the
    constraints.

    TESTS::

        sage: import sage.combinat.integer_list as integer_list
        sage: f = lambda x: lambda i: x
        sage: integer_list.upper_bound(0,4,f(0), f(1),-infinity,infinity)
        4
        sage: integer_list.upper_bound(0, infinity, f(0), f(1), -infinity, infinity)
        inf
        sage: integer_list.upper_bound(0, infinity, f(0), f(1), -infinity, -1)
        1
        sage: integer_list.upper_bound(0, infinity, f(0), f(5), -infinity, -1)
        15
        sage: integer_list.upper_bound(0, infinity, f(0), f(5), -infinity, -2)
        9
    """
    from sage.functions.all import floor as flr
    if max_length < float('inf'):
        return sum( [ ceiling(j) for j in range(max_length)] )
    elif max_slope < 0 and ceiling(1) < float('inf'):
        maxl = flr(-ceiling(1)/max_slope)
        return ceiling(1)*(maxl+1) + binomial(maxl+1,2)*max_slope
    #FIXME: only checking the first 10000 values, but that should generally
    #be enough
    elif [ceiling(j) for j in range(10000)] == [0]*10000:
        return 0
    else:
        return float('inf')

def Krawtchouk(n,q,l,i):
    """
    Compute ``K^{n,q}_l(i)``, the Krawtchouk polynomial:
    see :wikipedia:`Kravchuk_polynomials`.
    It is given by

    .. math::

        K^{n,q}_l(i)=\sum_{j=0}^l (-1)^j(q-1)^{(l-j)}{i \choose j}{n-i \choose l-j}

    EXAMPLES::

        sage: Krawtchouk(24,2,5,4)
        2224
        sage: Krawtchouk(12300,4,5,6)
        567785569973042442072

    """
    from sage.rings.arith import binomial
    # Use the expression in equation (55) of MacWilliams & Sloane, pg 151
    # We write jth term = some_factor * (j-1)th term
    kraw = jth_term = (q-1)**l * binomial(n, l) # j=0
    for j in range(1,l+1):
        jth_term *= -q*(l-j+1)*(i-j+1)/((q-1)*j*(n-j+1))
        kraw += jth_term
    return kraw

def loggam_binom(p,gam,z,n,M):
    r"""
    Returns the list of coefficients in the power series
    expansion (up to precision `M`) of `{\log_p(z)/\log_p(\gamma) \choose n}`

    INPUT:

        - ``p`` --  prime
        - ``gam`` -- topological generator e.g., `1+p`
        - ``z`` -- variable
        - ``n`` -- nonnegative integer
        - ``M`` -- precision

    OUTPUT:

    The list of coefficients in the power series expansion of 
    `{\log_p(z)/\log_p(\gamma) \choose n}`
    
    EXAMPLES:

        sage: R.<z> = QQ['z']  
        sage: loggam_binom(5,1+5,z,2,4)
        [0, -3/205, 651/84050, -223/42025]
        sage: loggam_binom(5,1+5,z,3,4)
        [0, 2/205, -223/42025, 95228/25845375]
    """
    L = logp(p,z,M)
    logpgam = L.substitute(z = (gam-1)) #log base p of gamma
    loggam = L/logpgam                  #log base gamma 
    return z.parent()(binomial(loggam,n)).truncate(M).list()