Combinatorial¶

This module implements various combinatorial functions.

bell¶

class sympy.functions.combinatorial.numbers.bell[source]¶

Bell numbers / Bell polynomials

The Bell numbers satisfy \(B_0 = 1\) and

\[B_n = \sum_{k=0}^{n-1} \binom{n-1}{k} B_k.\]

They are also given by:

\[B_n = \frac{1}{e} \sum_{k=0}^{\infty} \frac{k^n}{k!}.\]

The Bell polynomials are given by \(B_0(x) = 1\) and

\[B_n(x) = x \sum_{k=1}^{n-1} \binom{n-1}{k-1} B_{k-1}(x).\]

The second kind of Bell polynomials (are sometimes called “partial” Bell polynomials or incomplete Bell polynomials) are defined as

\[B_{n,k}(x_1, x_2,\dotsc x_{n-k+1}) = \sum_{j_1+j_2+j_2+\dotsb=k \atop j_1+2j_2+3j_2+\dotsb=n} \frac{n!}{j_1!j_2!\dotsb j_{n-k+1}!} \left(\frac{x_1}{1!} \right)^{j_1} \left(\frac{x_2}{2!} \right)^{j_2} \dotsb \left(\frac{x_{n-k+1}}{(n-k+1)!} \right) ^{j_{n-k+1}}.\]

bell(n) gives the \(n^{th}\) Bell number, \(B_n\).
bell(n, x) gives the \(n^{th}\) Bell polynomial, \(B_n(x)\).
bell(n, k, (x1, x2, ...)) gives Bell polynomials of the second kind, \(B_{n,k}(x_1, x_2, \dotsc, x_{n-k+1})\).

Notes

Not to be confused with Bernoulli numbers and Bernoulli polynomials, which use the same notation.

Examples

>>> from sympy import bell, Symbol, symbols

>>> [bell(n) for n in range(11)]
[1, 1, 2, 5, 15, 52, 203, 877, 4140, 21147, 115975]
>>> bell(30)
846749014511809332450147
>>> bell(4, Symbol('t'))
t**4 + 6*t**3 + 7*t**2 + t
>>> bell(6, 2, symbols('x:6')[1:])
6*x1*x5 + 15*x2*x4 + 10*x3**2

References

R138: https://en.wikipedia.org/wiki/Bell_number
R139: http://mathworld.wolfram.com/BellNumber.html
R140: http://mathworld.wolfram.com/BellPolynomial.html

bernoulli¶

class sympy.functions.combinatorial.numbers.bernoulli[source]¶

Bernoulli numbers / Bernoulli polynomials

The Bernoulli numbers are a sequence of rational numbers defined by \(B_0 = 1\) and the recursive relation (\(n > 0\)):

\[0 = \sum_{k=0}^n \binom{n+1}{k} B_k\]

They are also commonly defined by their exponential generating function, which is \(\frac{x}{e^x - 1}\). For odd indices > 1, the Bernoulli numbers are zero.

The Bernoulli polynomials satisfy the analogous formula:

\[B_n(x) = \sum_{k=0}^n \binom{n}{k} B_k x^{n-k}\]

Bernoulli numbers and Bernoulli polynomials are related as \(B_n(0) = B_n\).

We compute Bernoulli numbers using Ramanujan’s formula:

\[B_n = \frac{A(n) - S(n)}{\binom{n+3}{n}}\]

where:

\[\begin{split}A(n) = \begin{cases} \frac{n+3}{3} & n \equiv 0\ \text{or}\ 2 \pmod{6} \\ -\frac{n+3}{6} & n \equiv 4 \pmod{6} \end{cases}\end{split}\]

and:

\[S(n) = \sum_{k=1}^{[n/6]} \binom{n+3}{n-6k} B_{n-6k}\]

This formula is similar to the sum given in the definition, but cuts 2/3 of the terms. For Bernoulli polynomials, we use the formula in the definition.

bernoulli(n) gives the nth Bernoulli number, \(B_n\)
bernoulli(n, x) gives the nth Bernoulli polynomial in \(x\), \(B_n(x)\)

