<s> from __future__ import print_function , division <EOL> from sympy . core import S , sympify , Dummy , Mod <EOL> from sympy . core . function import Function , ArgumentIndexError <EOL> from sympy . core . logic import fuzzy_and <EOL> from sympy . core . numbers import Integer , pi <EOL> from sympy . core . relational import Eq <EOL> from sympy . ntheory import sieve <EOL> from math import sqrt as _sqrt <EOL> from sympy . core . compatibility import reduce , range <EOL> from sympy . core . cache import cacheit <EOL> from sympy . polys . polytools import poly_from_expr <EOL> from sympy . polys . polyerrors import PolificationFailed <EOL> class CombinatorialFunction ( Function ) : <EOL> """<STR_LIT>""" <EOL> def _eval_simplify ( self , ratio , measure ) : <EOL> from sympy . simplify . simplify import combsimp <EOL> expr = combsimp ( self ) <EOL> if measure ( expr ) <= ratio * measure ( self ) : <EOL> return expr <EOL> return self <EOL> class factorial ( CombinatorialFunction ) : <EOL> """<STR_LIT>""" <EOL> def fdiff ( self , argindex = <NUM_LIT:1> ) : <EOL> from sympy import gamma , polygamma <EOL> if argindex == <NUM_LIT:1> : <EOL> return gamma ( self . args [ <NUM_LIT:0> ] + <NUM_LIT:1> ) * polygamma ( <NUM_LIT:0> , self . args [ <NUM_LIT:0> ] + <NUM_LIT:1> ) <EOL> else : <EOL> raise ArgumentIndexError ( self , argindex ) <EOL> _small_swing = [ <EOL> <NUM_LIT:1> , <NUM_LIT:1> , <NUM_LIT:1> , <NUM_LIT:3> , <NUM_LIT:3> , <NUM_LIT:15> , <NUM_LIT:5> , <NUM_LIT> , <NUM_LIT> , <NUM_LIT> , <NUM_LIT> , <NUM_LIT> , <NUM_LIT> , <NUM_LIT> , <NUM_LIT> , <NUM_LIT> , <NUM_LIT> , <NUM_LIT> , <EOL> <NUM_LIT> , <NUM_LIT> , <NUM_LIT> , <NUM_LIT> , <NUM_LIT> , <NUM_LIT> , <NUM_LIT> , <NUM_LIT> , <NUM_LIT> , <EOL> <NUM_LIT> , <NUM_LIT> , <NUM_LIT> , <NUM_LIT> , <NUM_LIT> , <NUM_LIT> <EOL> ] <EOL> _small_factorials = [ ] <EOL> @ classmethod <EOL> def _swing ( cls , n ) : <EOL> if n < <NUM_LIT> : <EOL> return cls . _small_swing [ n ] <EOL> else : <EOL> N , primes = int ( _sqrt ( n ) ) , [ ] <EOL> for prime in sieve . primerange ( <NUM_LIT:3> , N + <NUM_LIT:1> ) : <EOL> p , q = <NUM_LIT:1> , n <EOL> while True : <EOL> q //= prime <EOL> if q > <NUM_LIT:0> : <EOL> if q & <NUM_LIT:1> == <NUM_LIT:1> : <EOL> p *= prime <EOL> else : <EOL> break <EOL> if p > <NUM_LIT:1> : <EOL> primes . append ( p ) <EOL> for prime in sieve . primerange ( N + <NUM_LIT:1> , n // <NUM_LIT:3> + <NUM_LIT:1> ) : <EOL> if ( n // prime ) & <NUM_LIT:1> == <NUM_LIT:1> : <EOL> primes . append ( prime ) <EOL> L_product = R_product = <NUM_LIT:1> <EOL> for prime in sieve . primerange ( n // <NUM_LIT:2> + <NUM_LIT:1> , n + <NUM_LIT:1> ) : <EOL> L_product *= prime <EOL> for prime in primes : <EOL> R_product *= prime <EOL> return L_product * R_product <EOL> @ classmethod <EOL> def _recursive ( cls , n ) : <EOL> if n < <NUM_LIT:2> : <EOL> return <NUM_LIT:1> <EOL> else : <EOL> return ( cls . _recursive ( n // <NUM_LIT:2> ) ** <NUM_LIT:2> ) * cls . _swing ( n ) <EOL> @ classmethod <EOL> def eval ( cls , n ) : <EOL> n = sympify ( n ) <EOL> if n . is_Number : <EOL> if n is S . Zero : <EOL> return S . One <EOL> elif n is S . Infinity : <EOL> return S . Infinity <EOL> elif n . is_Integer : <EOL> if n . is_negative : <EOL> return S . ComplexInfinity <EOL> else : <EOL> n = n . p <EOL> if n < <NUM_LIT:20> : <EOL> if not cls . _small_factorials : <EOL> result = <NUM_LIT:1> <EOL> for i in range ( <NUM_LIT:1> , <NUM_LIT:20> ) : <EOL> result *= i <EOL> cls . _small_factorials . append ( result ) <EOL> result = cls . _small_factorials [ n - <NUM_LIT:1> ] <EOL> else : <EOL> bits = bin ( n ) . count ( '<STR_LIT:1>' ) <EOL> result = cls . _recursive ( n ) * <NUM_LIT:2> ** ( n - bits ) <EOL> return Integer ( result ) <EOL> def _eval_rewrite_as_gamma ( self , n ) : <EOL> from sympy import gamma <EOL> return gamma ( n + <NUM_LIT:1> ) <EOL> def _eval_rewrite_as_Product ( self , n ) : <EOL> from sympy import Product <EOL> if n . is_nonnegative and n . is_integer : <EOL> i = Dummy ( '<STR_LIT:i>' , integer = True ) <EOL> return Product ( i , ( i , <NUM_LIT:1> , n ) ) <EOL> def _eval_is_integer ( self ) : <EOL> if self . args [ <NUM_LIT:0> ] . is_integer and self . args [ <NUM_LIT:0> ] . is_nonnegative : <EOL> return True <EOL> def _eval_is_positive ( self ) : <EOL> if self . args [ <NUM_LIT:0> ] . is_integer and self . args [ <NUM_LIT:0> ] . is_nonnegative : <EOL> return True <EOL> def _eval_is_composite ( self ) : <EOL> x = self . args [ <NUM_LIT:0> ] <EOL> if x . is_integer : <EOL> return ( x - <NUM_LIT:3> ) . is_nonnegative <EOL> def _eval_is_real ( self ) : <EOL> x = self . args [ <NUM_LIT:0> ] <EOL> if x . is_nonnegative or x . is_noninteger : <EOL> return True <EOL> class MultiFactorial ( CombinatorialFunction ) : <EOL> pass <EOL> class subfactorial ( CombinatorialFunction ) : <EOL> r"""<STR_LIT>""" <EOL> @ classmethod <EOL> @ cacheit <EOL> def _eval ( self , n ) : <EOL> if not n : <EOL> return S . One <EOL> elif n == <NUM_LIT:1> : <EOL> return S . Zero <EOL> return ( n - <NUM_LIT:1> ) * ( self . _eval ( n - <NUM_LIT:1> ) + self . _eval ( n - <NUM_LIT:2> ) ) <EOL> @ classmethod <EOL> def eval ( cls , arg ) : <EOL> if arg . is_Number : <EOL> if arg . is_Integer and arg . is_nonnegative : <EOL> return cls . _eval ( arg ) <EOL> elif arg is S . NaN : <EOL> return S . NaN <EOL> elif arg is S . Infinity : <EOL> return S . Infinity <EOL> def _eval_is_even ( self ) : <EOL> if self . args [ <NUM_LIT:0> ] . is_odd and self . args [ <NUM_LIT:0> ] . is_nonnegative : <EOL> return True <EOL> def _eval_is_integer ( self ) : <EOL> if self . args [ <NUM_LIT:0> ] . is_integer and self . args [ <NUM_LIT:0> ] . is_nonnegative : <EOL> return True <EOL> def _eval_rewrite_as_uppergamma ( self , arg ) : <EOL> from sympy import uppergamma <EOL> return uppergamma ( arg + <NUM_LIT:1> , - <NUM_LIT:1> ) / S . Exp1 <EOL> def _eval_is_nonnegative ( self ) : <EOL> if self . args [ <NUM_LIT:0> ] . is_integer and self . args [ <NUM_LIT:0> ] . is_nonnegative : <EOL> return True <EOL> def _eval_is_odd ( self ) : <EOL> if self . args [ <NUM_LIT:0> ] . is_even and self . args [ <NUM_LIT:0> ] . is_nonnegative : <EOL> return True <EOL> class factorial2 ( CombinatorialFunction ) : <EOL> """<STR_LIT>""" <EOL> @ classmethod <EOL> def eval ( cls , arg ) : <EOL> if arg . is_Number : <EOL> if not arg . is_Integer : <EOL> raise ValueError ( "<STR_LIT>" ) <EOL> if arg . is_nonnegative : <EOL> if arg . is_even : <EOL> k = arg / <NUM_LIT:2> <EOL> return <NUM_LIT:2> ** k * factorial ( k ) <EOL> return factorial ( arg ) / factorial2 ( arg - <NUM_LIT:1> ) <EOL> if arg . is_odd : <EOL> return arg * ( S . NegativeOne ) ** ( ( <NUM_LIT:1> - arg ) / <NUM_LIT:2> ) / factorial2 ( - arg ) <EOL> raise ValueError ( "<STR_LIT>" ) <EOL> def _eval_is_even ( self ) : <EOL> n = self . args [ <NUM_LIT:0> ] <EOL> if n . is_integer : <EOL> if n . is_odd : <EOL> return False <EOL> if n . is_even : <EOL> if n . is_positive : <EOL> return True <EOL> if n . is_zero : <EOL> return False <EOL> def _eval_is_integer ( self ) : <EOL> n = self . args [ <NUM_LIT:0> ] <EOL> if n . is_integer : <EOL> if ( n + <NUM_LIT:1> ) . is_nonnegative : <EOL> return True <EOL> if n . is_odd : <EOL> return ( n + <NUM_LIT:3> ) . is_nonnegative <EOL> def _eval_is_odd ( self ) : <EOL> n = self . args [ <NUM_LIT:0> ] <EOL> if n . is_odd : <EOL> return ( n + <NUM_LIT:3> ) . is_nonnegative <EOL> if n . is_even : <EOL> if n . is_positive : <EOL> return False <EOL> if n . is_zero : <EOL> return True <EOL> def _eval_is_positive ( self ) : <EOL> n = self . args [ <NUM_LIT:0> ] <EOL> if n . is_integer : <EOL> if ( n + <NUM_LIT:1> ) . is_nonnegative : <EOL> return True <EOL> if n . is_odd : <EOL> return ( ( n + <NUM_LIT:1> ) / <NUM_LIT:2> ) . is_even <EOL> def _eval_rewrite_as_gamma ( self , n ) : <EOL> from sympy import gamma , Piecewise , sqrt <EOL> return <NUM_LIT:2> ** ( n / <NUM_LIT:2> ) * gamma ( n / <NUM_LIT:2> + <NUM_LIT:1> ) * Piecewise ( ( <NUM_LIT:1> , Eq ( Mod ( n , <NUM_LIT:2> ) , <NUM_LIT:0> ) ) , ( sqrt ( <NUM_LIT:2> / pi ) , Eq ( Mod ( n , <NUM_LIT:2> ) , <NUM_LIT:1> ) ) ) <EOL> class RisingFactorial ( CombinatorialFunction ) : <EOL> """<STR_LIT>""" <EOL> @ classmethod <EOL> def eval ( cls , x , k ) : <EOL> x = sympify ( x ) <EOL> k = sympify ( k ) <EOL> if x is S . NaN or k is S . NaN : <EOL> return S . NaN <EOL> elif x is S . One : <EOL> return factorial ( k ) <EOL> elif k . is_Integer : <EOL> if k is S . Zero : <EOL> return S . One <EOL> else : <EOL> if k . is_positive : <EOL> if x is S . Infinity : <EOL> return S . Infinity <EOL> elif x is S . NegativeInfinity : <EOL> if k . is_odd : <EOL> return S . NegativeInfinity <EOL> else : <EOL> return S . Infinity <EOL> else : <EOL> try : <EOL> F , opt = poly_from_expr ( x ) <EOL> except PolificationFailed : <EOL> return reduce ( lambda r , i : r * ( x + i ) , range ( <NUM_LIT:0> , int ( k ) ) , <NUM_LIT:1> ) <EOL> if len ( opt . gens ) > <NUM_LIT:1> or F . degree ( ) <= <NUM_LIT:1> : <EOL> return reduce ( lambda r , i : r * ( x + i ) , range ( <NUM_LIT:0> , int ( k ) ) , <NUM_LIT:1> ) <EOL> else : <EOL> v = opt . gens [ <NUM_LIT:0> ] <EOL> return reduce ( lambda r , i : <EOL> r * ( F . subs ( v , v + i ) . expand ( ) ) , <EOL> range ( <NUM_LIT:0> , int ( k ) ) , <NUM_LIT:1> ) <EOL> else : <EOL> if x is S . Infinity : <EOL> return S . Infinity <EOL> elif x is S . NegativeInfinity : <EOL> return S . Infinity <EOL> else : <EOL> try : <EOL> F , opt = poly_from_expr ( x ) <EOL> except PolificationFailed : <EOL> return <NUM_LIT:1> / reduce ( lambda r , i : <EOL> r * ( x - i ) , <EOL> range ( <NUM_LIT:1> , abs ( int ( k ) ) + <NUM_LIT:1> ) , <NUM_LIT:1> ) <EOL> if len ( opt . gens ) > <NUM_LIT:1> or F . degree ( ) <= <NUM_LIT:1> : <EOL> return <NUM_LIT:1> / reduce ( lambda r , i : <EOL> r * ( x - i ) , <EOL> range ( <NUM_LIT:1> , abs ( int ( k ) ) + <NUM_LIT:1> ) , <NUM_LIT:1> ) <EOL> else : <EOL> v = opt . gens [ <NUM_LIT:0> ] <EOL> return <NUM_LIT:1> / reduce ( lambda r , i : <EOL> r * ( F . subs ( v , v - i ) . expand ( ) ) , <EOL> range ( <NUM_LIT:1> , abs ( int ( k ) ) + <NUM_LIT:1> ) , <NUM_LIT:1> ) <EOL> def _eval_rewrite_as_gamma ( self , x , k ) : <EOL> from sympy import gamma <EOL> return gamma ( x + k ) / gamma ( x ) <EOL> def _eval_rewrite_as_FallingFactorial ( self , x , k ) : <EOL> return FallingFactorial ( x + k - <NUM_LIT:1> , k ) <EOL> def _eval_rewrite_as_factorial ( self , x , k ) : <EOL> if x . is_integer and k . is_integer : <EOL> return factorial ( k + x - <NUM_LIT:1> ) / factorial ( x - <NUM_LIT:1> ) <EOL> def _eval_rewrite_as_binomial ( self , x , k ) : <EOL> if k . is_integer : <EOL> return factorial ( k ) * binomial ( x + k - <NUM_LIT:1> , k ) <EOL> def _eval_is_integer ( self ) : <EOL> return fuzzy_and ( ( self . args [ <NUM_LIT:0> ] . is_integer , self . args [ <NUM_LIT:1> ] . is_integer , <EOL> self . args [ <NUM_LIT:1> ] . is_nonnegative ) ) <EOL> def _sage_ ( self ) : <EOL> import sage . all as sage <EOL> return sage . rising_factorial ( self . args [ <NUM_LIT:0> ] . _sage_ ( ) , self . args [ <NUM_LIT:1> ] . _sage_ ( ) ) <EOL> class FallingFactorial ( CombinatorialFunction ) : <EOL> """<STR_LIT>""" <EOL> @ classmethod <EOL> def eval ( cls , x , k ) : <EOL> x = sympify ( x ) <EOL> k = sympify ( k ) <EOL> if x is S . NaN or k is S . NaN : <EOL> return S . NaN <EOL> elif k . is_integer and x == k : <EOL> return factorial ( x ) <EOL> elif k . is_Integer : <EOL> if k is S . Zero : <EOL> return S . One <EOL> else : <EOL> if k . is_positive : <EOL> if x is S . Infinity : <EOL> return S . Infinity <EOL> elif x is S . NegativeInfinity : <EOL> if k . is_odd : <EOL> return S . NegativeInfinity <EOL> else : <EOL> return S . Infinity <EOL> else : <EOL> try : <EOL> F , opt = poly_from_expr ( x ) <EOL> except PolificationFailed : <EOL> return reduce ( lambda r , i : r * ( x - i ) , <EOL> range ( <NUM_LIT:0> , int ( k ) ) , <NUM_LIT:1> ) <EOL> if len ( opt . gens ) > <NUM_LIT:1> or F . degree ( ) <= <NUM_LIT:1> : <EOL> return reduce ( lambda r , i : r * ( x - i ) , <EOL> range ( <NUM_LIT:0> , int ( k ) ) , <NUM_LIT:1> ) <EOL> else : <EOL> v = opt . gens [ <NUM_LIT:0> ] <EOL> return reduce ( lambda r , i : <EOL> r * ( F . subs ( v , v - i ) . expand ( ) ) , <EOL> range ( <NUM_LIT:0> , int ( k ) ) , <NUM_LIT:1> ) <EOL> else : <EOL> if x is S . Infinity : <EOL> return S . Infinity <EOL> elif x is S . NegativeInfinity : <EOL> return S . Infinity <EOL> else : <EOL> try : <EOL> F , opt = poly_from_expr ( x ) <EOL> except PolificationFailed : <EOL> return <NUM_LIT:1> / reduce ( lambda r , i : r * ( x + i ) , <EOL> range ( <NUM_LIT:1> , abs ( int ( k ) ) + <NUM_LIT:1> ) , <NUM_LIT:1> ) <EOL> if len ( opt . gens ) > <NUM_LIT:1> or F . degree ( ) <= <NUM_LIT:1> : <EOL> return <NUM_LIT:1> / reduce ( lambda r , i : r * ( x + i ) , <EOL> range ( <NUM_LIT:1> , abs ( int ( k ) ) + <NUM_LIT:1> ) , <NUM_LIT:1> ) <EOL> else : <EOL> v = opt . gens [ <NUM_LIT:0> ] <EOL> return <NUM_LIT:1> / reduce ( lambda r , i : <EOL> r * ( F . subs ( v , v + i ) . expand ( ) ) , <EOL> range ( <NUM_LIT:1> , abs ( int ( k ) ) + <NUM_LIT:1> ) , <NUM_LIT:1> ) <EOL> def _eval_rewrite_as_gamma ( self , x , k ) : <EOL> from sympy import gamma <EOL> return ( - <NUM_LIT:1> ) ** k * gamma ( k - x ) / gamma ( - x ) <EOL> def _eval_rewrite_as_RisingFactorial ( self , x , k ) : <EOL> return rf ( x - k + <NUM_LIT:1> , k ) <EOL> def _eval_rewrite_as_binomial ( self , x , k ) : <EOL> if k . is_integer : <EOL> return factorial ( k ) * binomial ( x , k ) <EOL> def _eval_rewrite_as_factorial ( self , x , k ) : <EOL> if x . is_integer and k . is_integer : <EOL> return factorial ( x ) / factorial ( x - k ) <EOL> def _eval_is_integer ( self ) : <EOL> return fuzzy_and ( ( self . args [ <NUM_LIT:0> ] . is_integer , self . args [ <NUM_LIT:1> ] . is_integer , <EOL> self . args [ <NUM_LIT:1> ] . is_nonnegative ) ) <EOL> def _sage_ ( self ) : <EOL> import sage . all as sage <EOL> return sage . falling_factorial ( self . args [ <NUM_LIT:0> ] . _sage_ ( ) , <EOL> self . args [ <NUM_LIT:1> ] . _sage_ ( ) ) <EOL> rf = RisingFactorial <EOL> ff = FallingFactorial <EOL> class binomial ( CombinatorialFunction ) : <EOL> """<STR_LIT>""" <EOL> def fdiff ( self , argindex = <NUM_LIT:1> ) : <EOL> from sympy import polygamma <EOL> if argindex == <NUM_LIT:1> : <EOL> n , k = self . args <EOL> return binomial ( n , k ) * ( polygamma ( <NUM_LIT:0> , n + <NUM_LIT:1> ) - polygamma ( <NUM_LIT:0> , n - k + <NUM_LIT:1> ) ) <EOL> elif argindex == <NUM_LIT:2> : <EOL> n , k = self . args <EOL> return binomial ( n , k ) * ( polygamma ( <NUM_LIT:0> , n - k + <NUM_LIT:1> ) - polygamma ( <NUM_LIT:0> , k + <NUM_LIT:1> ) ) <EOL> else : <EOL> raise ArgumentIndexError ( self , argindex ) <EOL> @ classmethod <EOL> def _eval ( self , n , k ) : <EOL> if k . is_Integer : <EOL> if n . is_Integer and n >= <NUM_LIT:0> : <EOL> n , k = int ( n ) , int ( k ) <EOL> if k > n : <EOL> return S . Zero <EOL> elif k > n // <NUM_LIT:2> : <EOL> k = n - k <EOL> M , result = int ( _sqrt ( n ) ) , <NUM_LIT:1> <EOL> for prime in sieve . primerange ( <NUM_LIT:2> , n + <NUM_LIT:1> ) : <EOL> if prime > n - k : <EOL> result *= prime <EOL> elif prime > n // <NUM_LIT:2> : <EOL> continue <EOL> elif prime > M : <EOL> if n % prime < k % prime : <EOL> result *= prime <EOL> else : <EOL> N , K = n , k <EOL> exp = a = <NUM_LIT:0> <EOL> while N > <NUM_LIT:0> : <EOL> a = int ( ( N % prime ) < ( K % prime + a ) ) <EOL> N , K = N // prime , K // prime <EOL> exp = a + exp <EOL> if exp > <NUM_LIT:0> : <EOL> result *= prime ** exp <EOL> return Integer ( result ) <EOL> else : <EOL> d = result = n - k + <NUM_LIT:1> <EOL> for i in range ( <NUM_LIT:2> , k + <NUM_LIT:1> ) : <EOL> d += <NUM_LIT:1> <EOL> result *= d <EOL> result /= i <EOL> return result <EOL> @ classmethod <EOL> def eval ( cls , n , k ) : <EOL> n , k = map ( sympify , ( n , k ) ) <EOL> d = n - k <EOL> if d . is_zero or k . is_zero : <EOL> return S . One <EOL> elif d . is_zero is False : <EOL> if ( k - <NUM_LIT:1> ) . is_zero : <EOL> return n <EOL> elif k . is_negative : <EOL> return S . Zero <EOL> elif n . is_integer and n . is_nonnegative and d . is_negative : <EOL> return S . Zero <EOL> if k . is_Integer and k > <NUM_LIT:0> and n . is_Number : <EOL> return cls . _eval ( n , k ) <EOL> def _eval_expand_func ( self , ** hints ) : <EOL> """<STR_LIT>""" <EOL> n = self . args [ <NUM_LIT:0> ] <EOL> if n . is_Number : <EOL> return binomial ( * self . args ) <EOL> k = self . args [ <NUM_LIT:1> ] <EOL> if k . is_Add and n in k . args : <EOL> k = n - k <EOL> if k . is_Integer : <EOL> if k == S . Zero : <EOL> return S . One <EOL> elif k < <NUM_LIT:0> : <EOL> return S . Zero <EOL> else : <EOL> n = self . args [ <NUM_LIT:0> ] <EOL> result = n - k + <NUM_LIT:1> <EOL> for i in range ( <NUM_LIT:2> , k + <NUM_LIT:1> ) : <EOL> result *= n - k + i <EOL> result /= i <EOL> return result <EOL> else : <EOL> return binomial ( * self . args ) <EOL> def _eval_rewrite_as_factorial ( self , n , k ) : <EOL> return factorial ( n ) / ( factorial ( k ) * factorial ( n - k ) ) <EOL> def _eval_rewrite_as_gamma ( self , n , k ) : <EOL> from sympy import gamma <EOL> return gamma ( n + <NUM_LIT:1> ) / ( gamma ( k + <NUM_LIT:1> ) * gamma ( n - k + <NUM_LIT:1> ) ) <EOL> def _eval_rewrite_as_tractable ( self , n , k ) : <EOL> return self . _eval_rewrite_as_gamma ( n , k ) . rewrite ( '<STR_LIT>' ) <EOL> def _eval_rewrite_as_FallingFactorial ( self , n , k ) : <EOL> if k . is_integer : <EOL> return ff ( n , k ) / factorial ( k ) <EOL> def _eval_is_integer ( self ) : <EOL> n , k = self . args <EOL> if n . is_integer and k . is_integer : <EOL> return True <EOL> elif k . is_integer is False : <EOL> return False </s>
