a RG5d@sdZddlmZddlmZddlmZmZm Z ddl m Z m Z mZmZddlmZddlmZddlmZdd lmZmZdd lmZdd lmZdd lmZmZm Z m!Z!m"Z"dd l#m$Z$m%Z%ddl&m'Z'ddl(m)Z)m*Z*m+Z+ddl,m-Z-m.Z.ddl/m0Z0ddl1m2Z2ddl3m4Z4m5Z5m6Z6ddl7m8Z8ddl9m:Z:ddl;mm?Z@ddZAe dZBddZCGdddeZDGdd d eZEGd!d"d"eZFGd#d$d$eZGGd%d&d&eZHGd'd(d(eZIGd)d*d*eZJGd+d,d,eZKGd-d.d.eZLGd/d0d0eZMd1d1gZNGd2d3d3eZOGd4d5d5ePZQd6ZRd7ZSeTd8eSZUd9d:ZVdXdd?ZXed@dAZYdZdBdCZZdDdEZ[edFdGZ\dHdIZ]edJdKZ^d[dMdNZ_edOdPZ`d\dQdRZaGdSdTdTeZbd]d8d8dUdVdWZcd8S)^a6 This module implements some special functions that commonly appear in combinatorial contexts (e.g. in power series); in particular, sequences of rational numbers such as Bernoulli and Fibonacci numbers. Factorials, binomial coefficients and related functions are located in the separate 'factorials' module. )prod) defaultdict)CallableDictTuple)SSymbolAddDummy)cacheit) pure_complex)Expr)Function expand_mul) fuzzy_not)Mul)EpiooRationalInteger)is_leis_gt) SYMPY_INTS)binomial factorial subfactorial)isprime is_square)MultisetPartitionTraversersympy_deprecation_warning)multisetmultiset_derangementsiterable)recurrence_memo)as_int)bernfracworkprec)ifibcCstt||dSN)rrange)abr/a/var/www/html/django/DPS/env/lib/python3.9/site-packages/sympy/functions/combinatorial/numbers.py_product$sr1xcCs ||dkS)Nrr/pnr/r/r0_divides2sr6c@sXeZdZdZeddZeddZeddZedd Zed d Z ed d Z dS) carmichaela0 Carmichael Numbers: Certain cryptographic algorithms make use of big prime numbers. However, checking whether a big number is prime is not so easy. Randomized prime number checking tests exist that offer a high degree of confidence of accurate determination at low cost, such as the Fermat test. Let 'a' be a random number between $2$ and $n - 1$, where $n$ is the number whose primality we are testing. Then, $n$ is probably prime if it satisfies the modular arithmetic congruence relation: .. math :: a^{n-1} = 1 \pmod{n} (where mod refers to the modulo operation) If a number passes the Fermat test several times, then it is prime with a high probability. Unfortunately, certain composite numbers (non-primes) still pass the Fermat test with every number smaller than themselves. These numbers are called Carmichael numbers. A Carmichael number will pass a Fermat primality test to every base $b$ relatively prime to the number, even though it is not actually prime. This makes tests based on Fermat's Little Theorem less effective than strong probable prime tests such as the Baillie-PSW primality test and the Miller-Rabin primality test. Examples ======== >>> from sympy import carmichael >>> carmichael.find_first_n_carmichaels(5) [561, 1105, 1729, 2465, 2821] >>> carmichael.find_carmichael_numbers_in_range(0, 562) [561] >>> carmichael.find_carmichael_numbers_in_range(0,1000) [561] >>> carmichael.find_carmichael_numbers_in_range(0,2000) [561, 1105, 1729] References ========== .. [1] https://en.wikipedia.org/wiki/Carmichael_number .. [2] https://en.wikipedia.org/wiki/Fermat_primality_test .. [3] https://www.jstor.org/stable/23248683?seq=1#metadata_info_tab_contents cCstddddt|S)Nzt is_perfect_square is just a wrapper around sympy.ntheory.primetest.is_square so use that directly instead. 1.11$deprecated-carmichael-static-methodsdeprecated_since_versionactive_deprecations_target)r!rr5r/r/r0is_perfect_squarehs zcarmichael.is_perfect_squarecCstdddd||dkS)NzI divides can be replaced by directly testing n % p == 0. r8r9r:rr r3r/r/r0dividests zcarmichael.dividescCstddddt|S)Nzi is_prime is just a wrapper around sympy.ntheory.primetest.isprime so use that directly instead. r8r9r:)r!rr=r/r/r0is_primes zcarmichael.is_primecsdkrdks$ts$ddkr(dSdg}|fddtddddD|D]:}t|rv|dkrvdSt|r\t|dds\dSq\dStd dS) Nrr+Fcsg|]}|dkr|qSrr/.0ir=r/r0 z,carmichael.is_carmichael..Tz6The provided number must be greater than or equal to 0)rextendr,rr6 ValueError)r5divisorsrEr/r=r0 is_carmichaels(zcarmichael.is_carmichaelcCsbd|kr|krVnn>|ddkr>ddt|d|dDSddt||dDSntddS)NrrAcSsg|]}t|r|qSr/r7rLrCr/r/r0rFrGz?carmichael.find_carmichael_numbers_in_range..r+cSsg|]}t|r|qSr/rMrCr/r/r0rFrGzQThe provided range is not valid. x and y must be non-negative integers and x <= y)r,rJ)r2yr/r/r0 find_carmichael_numbers_in_ranges  z+carmichael.find_carmichael_numbers_in_rangecCs8d}t}t||kr4t|r*|||d7}q |SNr+rA)listlenr7rLappend)r5rEZ carmichaelsr/r/r0find_first_n_carmichaelss    z#carmichael.find_first_n_carmichaelsN) __name__ __module__ __qualname____doc__ staticmethodr>r?r@rLrOrTr/r/r/r0r75s2     r7c@sVeZdZdZeddZeedeje gddZ e d ddZ d d Z d d ZdS) fibonaccia 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 .. math :: 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 See Also ======== bell, bernoulli, catalan, euler, harmonic, lucas, genocchi, partition, tribonacci References ========== .. [1] https://en.wikipedia.org/wiki/Fibonacci_number .. [2] http://mathworld.wolfram.com/FibonacciNumber.html cCst|SN)_ifibr=r/r/r0_fibszfibonacci._fibNcCs|dt|dS)N_symexpandr5prevr/r/r0_fibpolyszfibonacci._fibpolycCs||tjurtjS|jrx|durVt|}|dkrFtj|dt| St||Sn"|dkrftd| | t |SdS)Nrr+zDFibonacci polynomials are defined only for positive integer indices.) rInfinity is_Integerint NegativeOnerZrr]rJresubsraclsr5symr/r/r0evals zfibonacci.evalcKsDddlm}d| |dd|d||d d|dS)NrsqrtrAr+(sympy.functions.elementary.miscellaneousrpselfr5kwargsrpr/r/r0_eval_rewrite_as_sqrts zfibonacci._eval_rewrite_as_sqrtcKs(tj|dtj |dtjdSrP)r GoldenRatiorur5rvr/r/r0_eval_rewrite_as_GoldenRatio sz&fibonacci._eval_rewrite_as_GoldenRatio)N)rUrVrWrXrYr]r%rOnerare classmethodrnrwrzr/r/r/r0rZs)   rZc@s$eZdZdZeddZddZdS)lucasa 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] See Also ======== bell, bernoulli, catalan, euler, fibonacci, harmonic, genocchi, partition, tribonacci References ========== .. [1] https://en.wikipedia.org/wiki/Lucas_number .. [2] http://mathworld.wolfram.com/LucasNumber.html cCs2|tjurtjS|jr.t|dt|dSdSr*)rrfrgrZrlr5r/r/r0rn5s z lucas.evalcKs8ddlm}d| d|d||d d|S)NrrorAr+rqrrrtr/r/r0rw=s zlucas._eval_rewrite_as_sqrtN)rUrVrWrXr|rnrwr/r/r/r0r}s r}c@speZdZdZeeejejejgddZ eeejeje dgddZ e ddd Z d d Zd d ZdS) tribonaccia 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 See Also ======== bell, bernoulli, catalan, euler, fibonacci, harmonic, lucas, genocchi, partition References ========== .. [1] https://en.wikipedia.org/wiki/Generalizations_of_Fibonacci_numbers#Tribonacci_numbers .. [2] http://mathworld.wolfram.com/TribonacciNumber.html .. [3] https://oeis.org/A000073 cCs|d|d|dS)Nr^r_r/rcr/r/r0_tribpsztribonacci._tribrAcCs(|dt|dtd|dS)Nrr^rAr_r`rcr/r/r0 _tribpolyusztribonacci._tribpolyNcCsZ|tjurtjS|jrVt|}|dkr.td|durDt||S||t |SdS)NrzITribonacci polynomials are defined only for non-negative integer indices.) rrfrgrhrJrrrrjrarkr/r/r0rnzs ztribonacci.evalc Ks&ddlm}m}dtj|dd}d|dd|d|dd|dd}d||dd|d|d|dd|dd}d|d|dd|d||dd|dd}||d||||||d||||||d||||} | S) Nrcbrtrpr_rHrAr+!)rsrrpr ImaginaryUnit) rur5rvrrpwr-r.cTnr/r/r0rws0<<z tribonacci._eval_rewrite_as_sqrtcKsdddlm}ddlm}m}|dd|d}d|tj||dd|d }||tjS) NrfloorriJfrrHrA)#sympy.functions.elementary.integersrrsrrprTribonacciConstantHalf)rur5rvrrrpr.rr/r/r0#_eval_rewrite_as_TribonacciConstants  &z.tribonacci._eval_rewrite_as_TribonacciConstant)N)rUrVrWrXrYr%rZeror{rrarr|rnrwrr/r/r/r0rIs&    rc@s^eZdZUdZeeed<eddZe j e dde ddd Z d d d d Z edddZd S) bernoulliaI Bernoulli numbers / Bernoulli polynomials The Bernoulli numbers are a sequence of rational numbers defined by `B_0 = 1` and the recursive relation (`n > 0`): .. math :: 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: .. math :: 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: .. math :: B_n = \frac{A(n) - S(n)}{\binom{n+3}{n}} where: .. math :: 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} and: .. math :: 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 See Also ======== bell, catalan, euler, fibonacci, harmonic, lucas, genocchi, partition, tribonacci References ========== .. [1] https://en.wikipedia.org/wiki/Bernoulli_number .. [2] https://en.wikipedia.org/wiki/Bernoulli_polynomial .. [3] http://mathworld.wolfram.com/BernoulliNumber.html .. [4] http://mathworld.wolfram.com/BernoulliPolynomial.html argscCsd}tt|d|d}td|ddD]`}||t|d|7}|t|dd|d|d|9}|td|dd|d}q,|ddkrt|dd |}nt|dd|}|t|d|S)NrrHr+r )rhrr,rr1r)r5sr-jr/r/r0_calc_bernoullis&  zbernoulli._calc_bernoullir+rr_)rrArrrArNc sNjr&jrjrjr$tjStjurLdur@tddStjSnЈdurjr`tj St dkrt \}}tt |t |Sd}j |}|krj St|dddD]"}|}|j |<|j |<q|St fddtdD} t| SntddurJjrJdjrJtj SdS) Nr_rAircs*g|]"}t|||qSr/rrDkrkr/r0rFrGz"bernoulli.eval..r+zCBernoulli numbers are defined only for nonnegative integer indices.) is_Numberrgis_nonnegativeis_zerorr{rris_oddrrhr'_highest_cacher,rr rJ is_positive) rlr5rmr4qcaseZhighest_cachedrEr.resultr/rkr0rns>           zbernoulli.eval)N)rUrVrWrXtTupler__annotations__rYrrr{rrrr|rnr/r/r/r0rs A   rc@sfeZdZdZeeddgddZeeeje gddZ eddZ e dd d Z dd d Zd S)bellax Bell numbers / Bell polynomials The Bell numbers satisfy `B_0 = 1` and .. math:: B_n = \sum_{k=0}^{n-1} \binom{n-1}{k} B_k. They are also given by: .. math:: B_n = \frac{1}{e} \sum_{k=0}^{\infty} \frac{k^n}{k!}. The Bell polynomials are given by `B_0(x) = 1` and .. math:: 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 .. math:: 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 See Also ======== bernoulli, catalan, euler, fibonacci, harmonic, lucas, genocchi, partition, tribonacci References ========== .. [1] https://en.wikipedia.org/wiki/Bell_number .. [2] http://mathworld.wolfram.com/BellNumber.html .. [3] http://mathworld.wolfram.com/BellPolynomial.html r+cCs<d}d}td|D]$}||||}||||7}q|Sr*)r,r5rdrr-rr/r/r0_bellrs z bell._bellcCsTd}d}td|dD]0}|||d|d}||||d7}qtt|SrP)r,rrarr/r/r0 _bell_poly|s zbell._bell_polycCs|dkr|dkrtjS|dks&|dkr,tjStj}tj}td||dD]>}||t|||d|||d7}||||}qJt|S)a The second kind of Bell polynomials (incomplete Bell polynomials). Calculated by recurrence formula: .. math:: B_{n,k}(x_1, x_2, \dotsc, x_{n-k+1}) = \sum_{m=1}^{n-k+1} \x_m \binom{n-1}{m-1} B_{n-m,k-1}(x_1, x_2, \dotsc, x_{n-m-k}) where `B_{0,0} = 1;` `B_{n,0} = 0; for n \ge 1` `B_{0,k} = 0; for k \ge 1` rr+rA)rr{rr,r_bell_incomplete_polyr)r5rsymbolsrr-mr/r/r0rs zbell._bell_incomplete_polyNcCs|tjur |durtjStd|js0|jdur8td|jr|jr|dur^t|t |S|dur|| t | t |S| t |t ||}|SdS)NzBell polynomial is not definedFza non-negative integer expected)rrfrJ is_negative is_integerrgrrrrhrrjrar)rlr5k_symrrr/r/r0rns  z bell.evalcKs^ddlm}|dus|dur |S|js*|Stdddd}dt|||t||dtjfS)NrSumrT)integer nonnegativer+)sympy.concrete.summationsrrr rrrrf)rur5rrrvrrr/r/r0_eval_rewrite_as_Sums zbell._eval_rewrite_as_Sum)NN)NN)rUrVrWrXrYr%rrr{rarrr|rnrr/r/r/r0r1s@      rc@sdeZdZdZiZedddZdddZddd Zdd d Z dd d Z ddZ dddZ ddZ dS)harmonica+ 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: .. math:: \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(2*pi/7)**(2*cos(16*pi/7))/(14*sin(pi/7)**(2*cos(pi/7))*cos(pi/14)**(2*sin(pi/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 cannot compute the general relation yet: >>> limit(harmonic(n, m), n, oo) harmonic(oo, m) which equals ``zeta(m)`` for ``m > 1``. See Also ======== bell, bernoulli, catalan, euler, fibonacci, lucas, genocchi, partition, tribonacci References ========== .. [1] https://en.wikipedia.org/wiki/Harmonic_number .. [2] http://functions.wolfram.com/GammaBetaErf/HarmonicNumber/ .. [3] http://functions.wolfram.com/GammaBetaErf/HarmonicNumber2/ Ncsddlm}tjur||Sdur,tjjr6|S|tjurvjrLtjSttjr^tjSt tjrr|SdS|dkrtj S|j r|j r҈j r҈|j vrtdgfdd}||j <|j t|SdS)Nr)zetacs|dtj|S)Nr_rr{rcrr/r0ffszharmonic.eval..f)&sympy.functions.special.zeta_functionsrrr{rrfrNaNrrrrgr _functionsr%rh)rlr5rrrr/rr0rnLs.       z harmonic.evalr+cKsBddlm}tj|t|d||dd||d|dS)Nr polygammar+)'sympy.functions.special.gamma_functionsrrrirrur5rrvrr/r/r0_eval_rewrite_as_polygammals z#harmonic._eval_rewrite_as_polygammacKsddlm}||SNrrrrrewriterr/r/r0_eval_rewrite_as_digammaps z!harmonic._eval_rewrite_as_digammacKsddlm}||Srrrr/r/r0_eval_rewrite_as_trigammats z"harmonic._eval_rewrite_as_trigammacKs<ddlm}tddd}|dur&tj}||| |d|fS)NrrrTrr+)rrr rr{)rur5rrvrrr/r/r0rxs   zharmonic._eval_rewrite_as_Sumc sddlm}|jd}t|jdkr.|jdnd}|tjkr|jr|jd}|||jr|jrfddt |ddDt g}t |S|jr|j rƇfddt d|dDt g}t |S|j r|\}}||} || |}| jr|jr|jr||krdd lm} dd lm} dd lm} m} m}td }||d||||d| f}d|| dt||t|| | t|t||d| |dtdf}td|t||| d|}|||S|S) NrrrAr+csg|]}tj|qSr/rrCnnewr/r0rFrGz.harmonic._eval_expand_func..r_csg|]}tj |qSr/rrCrr/r0rFrGlogr)sincoscotr)rrrrRrr{is_Addrgrr,rr r is_Rationalas_numer_denomr&sympy.functions.elementary.exponentialrrr(sympy.functions.elementary.trigonometricrrrr r)ruhintsrr5roffrr4rurrrrrrt1t2t3r/rr0_eval_expand_funcs>     $ $  "   $ zharmonic._eval_expand_funccKs ddlm}||jdddS)Nrr tractableT)deepr)rur5rlimitvarrvrr/r/r0_eval_rewrite_as_tractables z#harmonic._eval_rewrite_as_tractablecCs4ddlm}tdd|jDr0|||SdS)Nrrcss|] }|jVqdSr[) is_numberrCr/r/r0 rGz'harmonic._eval_evalf..)rrallrr _eval_evalf)ruprecrr/r/r0rs zharmonic._eval_evalf)N)r+)r+)r+)N)r+N)rUrVrWrXrr|rnrrrrrrrr/r/r/r0rs     $ rc@s0eZdZdZed ddZd ddZddZdS) eulera Euler numbers / Euler polynomials The Euler numbers are given by: .. math:: 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} .. math:: E_{2n+1} = 0 Euler numbers and Euler polynomials are related by .. math:: E_n = 2^n E_n\left(\frac{1}{2}\right). We compute symbolic Euler polynomials using [5]_ .. math:: 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 euler, Symbol, S >>> [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 See Also ======== bell, bernoulli, catalan, fibonacci, harmonic, lucas, genocchi, partition, tribonacci References ========== .. [1] https://en.wikipedia.org/wiki/Euler_numbers .. [2] http://mathworld.wolfram.com/EulerNumber.html .. [3] https://en.wikipedia.org/wiki/Alternating_permutation .. [4] http://mathworld.wolfram.com/AlternatingPermutation.html .. [5] http://dlmf.nist.gov/24.2#ii NcsTjr0jr(jr(durZjr,tjSddlm}|j |j dd}t |St dd}|rt dd|Drtdd|Drddlm}ttd d |D}t||}Wdn1s0Yt||Stfd d td D}t|Sntd durPjrPjrPtjSdS)NrmpT)exact)or_realcss|]}|jp|jVqdSr[)is_FloatrgrDr-r/r/r0r rGzeuler.eval..css|] }|jVqdSr[)rrr/r/r0r rGcSsg|]}|jr|jqSr/)r_precrr/r/r0rFrGzeuler.eval..cs8g|]0}t||d|tj|qS)rA)rrrrrlrrmr/r0rFrGr+z?Euler numbers are defined only for nonnegative integer indices.)rrgrrrrmpmathr _to_mpmathreulernumrr ranyrhminr( eulerpolyr _from_mpmathr,r rbrJr)rlrrmrresZreimrrr/rr0rns4     *  z euler.evalcKsddlm}|dur|jrtddd}tddd}|d}tj||t||tj||d|d|dd|tj|||d|f|dd|df}|S|rtddd}|t||t|d||tj |||d|fSdS) NrrrTrrrAr+) rris_evenr rrrrirr)rur5r2rvrrrZEmr/r/r0r s$    zeuler._eval_rewrite_as_SumcCst|jdkr|jddfn|j\}}|dur|jr|jrddlm}||}t|||}Wdn1sx0Yt ||S|r |j r |jr |jr ddlm}t |}||}t|| ||}Wdn1s0Yt ||SdS)Nr+rr)rRrrgrrrrr(rr rrrhr)rurrr2rrr/r/r0r.s&   (    *zeuler._eval_evalf)N)N)rUrVrWrXr|rnrrr/r/r/r0rs G # rc@speZdZdZeddZdddZddZd d Zdd d Z ddZ ddZ ddZ ddZ ddZddZdS)catalana; Catalan numbers The `n^{th}` catalan number is given by: .. math :: 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, ... 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 See Also ======== bell, bernoulli, euler, fibonacci, harmonic, lucas, genocchi, partition, tribonacci sympy.functions.combinatorial.factorials.binomial References ========== .. [1] https://en.wikipedia.org/wiki/Catalan_number .. [2] http://mathworld.wolfram.com/CatalanNumber.html .. [3] http://functions.wolfram.com/GammaBetaErf/CatalanNumber/ .. [4] http://geometer.org/mathcircles/catalan.pdf cCsddlm}|jr|js$|jrP|jrPd|||tj|tj||dS|jr|jr|djrltj S|dj rt ddSdS)NrgammarrAr+r_) rrrgr is_nonintegerrrrrrrr)rlr5rr/r/r0rns  ,   z catalan.evalr+cCsPddlm}ddlm}|jd}t||d|tj|d|d|dS)NrrrrAr)rrrrrrrr)ruargindexrrr5r/r/r0fdiffs   z catalan.fdiffcKstd|||dSNrAr+rryr/r/r0_eval_rewrite_as_binomialsz!catalan._eval_rewrite_as_binomialcKs td|t|dt|Srrryr/r/r0_eval_rewrite_as_factorialsz"catalan._eval_rewrite_as_factorialTcKs8ddlm}d|||tj|tj||dS)NrrrrA)rrrr)rur5 piecewiservrr/r/r0_eval_rewrite_as_gammas zcatalan._eval_rewrite_as_gammacKs$ddlm}|d|| gdgdS)Nr)hyperr+rA)sympy.functions.special.hyperr )rur5rvr r/r/r0_eval_rewrite_as_hypers zcatalan._eval_rewrite_as_hypercKsBddlm}|jr|js|Stdddd}|||||d|fS)Nr)ProductrT)rpositiverA)sympy.concrete.productsr rrr )rur5rvr rr/r/r0_eval_rewrite_as_Products   z catalan._eval_rewrite_as_ProductcCs |jdjr|jdjrdSdSNrT)rrrrur/r/r0_eval_is_integerszcatalan._eval_is_integercCs|jdjrdSdSr)rrrr/r/r0_eval_is_positives zcatalan._eval_is_positivecCs$|jdjr |jddjr dSdS)NrrHTrrrrr/r/r0_eval_is_compositeszcatalan._eval_is_compositecCs,ddlm}|jdjr(|||SdS)Nrr)rrrrrr)rurrr/r/r0rs  zcatalan._eval_evalfN)r+)T)rUrVrWrXr|rnrrrrr rrrrrr/r/r/r0rFsO  rc@sTeZdZdZeddZddZddZdd Zd d Z d d Z ddZ ddZ dS)genocchia Genocchi numbers The Genocchi numbers are a sequence of integers `G_n` that satisfy the relation: .. math:: \frac{2t}{e^t + 1} = \sum_{n=1}^\infty \frac{G_n t^n}{n!} Examples ======== >>> from sympy import genocchi, Symbol >>> [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 See Also ======== bell, bernoulli, catalan, euler, fibonacci, harmonic, lucas, partition, tribonacci References ========== .. [1] https://en.wikipedia.org/wiki/Genocchi_number .. [2] http://mathworld.wolfram.com/GenocchiNumber.html cCs`|jr6|jr|jrtdddtd|t|S|jrL|djrLtjS|dj r\tj SdS)Nz7Genocchi numbers are defined only for positive integersrAr+) rrgis_nonpositiverJrrrrrrr{r~r/r/r0rns  z genocchi.evalcKs,|jr(|jr(dtd|t|dSdSrP)rrrrryr/r/r0_eval_rewrite_as_bernoullis z#genocchi._eval_rewrite_as_bernoullicCs |jdjr|jdjrdSdSrrrr/r/r0r szgenocchi._eval_is_integercCs.|jd}|jr*|jr*|jr dS|djSdS)NrFrA)rrrrrur5r/r/r0_eval_is_negatives   zgenocchi._eval_is_negativecCs8|jd}|jr4|jr4|jr*t|djS|djSdSNrr+rA)rrrrrrrr/r/r0rs   zgenocchi._eval_is_positivecCs.|jd}|jr*|jr*|jr dS|djSdS)NrFr+)rrrrrr/r/r0 _eval_is_evens   zgenocchi._eval_is_evencCs2|jd}|jr.|jr.|jr dSt|djSdS)NrTr+)rrrrrrr/r/r0 _eval_is_odd$s   zgenocchi._eval_is_oddcCs|jd}|djS)Nr)rrrr/r/r0_eval_is_prime+s zgenocchi._eval_is_primeN) rUrVrWrXr|rnrrrrrrrr/r/r/r0rs rr+c@s@eZdZdZeddZeddZddZdd Z d d Z d S) partitionaN 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: .. math:: \sum_{n=0}^\infty p_n x^n = \prod_{k=1}^\infty (1 - x^k)^{-1} Examples ======== >>> from sympy import partition, Symbol >>> [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 See Also ======== bell, bernoulli, catalan, euler, fibonacci, harmonic, lucas, genocchi, tribonacci References ========== .. [1] https://en.wikipedia.org/wiki/Partition_(number_theory%29 .. [2] https://en.wikipedia.org/wiki/Pentagonal_number_theorem cCstt}||krt|St||dD]}d\}}}d}|d|d7}||kr`|t||7}|d7}||}||kr|t||7}|dkrqq4||ddkr|n| 7}q4t|q&|S)Nr+)rrrrrHrA)rR _npartitionr,rS)r5L_nvr4rErgpr/r/r0 _partition[s$  zpartition._partitioncCsV|j}|dkrtdn:|rR|jr(tjS|js8|djr>tjS|jrRt| |SdS)NFz/Partition numbers are defined only for integersr+) rrJrrrrr{rgrr&)rlr5is_intr/r/r0rnss zpartition.evalcCs|jdjrdSdSrrrrr/r/r0rs zpartition._eval_is_integercCs|jdjrdSdS)NrFr(rr/r/r0rs zpartition._eval_is_negativecCs|jd}|jr|jrdSdSr)rrrrr/r/r0rs  zpartition._eval_is_positiveN) rUrVrWrXrYr&r|rnrrrr/r/r/r0r :s   r c@s eZdZdS)_MultisetHistogramN)rUrVrWr/r/r/r0r)sr)r_r^Nc sttrdtddDs$tt}tfddD}tfddD||gStt}t |}t }||krdg|||gtStt |t |}tt t |d|}D]}|||d7<qt |SdS) a6Return tuple used in permutation and combination counting. Input is a dictionary giving items with counts as values or a sequence of items (which need not be sorted). The data is stored in a class deriving from tuple so it is easily recognized and so it can be converted easily to a list. css |]}t|to|dkVqdS)rN) isinstancerh)rDr$r/r/r0rrGz&_multiset_histogram..c3s|]}|dkrdVqdS)rr+Nr/rr=r/r0rrGcs g|]}|dkr|qSrBr/rr=r/r0rFrGz'_multiset_histogram..r+rBN) r*dictrvaluesrJsumr)rQsetrRzipr,_multiset_histogram) r5totitemsrlenslennrdrEr/r=r0r0s$  r0FcCsDz t|}Wn&ty2ttt|||YS0tt|||S)a`Return the number of permutations of ``n`` items taken ``k`` at a time. Possible values for ``n``: integer - set of length ``n`` sequence - converted to a multiset internally multiset - {element: multiplicity} If ``k`` is None then the total of all permutations of length 0 through the number of items represented by ``n`` will be returned. If ``replacement`` is True then a given item can appear more than once in the ``k`` items. (For example, for 'ab' permutations of 2 would include 'aa', 'ab', 'ba' and 'bb'.) The multiplicity of elements in ``n`` is ignored when ``replacement`` is True but the total number of elements is considered since no element can appear more times than the number of elements in ``n``. Examples ======== >>> from sympy.functions.combinatorial.numbers import nP >>> from sympy.utilities.iterables import multiset_permutations, multiset >>> nP(3, 2) 6 >>> nP('abc', 2) == nP(multiset('abc'), 2) == 6 True >>> nP('aab', 2) 3 >>> nP([1, 2, 2], 2) 3 >>> [nP(3, i) for i in range(4)] [1, 3, 6, 6] >>> nP(3) == sum(_) True When ``replacement`` is True, each item can have multiplicity equal to the length represented by ``n``: >>> nP('aabc', replacement=True) 121 >>> [len(list(multiset_permutations('aaaabbbbcccc', i))) for i in range(5)] [1, 3, 9, 27, 81] >>> sum(_) 121 See Also ======== sympy.utilities.iterables.multiset_permutations References ========== .. [1] https://en.wikipedia.org/wiki/Permutation )r&rJr_nPr0)r5r replacementr/r/r0nPs ;  r8cs|dkr dSttr|dur>tfddtdDSrJ|S|krVdS|krft|S|dkrrSt|dSnzttr|durtfddttdDSrЈt|S|tkrt|t ddt DS|tkrdS|dkr tSd}t tt t D]}|sNq<td8<|dkrd|<td8<|t t|d7}td7<d|<n6|d8<|t t|d7}|d7<td7<q<|SdS)Nrr+c3s|]}t|VqdSr[r6rCr5r7r/r0r rGz_nP..c3s|]}t|VqdSr[r9rCr:r/r0rrGcSsg|]}|dkrt|qS)r+rrCr/r/r0rFrGz_nP..)r*rr-r,rr1r)_N_ITEMSr_MrQrRr6)r5rr7r1rEr/r:r0r6sT   $    r6c Cs*t|}t|}|dd}dg|d}|d|t||d|}|r|}|d}|dd}tdt|t|D]}||||d7<qt||D](}||||d|||7<qqPtt|}|dr||}n ||dd<tt } t |D]\}} | | |<q| S)afor n = (m1, m2, .., mk) return the coefficients of the polynomial, prod(sum(x**i for i in range(nj + 1)) for nj in n); i.e. the coefficients of the product of AOPs (all-one polynomials) or order given in n. The resulting coefficient corresponding to x**r is the number of r-length combinations of sum(n) elements with multiplicities given in n. The coefficients are given as a default dictionary (so if a query is made for a key that is not present, 0 will be returned). Examples ======== >>> from sympy.functions.combinatorial.numbers import _AOP_product >>> from sympy.abc import x >>> n = (2, 2, 3) # e.g. aabbccc >>> prod = ((x**2 + x + 1)*(x**2 + x + 1)*(x**3 + x**2 + x + 1)).expand() >>> c = _AOP_product(n); dict(c) {0: 1, 1: 3, 2: 6, 3: 8, 4: 8, 5: 6, 6: 3, 7: 1} >>> [c[i] for i in range(8)] == [prod.coeff(x, i) for i in range(8)] True The generating poly used here is the same as that listed in http://tinyurl.com/cep849r, but in a refactored form. rAr+rBNr_) rQr-poprIrRr,rreversedrrh enumerate) r5ordneedrvniNwasrErevr5rr/r/r0 _AOP_product6s,   (    rHcs*ttrn|dur>sdStfddtdDS|dkrNtdrdt|d|St|Sttrt}|durĈstddt DStfd dt|dDSrt t |S|d|dfvrt S|d|fvrdSt t t |St t|SdS) aReturn the number of combinations of ``n`` items taken ``k`` at a time. Possible values for ``n``: integer - set of length ``n`` sequence - converted to a multiset internally multiset - {element: multiplicity} If ``k`` is None then the total of all combinations of length 0 through the number of items represented in ``n`` will be returned. If ``replacement`` is True then a given item can appear more than once in the ``k`` items. (For example, for 'ab' sets of 2 would include 'aa', 'ab', and 'bb'.) The multiplicity of elements in ``n`` is ignored when ``replacement`` is True but the total number of elements is considered since no element can appear more times than the number of elements in ``n``. Examples ======== >>> from sympy.functions.combinatorial.numbers import nC >>> from sympy.utilities.iterables import multiset_combinations >>> nC(3, 2) 3 >>> nC('abc', 2) 3 >>> nC('aab', 2) 2 When ``replacement`` is True, each item can have multiplicity equal to the length represented by ``n``: >>> nC('aabc', replacement=True) 35 >>> [len(list(multiset_combinations('aaaabbbbcccc', i))) for i in range(5)] [1, 3, 6, 10, 15] >>> sum(_) 35 If there are ``k`` items with multiplicities ``m_1, m_2, ..., m_k`` then the total of all combinations of length 0 through ``k`` is the product, ``(m_1 + 1)*(m_2 + 1)*...*(m_k + 1)``. When the multiplicity of each item is 1 (i.e., k unique items) then there are 2**k combinations. For example, if there are 4 unique items, the total number of combinations is 16: >>> sum(nC(4, i) for i in range(5)) 16 See Also ======== sympy.utilities.iterables.multiset_combinations References ========== .. [1] https://en.wikipedia.org/wiki/Combination .. [2] http://tinyurl.com/cep849r NrAc3s|]}t|VqdSr[nCrCr:r/r0rrGznC..r+rzk cannot be negativecss|]}|dVqdSr+Nr/)rDrr/r/r0rrGc3s|]}t|VqdSr[rIrCr:r/r0rrG)r*rr-r,rJrr)r;rr=rJr<rHtupler0)r5rr7rEr/r:r0rJjs0B     rJcCs||krdkrnntjSd||fvr0tjS||kr>tjS||dkrTt|dS||dkrzd|dt|ddS||dkrt|dt|dSt||S)Nrr+rArHr)rr{rr _stirling1r5rr/r/r0_eval_stirling1s     rOcCsnddgdg|d}td|dD]<}tt||ddD]$}|d||||d||<q:q$t||SNrr+rAr_r,rrr5rrowrErr/r/r0rMs $rMcCs||krdkrnntjSd||fvr0tjS||kr>tjS||dkrTt|dS|dkrbtjS|dkr~td|ddSt||Sr)rr{rrr _stirling2rNr/r/r0_eval_stirling2s   rUcCsjddgdg|d}td|dD]8}tt||ddD] }|||||d||<q:q$t||SrPrQrRr/r/r0rTs  rTrAcCst|}t|}|dkr td||kr.tjS|rLt||d||dS|rhtj||t||S|dkrzt||S|dkrt||Std|dS)a 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: .. math :: {0 \brace 0} = 1; {n \brace 0} = {0 \brace k} = 0; .. math :: {{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 \ge k \ge 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 See Also ======== sympy.utilities.iterables.multiset_partitions References ========== .. [1] https://en.wikipedia.org/wiki/Stirling_numbers_of_the_first_kind .. [2] https://en.wikipedia.org/wiki/Stirling_numbers_of_the_second_kind rzn must be nonnegativer+rAzkind must be 1 or 2, not %sN)r&rJrrrUrirO)r5rr5kindsignedr/r/r0stirlingsV  rXcCs||ks|dkrdS|d|fvr$dS|dkr0dS|dkr@|dS||}|dkrT|Sd||krt|}||dkr|ttd||8}|Sdg|}td|dD]8}t|d|D]}|||||7<q|d8}qdt|d|dS)zxReturn the partitions of ``n`` items into ``k`` parts. This is used by ``nT`` for the case when ``n`` is an integer.rr+rArHN)r r&r-r!r,)r5rr5r1r4rErr/r/r0_nTls,      rYc sttrB|durtSt|trBtt|}tt|Sttsz@tt}|dkrrt t|WS|tkrt |t Wnt yt Yn0t }|dur|dkrdS|d|fvrdS|dks|dkrL|durLt|d\}}tfddt d|dD}|s6|t|d8}|durH|d7}|S|tkrv|durlt|St||St}|dur|tSd}|t|d|D]}|d7}q|S)aIReturn the number of ``k``-sized partitions of ``n`` items. Possible values for ``n``: integer - ``n`` identical items sequence - converted to a multiset internally multiset - {element: multiplicity} Note: the convention for ``nT`` is different than that of ``nC`` and ``nP`` in that here an integer indicates ``n`` *identical* items instead of a set of length ``n``; this is in keeping with the ``partitions`` function which treats its integer-``n`` input like a list of ``n`` 1s. One can use ``range(n)`` for ``n`` to indicate ``n`` distinct items. If ``k`` is None then the total number of ways to partition the elements represented in ``n`` will be returned. Examples ======== >>> from sympy.functions.combinatorial.numbers import nT Partitions of the given multiset: >>> [nT('aabbc', i) for i in range(1, 7)] [1, 8, 11, 5, 1, 0] >>> nT('aabbc') == sum(_) True >>> [nT("mississippi", i) for i in range(1, 12)] [1, 74, 609, 1521, 1768, 1224, 579, 197, 50, 9, 1] Partitions when all items are identical: >>> [nT(5, i) for i in range(1, 6)] [1, 2, 2, 1, 1] >>> nT('1'*5) == sum(_) True When all items are different: >>> [nT(range(5), i) for i in range(1, 6)] [1, 15, 25, 10, 1] >>> nT(range(5)) == sum(_) True Partitions of an integer expressed as a sum of positive integers: >>> from sympy import partition >>> partition(4) 5 >>> nT(4, 1) + nT(4, 2) + nT(4, 3) + nT(4, 4) 5 >>> nT('1'*4) 5 See Also ======== sympy.utilities.iterables.partitions sympy.utilities.iterables.multiset_partitions sympy.functions.combinatorial.numbers.partition References ========== .. [1] http://undergraduate.csse.uwa.edu.au/units/CITS7209/partition.pdf Nr+rAc3s|]}t|VqdSr[rIrCr=r/r0rrGznT..r)r*rr r&rrYr)rRr.nTr, TypeErrorr0r;divmodr-rJr<rrXrcount_partitionsr= enum_range) r5rrrErrrCr1discardr/r=r0rZsRI             rZc@s\eZdZdZeddZeddZeddZeee j e j gdd Z e d d Z d S) motzkina The nth Motzkin number is the number of ways of drawing non-intersecting chords between n points on a circle (not necessarily touching every point by a chord). The Motzkin numbers are named after Theodore Motzkin and have diverse applications in geometry, combinatorics and number theory. Motzkin numbers are the integer sequence defined by the initial terms `M_0 = 1`, `M_1 = 1` and the two-term recurrence relation `M_n = rac{2*n + 1}{n + 2} * M_{n-1} + rac{3n - 3}{n + 2} * M_{n-2}`. Examples ======== >>> from sympy import motzkin >>> motzkin.is_motzkin(5) False >>> motzkin.find_motzkin_numbers_in_range(2,300) [2, 4, 9, 21, 51, 127] >>> motzkin.find_motzkin_numbers_in_range(2,900) [2, 4, 9, 21, 51, 127, 323, 835] >>> motzkin.find_first_n_motzkins(10) [1, 1, 2, 4, 9, 21, 51, 127, 323, 835] References ========== .. [1] https://en.wikipedia.org/wiki/Motzkin_number .. [2] https://mathworld.wolfram.com/MotzkinNumber.html cCsz t|}Wnty YdS0|dkr|dvr6dSd}d}d}||krd|d|d|d||d}|d7}|}|}qB||krdSdSndSdS)NFr)r+rATr+rArH)r&rJ)r5tn1tnrEr-r/r/r0 is_motzkin>s&  (zmotzkin.is_motzkincCsd|kr|krnnt}|dkr2|kr@nn |dd}d}d}||kr||krf||d|d|d|d||d}|d7}|}t|}qL|StddS)Nrr+rArHzEThe provided range is not valid. This condition should satisfy x <= y)rQrSrhrJ)r2rNmotzkinsrarbrEr-r/r/r0find_motzkin_numbers_in_rangeYs   ( z%motzkin.find_motzkin_numbers_in_rangecCsz t|}Wnty&tdYn0|dkr8tddg}|dkrP|dd}d}d}||kr||d|d|d|d||d}|d7}|}t|}q\|S)N.The provided number must be a positive integerrr+rArH)r&rJrSrh)r5rdrarbrEr-r/r/r0find_first_n_motzkinsos&    ( zmotzkin.find_first_n_motzkinscCs0d|d|dd|d|d|dS)NrAr+r_rHr^r/rcr/r/r0_motzkinszmotzkin._motzkincCsJz t|}Wnty&tdYn0|dkr8tdt||dS)Nrfrr+)r&rJrrhr~r/r/r0rns  z motzkin.evalN)rUrVrWrXrYrcrergr%rr{rhr|rnr/r/r/r0r`s$    r`)r5rcsddlm}ddlmddlmdd|||fddkrHtd |durt|t rbt d |slt j St |turt|}t|}n>> from sympy.utilities.iterables import generate_derangements as enum >>> from sympy.functions.combinatorial.numbers import nD A derangement ``d`` of sequence ``s`` has all ``d[i] != s[i]``: >>> set([''.join(i) for i in enum('abc')]) {'bca', 'cab'} >>> nD('abc') 2 Input as iterable or dictionary (multiset form) is accepted: >>> assert nD([1, 2, 2, 3, 3, 3]) == nD({1: 1, 2: 2, 3: 3}) By default, a brute-force enumeration and count of multiset permutations is only done if there are fewer than 9 elements. There may be cases when there is high multiplicty with few unique elements that will benefit from a brute-force enumeration, too. For this reason, the `brute` keyword (default None) is provided. When False, the brute-force enumeration will never be used. When True, it will always be used. >>> nD('1111222233', brute=True) 44 For convenience, one may specify ``n`` distinct items using the ``n`` keyword: >>> assert nD(n=3) == nD('abc') == 2 Since the number of derangments depends on the multiplicity of the elements and not the elements themselves, it may be more convenient to give a list or multiset of multiplicities using keyword ``m``: >>> assert nD('abc') == nD(m=(1,1,1)) == nD(m={1:3}) == 2 r) integrate)laguerrer2cSs&t|tstd|dkr"tddS)Nzexpecting integer valuesrzvalue must not be negativeT)r*rr[rJrkr/r/r0oks  znD..okNrAzenter only 1 of i, n, or mz"items must be a list or dictionaryc3s|]}|VqdSr[r/)rD_rlr/r0rrGznD..cSsi|]\}}|r||qSr/r/rDrr$r/r/r0 rGznD..c3s"|]\}}|o|VqdSr[r/)rDrErrnr/r0rrGcSsi|]\}}||r||qSr/r/ror/r/r0rprGc3s|]}|VqdSr[r/rCrnr/r0rrGcSsg|] }|r|qSr/r/rCr/r/r0rFrGznD..zexpecting iterablecss|]\}}||VqdSr[r/ror/r/r0rrGr+rcss|] }dVqdSrKr/rCr/r/r0r rG)expcsg|]\}}||qSr/r/)rDrEr)rjr2r/r0rF s)&sympy.integrals.integralsri#sympy.functions.special.polynomialsrjZ sympy.abcr2countrJr*rr[rrtyper+rQr"rr,r2r-rRrr$strrhmaxr{rr,rIrr#rrqabsrr)rEZbruter5rrirmsrEcountsZnkeybigZnvalr$rrqr/)rjrlr2r0nDs,                   r|)NF)NF)NF)NrAF)N)NN)drXmathr collectionsrtypingrrtDictrr sympy.corerrr r sympy.core.cacher Zsympy.core.evalfr sympy.core.exprr sympy.core.functionrrsympy.core.logicrZsympy.core.mulrsympy.core.numbersrrrrrsympy.core.relationalrrsympy.external.gmpyr(sympy.functions.combinatorial.factorialsrrrsympy.ntheory.primetestrrsympy.utilities.enumerativersympy.utilities.exceptionsr!sympy.utilities.iterablesr"r#r$Zsympy.utilities.memoizationr%sympy.utilities.miscr&rr'r(Z mpmath.libmpr)r\r1rar6r7rZr}rrrrrrrr!r rLr)r;r<slicer=r0r8r6rHrJrOrMrUrTrXrYrZr`r|r/r/r/r0s               T4Zk`_  B 4 3 ^   l + }