Examples

>>> from sympy import bernoulli

>>> [bernoulli(n) for n in range(11)]
[1, -1/2, 1/6, 0, -1/30, 0, 1/42, 0, -1/30, 0, 5/66]
>>> bernoulli(1000001)
0

References

R141: https://en.wikipedia.org/wiki/Bernoulli_number
R142: https://en.wikipedia.org/wiki/Bernoulli_polynomial
R143: http://mathworld.wolfram.com/BernoulliNumber.html
R144: http://mathworld.wolfram.com/BernoulliPolynomial.html

binomial¶

class sympy.functions.combinatorial.factorials.binomial[source]¶

Implementation of the binomial coefficient. It can be defined in two ways depending on its desired interpretation:

\[\binom{n}{k} = \frac{n!}{k!(n-k)!}\ \text{or}\ \binom{n}{k} = \frac{ff(n, k)}{k!}\]

First, in a strict combinatorial sense it defines the number of ways we can choose \(k\) elements from a set of \(n\) elements. In this case both arguments are nonnegative integers and binomial is computed using an efficient algorithm based on prime factorization.

The other definition is generalization for arbitrary \(n\), however \(k\) must also be nonnegative. This case is very useful when evaluating summations.

For the sake of convenience for negative integer \(k\) this function will return zero no matter what valued is the other argument.

To expand the binomial when \(n\) is a symbol, use either expand_func() or expand(func=True). The former will keep the polynomial in factored form while the latter will expand the polynomial itself. See examples for details.

Examples

>>> from sympy import Symbol, Rational, binomial, expand_func
>>> n = Symbol('n', integer=True, positive=True)

>>> binomial(15, 8)
6435

>>> binomial(n, -1)
0

Rows of Pascal’s triangle can be generated with the binomial function:

>>> for N in range(8):
...     print([binomial(N, i) for i in range(N + 1)])
...
[1]
[1, 1]
[1, 2, 1]
[1, 3, 3, 1]
[1, 4, 6, 4, 1]
[1, 5, 10, 10, 5, 1]
[1, 6, 15, 20, 15, 6, 1]
[1, 7, 21, 35, 35, 21, 7, 1]

As can a given diagonal, e.g. the 4th diagonal:

>>> N = -4
>>> [binomial(N, i) for i in range(1 - N)]
[1, -4, 10, -20, 35]

>>> binomial(Rational(5, 4), 3)
-5/128
>>> binomial(Rational(-5, 4), 3)
-195/128

>>> binomial(n, 3)
binomial(n, 3)

>>> binomial(n, 3).expand(func=True)
n**3/6 - n**2/2 + n/3

>>> expand_func(binomial(n, 3))
n*(n - 2)*(n - 1)/6

References

R145: https://www.johndcook.com/blog/binomial_coefficients/

catalan¶

class sympy.functions.combinatorial.numbers.catalan[source]¶

Catalan numbers

The \(n^{th}\) catalan number is given by:

\[C_n = \frac{1}{n+1} \binom{2n}{n}\]

catalan(n) gives the \(n^{th}\) Catalan number, \(C_n\)

Examples

>>> from sympy import (Symbol, binomial, gamma, hyper, polygamma,
...             catalan, diff, combsimp, Rational, I)

>>> [catalan(i) for i in range(1,10)]
[1, 2, 5, 14, 42, 132, 429, 1430, 4862]

>>> n = Symbol("n", integer=True)

>>> catalan(n)
catalan(n)

Catalan numbers can be transformed into several other, identical expressions involving other mathematical functions

>>> catalan(n).rewrite(binomial)
binomial(2*n, n)/(n + 1)

>>> catalan(n).rewrite(gamma)
4**n*gamma(n + 1/2)/(sqrt(pi)*gamma(n + 2))

>>> catalan(n).rewrite(hyper)
hyper((1 - n, -n), (2,), 1)

For some non-integer values of n we can get closed form expressions by rewriting in terms of gamma functions:

>>> catalan(Rational(1,2)).rewrite(gamma)
8/(3*pi)

We can differentiate the Catalan numbers C(n) interpreted as a continuous real function in n:

>>> diff(catalan(n), n)
(polygamma(0, n + 1/2) - polygamma(0, n + 2) + log(4))*catalan(n)

As a more advanced example consider the following ratio between consecutive numbers:

>>> combsimp((catalan(n + 1)/catalan(n)).rewrite(binomial))
2*(2*n + 1)/(n + 2)

The Catalan numbers can be generalized to complex numbers:

>>> catalan(I).rewrite(gamma)
4**I*gamma(1/2 + I)/(sqrt(pi)*gamma(2 + I))

and evaluated with arbitrary precision:

>>> catalan(I).evalf(20)
0.39764993382373624267 - 0.020884341620842555705*I

References

R146: https://en.wikipedia.org/wiki/Catalan_number
R147: http://mathworld.wolfram.com/CatalanNumber.html
R148: http://functions.wolfram.com/GammaBetaErf/CatalanNumber/
R149: http://geometer.org/mathcircles/catalan.pdf

euler¶

class sympy.functions.combinatorial.numbers.euler[source]¶

Euler numbers / Euler polynomials

The Euler numbers are given by:

\[E_{2n} = I \sum_{k=1}^{2n+1} \sum_{j=0}^k \binom{k}{j} \frac{(-1)^j (k-2j)^{2n+1}}{2^k I^k k}\]

\[E_{2n+1} = 0\]

Euler numbers and Euler polynomials are related by

\[E_n = 2^n E_n\left(\frac{1}{2}\right).\]

We compute symbolic Euler polynomials using [R154]

\[E_n(x) = \sum_{k=0}^n \binom{n}{k} \frac{E_k}{2^k} \left(x - \frac{1}{2}\right)^{n-k}.\]

However, numerical evaluation of the Euler polynomial is computed more efficiently (and more accurately) using the mpmath library.

euler(n) gives the \(n^{th}\) Euler number, \(E_n\).
euler(n, x) gives the \(n^{th}\) Euler polynomial, \(E_n(x)\).

Examples

>>> from sympy import Symbol, S
>>> from sympy.functions import euler
>>> [euler(n) for n in range(10)]
[1, 0, -1, 0, 5, 0, -61, 0, 1385, 0]
>>> n = Symbol("n")
>>> euler(n + 2*n)
euler(3*n)

>>> x = Symbol("x")
>>> euler(n, x)
euler(n, x)

>>> euler(0, x)
1
>>> euler(1, x)
x - 1/2
>>> euler(2, x)
x**2 - x
>>> euler(3, x)
x**3 - 3*x**2/2 + 1/4
>>> euler(4, x)
x**4 - 2*x**3 + x

>>> euler(12, S.Half)
2702765/4096
>>> euler(12)
2702765

References

R150: https://en.wikipedia.org/wiki/Euler_numbers
R151: http://mathworld.wolfram.com/EulerNumber.html
R152: https://en.wikipedia.org/wiki/Alternating_permutation
R153: http://mathworld.wolfram.com/AlternatingPermutation.html
R154(1,2): http://dlmf.nist.gov/24.2#ii

factorial¶

class sympy.functions.combinatorial.factorials.factorial[source]¶

Implementation of factorial function over nonnegative integers. By convention (consistent with the gamma function and the binomial coefficients), factorial of a negative integer is complex infinity.

The factorial is very important in combinatorics where it gives the number of ways in which \(n\) objects can be permuted. It also arises in calculus, probability, number theory, etc.

There is strict relation of factorial with gamma function. In fact \(n! = gamma(n+1)\) for nonnegative integers. Rewrite of this kind is very useful in case of combinatorial simplification.

Computation of the factorial is done using two algorithms. For small arguments a precomputed look up table is used. However for bigger input algorithm Prime-Swing is used. It is the fastest algorithm known and computes \(n!\) via prime factorization of special class of numbers, called here the ‘Swing Numbers’.

Examples

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

>>> factorial(0)
1

>>> factorial(7)
5040

>>> factorial(-2)
zoo

>>> factorial(n)
factorial(n)

>>> factorial(2*n)
factorial(2*n)

>>> factorial(S(1)/2)
factorial(1/2)

subfactorial¶

class sympy.functions.combinatorial.factorials.subfactorial[source]¶

The subfactorial counts the derangements of n items and is defined for non-negative integers as:

\[\begin{split}!n = \begin{cases} 1 & n = 0 \\ 0 & n = 1 \\ (n-1)(!(n-1) + !(n-2)) & n > 1 \end{cases}\end{split}\]

It can also be written as int(round(n!/exp(1))) but the recursive definition with caching is implemented for this function.

An interesting analytic expression is the following [R156]

\[!x = \Gamma(x + 1, -1)/e\]

which is valid for non-negative integers \(x\). The above formula is not very useful incase of non-integers. \(\Gamma(x + 1, -1)\) is single-valued only for integral arguments \(x\), elsewhere on the positive real axis it has an infinite number of branches none of which are real.

Examples

>>> from sympy import subfactorial
>>> from sympy.abc import n
>>> subfactorial(n + 1)
subfactorial(n + 1)
>>> subfactorial(5)
44

References

R155: https://en.wikipedia.org/wiki/Subfactorial
R156(1,2): http://mathworld.wolfram.com/Subfactorial.html

factorial2 / double factorial¶

class sympy.functions.combinatorial.factorials.factorial2[source]¶

The double factorial \(n!!\), not to be confused with \((n!)!\)

The double factorial is defined for nonnegative integers and for odd negative integers as:

\[\begin{split}n!! = \begin{cases} 1 & n = 0 \\ n(n-2)(n-4) \cdots 1 & n\ \text{positive odd} \\ n(n-2)(n-4) \cdots 2 & n\ \text{positive even} \\ (n+2)!!/(n+2) & n\ \text{negative odd} \end{cases}\end{split}\]

Examples

>>> from sympy import factorial2, var
>>> var('n')
n
>>> factorial2(n + 1)
factorial2(n + 1)
>>> factorial2(5)
15
>>> factorial2(-1)
1
>>> factorial2(-5)
1/3

References

R157: https://en.wikipedia.org/wiki/Double_factorial

FallingFactorial¶

class sympy.functions.combinatorial.factorials.FallingFactorial[source]¶

Falling factorial (related to rising factorial) is a double valued function arising in concrete mathematics, hypergeometric functions and series expansions. It is defined by

\[ff(x,k) = x \cdot (x-1) \cdots (x-k+1)\]

where \(x\) can be arbitrary expression and \(k\) is an integer. For more information check “Concrete mathematics” by Graham, pp. 66 or visit http://mathworld.wolfram.com/FallingFactorial.html page.

When \(x\) is a Poly instance of degree >= 1 with single variable, \(ff(x,k) = x(y) \cdot x(y-1) \cdots x(y-k+1)\), where \(y\) is the variable of \(x\). This is as described in Peter Paule, “Greatest Factorial Factorization and Symbolic Summation”, Journal of Symbolic Computation, vol. 20, pp. 235-268, 1995.

>>> from sympy import ff, factorial, rf, gamma, polygamma, binomial, symbols, Poly
>>> from sympy.abc import x, k
>>> n, m = symbols('n m', integer=True)
>>> ff(x, 0)
1
>>> ff(5, 5)
120
>>> ff(x, 5) == x*(x-1)*(x-2)*(x-3)*(x-4)
True
>>> ff(Poly(x**2, x), 2)
Poly(x**4 - 2*x**3 + x**2, x, domain='ZZ')
>>> ff(n, n)
factorial(n)

Rewrite

>>> ff(x, k).rewrite(gamma)
(-1)**k*gamma(k - x)/gamma(-x)
>>> ff(x, k).rewrite(rf)
RisingFactorial(-k + x + 1, k)
>>> ff(x, m).rewrite(binomial)
binomial(x, m)*factorial(m)
>>> ff(n, m).rewrite(factorial)
factorial(n)/factorial(-m + n)

References

R158: http://mathworld.wolfram.com/FallingFactorial.html

fibonacci¶

class sympy.functions.combinatorial.numbers.fibonacci[source]¶

Fibonacci numbers / Fibonacci polynomials

The Fibonacci numbers are the integer sequence defined by the initial terms \(F_0 = 0\), \(F_1 = 1\) and the two-term recurrence relation \(F_n = F_{n-1} + F_{n-2}\). This definition extended to arbitrary real and complex arguments using the formula

\[F_z = \frac{\phi^z - \cos(\pi z) \phi^{-z}}{\sqrt 5}\]

The Fibonacci polynomials are defined by \(F_1(x) = 1\), \(F_2(x) = x\), and \(F_n(x) = x*F_{n-1}(x) + F_{n-2}(x)\) for \(n > 2\). For all positive integers \(n\), \(F_n(1) = F_n\).

fibonacci(n) gives the \(n^{th}\) Fibonacci number, \(F_n\)
fibonacci(n, x) gives the \(n^{th}\) Fibonacci polynomial in \(x\), \(F_n(x)\)

Examples

>>> from sympy import fibonacci, Symbol

>>> [fibonacci(x) for x in range(11)]
[0, 1, 1, 2, 3, 5, 8, 13, 21, 34, 55]
>>> fibonacci(5, Symbol('t'))
t**4 + 3*t**2 + 1

References

R159: https://en.wikipedia.org/wiki/Fibonacci_number
R160: http://mathworld.wolfram.com/FibonacciNumber.html

tribonacci¶

class sympy.functions.combinatorial.numbers.tribonacci[source]¶

Tribonacci numbers / Tribonacci polynomials

The Tribonacci numbers are the integer sequence defined by the initial terms \(T_0 = 0\), \(T_1 = 1\), \(T_2 = 1\) and the three-term recurrence relation \(T_n = T_{n-1} + T_{n-2} + T_{n-3}\).

The Tribonacci polynomials are defined by \(T_0(x) = 0\), \(T_1(x) = 1\), \(T_2(x) = x^2\), and \(T_n(x) = x^2 T_{n-1}(x) + x T_{n-2}(x) + T_{n-3}(x)\) for \(n > 2\). For all positive integers \(n\), \(T_n(1) = T_n\).

tribonacci(n) gives the \(n^{th}\) Tribonacci number, \(T_n\)
tribonacci(n, x) gives the \(n^{th}\) Tribonacci polynomial in \(x\), \(T_n(x)\)

Examples

>>> from sympy import tribonacci, Symbol

>>> [tribonacci(x) for x in range(11)]
[0, 1, 1, 2, 4, 7, 13, 24, 44, 81, 149]
>>> tribonacci(5, Symbol('t'))
t**8 + 3*t**5 + 3*t**2

References

R161: https://en.wikipedia.org/wiki/Generalizations_of_Fibonacci_numbers#Tribonacci_numbers
R162: http://mathworld.wolfram.com/TribonacciNumber.html
R163: https://oeis.org/A000073

harmonic¶

class sympy.functions.combinatorial.numbers.harmonic[source]¶

Harmonic numbers

The nth harmonic number is given by \(\operatorname{H}_{n} = 1 + \frac{1}{2} + \frac{1}{3} + \ldots + \frac{1}{n}\).

More generally:

\[\operatorname{H}_{n,m} = \sum_{k=1}^{n} \frac{1}{k^m}\]

As \(n \rightarrow \infty\), \(\operatorname{H}_{n,m} \rightarrow \zeta(m)\), the Riemann zeta function.

harmonic(n) gives the nth harmonic number, \(\operatorname{H}_n\)
harmonic(n, m) gives the nth generalized harmonic number of order \(m\), \(\operatorname{H}_{n,m}\), where harmonic(n) == harmonic(n, 1)

Examples

>>> from sympy import harmonic, oo

>>> [harmonic(n) for n in range(6)]
[0, 1, 3/2, 11/6, 25/12, 137/60]
>>> [harmonic(n, 2) for n in range(6)]
[0, 1, 5/4, 49/36, 205/144, 5269/3600]
>>> harmonic(oo, 2)
pi**2/6

>>> from sympy import Symbol, Sum
>>> n = Symbol("n")

>>> harmonic(n).rewrite(Sum)
Sum(1/_k, (_k, 1, n))

We can evaluate harmonic numbers for all integral and positive rational arguments:

>>> from sympy import S, expand_func, simplify
>>> harmonic(8)
761/280
>>> harmonic(11)
83711/27720

>>> H = harmonic(1/S(3))
>>> H
harmonic(1/3)
>>> He = expand_func(H)
>>> He
-log(6) - sqrt(3)*pi/6 + 2*Sum(log(sin(_k*pi/3))*cos(2*_k*pi/3), (_k, 1, 1))
                       + 3*Sum(1/(3*_k + 1), (_k, 0, 0))
>>> He.doit()
-log(6) - sqrt(3)*pi/6 - log(sqrt(3)/2) + 3
>>> H = harmonic(25/S(7))
>>> He = simplify(expand_func(H).doit())
>>> He
log(sin(pi/7)**(-2*cos(pi/7))*sin(2*pi/7)**(2*cos(16*pi/7))*cos(pi/14)**(-2*sin(pi/14))/14)
+ pi*tan(pi/14)/2 + 30247/9900
>>> He.n(40)
1.983697455232980674869851942390639915940
>>> harmonic(25/S(7)).n(40)
1.983697455232980674869851942390639915940

We can rewrite harmonic numbers in terms of polygamma functions:

>>> from sympy import digamma, polygamma
>>> m = Symbol("m")

>>> harmonic(n).rewrite(digamma)
polygamma(0, n + 1) + EulerGamma

>>> harmonic(n).rewrite(polygamma)
polygamma(0, n + 1) + EulerGamma

>>> harmonic(n,3).rewrite(polygamma)
polygamma(2, n + 1)/2 - polygamma(2, 1)/2

>>> harmonic(n,m).rewrite(polygamma)
(-1)**m*(polygamma(m - 1, 1) - polygamma(m - 1, n + 1))/factorial(m - 1)

Integer offsets in the argument can be pulled out:

>>> from sympy import expand_func

>>> expand_func(harmonic(n+4))
harmonic(n) + 1/(n + 4) + 1/(n + 3) + 1/(n + 2) + 1/(n + 1)

>>> expand_func(harmonic(n-4))
harmonic(n) - 1/(n - 1) - 1/(n - 2) - 1/(n - 3) - 1/n

Some limits can be computed as well:

>>> from sympy import limit, oo

>>> limit(harmonic(n), n, oo)
oo

>>> limit(harmonic(n, 2), n, oo)
pi**2/6

>>> limit(harmonic(n, 3), n, oo)
-polygamma(2, 1)/2

However we can not compute the general relation yet:

>>> limit(harmonic(n, m), n, oo)
harmonic(oo, m)

which equals zeta(m) for m > 1.

References

R164: https://en.wikipedia.org/wiki/Harmonic_number
R165: http://functions.wolfram.com/GammaBetaErf/HarmonicNumber/
R166: http://functions.wolfram.com/GammaBetaErf/HarmonicNumber2/

lucas¶

class sympy.functions.combinatorial.numbers.lucas[source]¶

Lucas numbers

Lucas numbers satisfy a recurrence relation similar to that of the Fibonacci sequence, in which each term is the sum of the preceding two. They are generated by choosing the initial values \(L_0 = 2\) and \(L_1 = 1\).

lucas(n) gives the \(n^{th}\) Lucas number

Examples

>>> from sympy import lucas

>>> [lucas(x) for x in range(11)]
[2, 1, 3, 4, 7, 11, 18, 29, 47, 76, 123]

References

R167: https://en.wikipedia.org/wiki/Lucas_number
R168: http://mathworld.wolfram.com/LucasNumber.html

genocchi¶

class sympy.functions.combinatorial.numbers.genocchi[source]¶

Genocchi numbers

The Genocchi numbers are a sequence of integers \(G_n\) that satisfy the relation:

\[\frac{2t}{e^t + 1} = \sum_{n=1}^\infty \frac{G_n t^n}{n!}\]

Examples

>>> from sympy import Symbol
>>> from sympy.functions import genocchi
>>> [genocchi(n) for n in range(1, 9)]
[1, -1, 0, 1, 0, -3, 0, 17]
>>> n = Symbol('n', integer=True, positive=True)
>>> genocchi(2*n + 1)
0

References

R169: https://en.wikipedia.org/wiki/Genocchi_number
R170: http://mathworld.wolfram.com/GenocchiNumber.html

partition¶

class sympy.functions.combinatorial.numbers.partition[source]¶

Partition numbers

The Partition numbers are a sequence of integers \(p_n\) that represent the number of distinct ways of representing \(n\) as a sum of natural numbers (with order irrelevant). The generating function for \(p_n\) is given by:

\[\sum_{n=0}^\infty p_n x^n = \prod_{k=1}^\infty (1 - x^k)^{-1}\]

Examples

>>> from sympy import Symbol
>>> from sympy.functions import partition
>>> [partition(n) for n in range(9)]
[1, 1, 2, 3, 5, 7, 11, 15, 22]
>>> n = Symbol('n', integer=True, negative=True)
>>> partition(n)
0

References

R171: https://en.wikipedia.org/wiki/Partition_(number_theory%29
R172: https://en.wikipedia.org/wiki/Pentagonal_number_theorem

MultiFactorial¶

class sympy.functions.combinatorial.factorials.MultiFactorial[source]¶

RisingFactorial¶

class sympy.functions.combinatorial.factorials.RisingFactorial[source]¶

Rising factorial (also called Pochhammer symbol) is a double valued function arising in concrete mathematics, hypergeometric functions and series expansions. It is defined by:

\[rf(x,k) = x \cdot (x+1) \cdots (x+k-1)\]

where \(x\) can be arbitrary expression and \(k\) is an integer. For more information check “Concrete mathematics” by Graham, pp. 66 or visit http://mathworld.wolfram.com/RisingFactorial.html page.

When \(x\) is a Poly instance of degree >= 1 with a single variable, \(rf(x,k) = x(y) \cdot x(y+1) \cdots x(y+k-1)\), where \(y\) is the variable of \(x\). This is as described in Peter Paule, “Greatest Factorial Factorization and Symbolic Summation”, Journal of Symbolic Computation, vol. 20, pp. 235-268, 1995.

Examples

>>> from sympy import rf, symbols, factorial, ff, binomial, Poly
>>> from sympy.abc import x
>>> n, k = symbols('n k', integer=True)
>>> rf(x, 0)
1
>>> rf(1, 5)
120
>>> rf(x, 5) == x*(1 + x)*(2 + x)*(3 + x)*(4 + x)
True
>>> rf(Poly(x**3, x), 2)
Poly(x**6 + 3*x**5 + 3*x**4 + x**3, x, domain='ZZ')

Rewrite

>>> rf(x, k).rewrite(ff)
FallingFactorial(k + x - 1, k)
>>> rf(x, k).rewrite(binomial)
binomial(k + x - 1, k)*factorial(k)
>>> rf(n, k).rewrite(factorial)
factorial(k + n - 1)/factorial(n - 1)

References

R173: https://en.wikipedia.org/wiki/Pochhammer_symbol

stirling¶

sympy.functions.combinatorial.numbers.stirling(n, k, d=None, kind=2, signed=False)[source]¶

Return Stirling number \(S(n, k)\) of the first or second (default) kind.

The sum of all Stirling numbers of the second kind for \(k = 1\) through \(n\) is bell(n). The recurrence relationship for these numbers is:

\[{0 \brace 0} = 1; {n \brace 0} = {0 \brace k} = 0;\]

\[{{n+1} \brace k} = j {n \brace k} + {n \brace {k-1}}\]

where \(j\) is:: \(n\) for Stirling numbers of the first kind \(-n\) for signed Stirling numbers of the first kind \(k\) for Stirling numbers of the second kind

The first kind of Stirling number counts the number of permutations of n distinct items that have k cycles; the second kind counts the ways in which n distinct items can be partitioned into k parts. If d is given, the “reduced Stirling number of the second kind” is returned: S^{d}(n, k) = S(n - d + 1, k - d + 1) with n >= k >= d. (This counts the ways to partition n consecutive integers into k groups with no pairwise difference less than d. See example below.)

To obtain the signed Stirling numbers of the first kind, use keyword signed=True. Using this keyword automatically sets kind to 1.

Examples

>>> from sympy.functions.combinatorial.numbers import stirling, bell
>>> from sympy.combinatorics import Permutation
>>> from sympy.utilities.iterables import multiset_partitions, permutations

First kind (unsigned by default):

>>> [stirling(6, i, kind=1) for i in range(7)]
[0, 120, 274, 225, 85, 15, 1]
>>> perms = list(permutations(range(4)))
>>> [sum(Permutation(p).cycles == i for p in perms) for i in range(5)]
[0, 6, 11, 6, 1]
>>> [stirling(4, i, kind=1) for i in range(5)]
[0, 6, 11, 6, 1]

First kind (signed):

>>> [stirling(4, i, signed=True) for i in range(5)]
[0, -6, 11, -6, 1]

Second kind:

>>> [stirling(10, i) for i in range(12)]
[0, 1, 511, 9330, 34105, 42525, 22827, 5880, 750, 45, 1, 0]
>>> sum(_) == bell(10)
True
>>> len(list(multiset_partitions(range(4), 2))) == stirling(4, 2)
True

Reduced second kind:

>>> from sympy import subsets, oo
>>> def delta(p):
...    if len(p) == 1:
...        return oo
...    return min(abs(i[0] - i[1]) for i in subsets(p, 2))
>>> parts = multiset_partitions(range(5), 3)
>>> d = 2
>>> sum(1 for p in parts if all(delta(i) >= d for i in p))
7
>>> stirling(5, 3, 2)
7

References

R174: https://en.wikipedia.org/wiki/Stirling_numbers_of_the_first_kind
R175: https://en.wikipedia.org/wiki/Stirling_numbers_of_the_second_kind

Enumeration¶

Three functions are available. Each of them attempts to efficiently compute a given combinatorial quantity for a given set or multiset which can be entered as an integer, sequence or multiset (dictionary with elements as keys and multiplicities as values). The k parameter indicates the number of elements to pick (or the number of partitions to make). When k is None, the sum of the enumeration for all k (from 0 through the number of items represented by n) is returned. A replacement parameter is recognized for combinations and permutations; this indicates that any item may appear with multiplicity as high as the number of items in the original set.

>>> from sympy.functions.combinatorial.numbers import nC, nP, nT
>>> items = 'baby'

nC¶

Calculate the number of combinations of length k.

>>> [nC(items, k) for k in range(len(items) + 1)], nC(items)
([1, 3, 4, 3, 1], 12)
>>> nC('aaa', 2)
1
>>> nC('abc', 2)
3
>>> nC(3, 2)
3

nP¶

Calculate the number of permutations of length k.

>>> [nP(items, k) for k in range(len(items) + 1)], nP(items)
([1, 3, 7, 12, 12], 35)
>>> nC('aaa', 2)
1
>>> nC('abc', 2)
3
>>> nC(3, 2)
3

nT¶

Calculate the number of partitions that have k parts.

>>> [nT(items, k) for k in range(len(items) + 1)], nT(items)
([0, 1, 5, 4, 1], 11)
>>> nT('aaa', 2)
1
>>> nT('abc', 2)
3
>>> nT(3, 2)
1

Note that the integer for n indicates identical items for nT but indicates n different items for nC and nP.

Combinatorial¶

bell¶

bernoulli¶

binomial¶

catalan¶

euler¶

factorial¶

subfactorial¶

factorial2 / double factorial¶

FallingFactorial¶

fibonacci¶

tribonacci¶

harmonic¶

lucas¶

genocchi¶

partition¶

MultiFactorial¶

RisingFactorial¶

stirling¶

Enumeration¶

nC¶

nP¶

nT¶

Table of Contents

Previous topic

Next topic

This Page