Concrete¶
Hypergeometric terms¶
The center stage, in recurrence solving and summations, play hypergeometric terms. Formally these are sequences annihilated by first order linear recurrence operators. In simple words if we are given term \(a(n)\) then it is hypergeometric if its consecutive term ratio is a rational function in \(n\).
To check if a sequence is of this type you can use the is_hypergeometric
method which is available in Basic class. Here is simple example involving a
polynomial:
>>> from sympy import *
>>> n, k = symbols('n,k')
>>> (n**2 + 1).is_hypergeometric(n)
True
Of course polynomials are hypergeometric but are there any more complicated sequences of this type? Here are some trivial examples:
>>> factorial(n).is_hypergeometric(n)
True
>>> binomial(n, k).is_hypergeometric(n)
True
>>> rf(n, k).is_hypergeometric(n)
True
>>> ff(n, k).is_hypergeometric(n)
True
>>> gamma(n).is_hypergeometric(n)
True
>>> (2**n).is_hypergeometric(n)
True
We see that all species used in summations and other parts of concrete mathematics are hypergeometric. Note also that binomial coefficients and both rising and falling factorials are hypergeometric in both their arguments:
>>> binomial(n, k).is_hypergeometric(k)
True
>>> rf(n, k).is_hypergeometric(k)
True
>>> ff(n, k).is_hypergeometric(k)
True
To say more, all previously shown examples are valid for integer linear arguments:
>>> factorial(2*n).is_hypergeometric(n)
True
>>> binomial(3*n+1, k).is_hypergeometric(n)
True
>>> rf(n+1, k1).is_hypergeometric(n)
True
>>> ff(n1, k+1).is_hypergeometric(n)
True
>>> gamma(5*n).is_hypergeometric(n)
True
>>> (2**(n7)).is_hypergeometric(n)
True
However nonlinear arguments make those sequences fail to be hypergeometric:
>>> factorial(n**2).is_hypergeometric(n)
False
>>> (2**(n**3 + 1)).is_hypergeometric(n)
False
If not only the knowledge of being hypergeometric or not is needed, you can use
hypersimp()
function. It will try to simplify combinatorial expression and
if the term given is hypergeometric it will return a quotient of polynomials of
minimal degree. Otherwise is will return \(None\) to say that sequence is not
hypergeometric:
>>> hypersimp(factorial(2*n), n)
2*(n + 1)*(2*n + 1)
>>> hypersimp(factorial(n**2), n)
Concrete Class Reference¶

class
sympy.concrete.summations.
Sum
[source]¶ Represents unevaluated summation.
Sum
represents a finite or infinite series, with the first argument being the general form of terms in the series, and the second argument being(dummy_variable, start, end)
, withdummy_variable
taking all integer values fromstart
throughend
. In accordance with longstanding mathematical convention, the end term is included in the summation.Finite Sums
For finite sums (and sums with symbolic limits assumed to be finite) we follow the summation convention described by Karr [1], especially definition 3 of section 1.4. The sum:
\[\sum_{m \leq i < n} f(i)\]has the obvious meaning for \(m < n\), namely:
\[\sum_{m \leq i < n} f(i) = f(m) + f(m+1) + \ldots + f(n2) + f(n1)\]with the upper limit value \(f(n)\) excluded. The sum over an empty set is zero if and only if \(m = n\):
\[\sum_{m \leq i < n} f(i) = 0 \quad \mathrm{for} \quad m = n\]Finally, for all other sums over empty sets we assume the following definition:
\[\sum_{m \leq i < n} f(i) =  \sum_{n \leq i < m} f(i) \quad \mathrm{for} \quad m > n\]It is important to note that Karr defines all sums with the upper limit being exclusive. This is in contrast to the usual mathematical notation, but does not affect the summation convention. Indeed we have:
\[\sum_{m \leq i < n} f(i) = \sum_{i = m}^{n  1} f(i)\]where the difference in notation is intentional to emphasize the meaning, with limits typeset on the top being inclusive.
Examples
>>> from sympy.abc import i, k, m, n, x >>> from sympy import Sum, factorial, oo, IndexedBase, Function >>> Sum(k, (k, 1, m)) Sum(k, (k, 1, m)) >>> Sum(k, (k, 1, m)).doit() m**2/2 + m/2 >>> Sum(k**2, (k, 1, m)) Sum(k**2, (k, 1, m)) >>> Sum(k**2, (k, 1, m)).doit() m**3/3 + m**2/2 + m/6 >>> Sum(x**k, (k, 0, oo)) Sum(x**k, (k, 0, oo)) >>> Sum(x**k, (k, 0, oo)).doit() Piecewise((1/(1  x), Abs(x) < 1), (Sum(x**k, (k, 0, oo)), True)) >>> Sum(x**k/factorial(k), (k, 0, oo)).doit() exp(x)
Here are examples to do summation with symbolic indices. You can use either Function of IndexedBase classes:
>>> f = Function('f') >>> Sum(f(n), (n, 0, 3)).doit() f(0) + f(1) + f(2) + f(3) >>> Sum(f(n), (n, 0, oo)).doit() Sum(f(n), (n, 0, oo)) >>> f = IndexedBase('f') >>> Sum(f[n]**2, (n, 0, 3)).doit() f[0]**2 + f[1]**2 + f[2]**2 + f[3]**2
An example showing that the symbolic result of a summation is still valid for seemingly nonsensical values of the limits. Then the Karr convention allows us to give a perfectly valid interpretation to those sums by interchanging the limits according to the above rules:
>>> S = Sum(i, (i, 1, n)).doit() >>> S n**2/2 + n/2 >>> S.subs(n, 4) 6 >>> Sum(i, (i, 1, 4)).doit() 6 >>> Sum(i, (i, 3, 0)).doit() 6
An explicit example of the Karr summation convention:
>>> S1 = Sum(i**2, (i, m, m+n1)).doit() >>> S1 m**2*n + m*n**2  m*n + n**3/3  n**2/2 + n/6 >>> S2 = Sum(i**2, (i, m+n, m1)).doit() >>> S2 m**2*n  m*n**2 + m*n  n**3/3 + n**2/2  n/6 >>> S1 + S2 0 >>> S3 = Sum(i, (i, m, m1)).doit() >>> S3 0
See also
summation
,Product
,product
References
 R67
Michael Karr, “Summation in Finite Terms”, Journal of the ACM, Volume 28 Issue 2, April 1981, Pages 305350 http://dl.acm.org/citation.cfm?doid=322248.322255
 R68
https://en.wikipedia.org/wiki/Summation#Capitalsigma_notation
 R69

euler_maclaurin
(m=0, n=0, eps=0, eval_integral=True)[source]¶ Return an EulerMaclaurin approximation of self, where m is the number of leading terms to sum directly and n is the number of terms in the tail.
With m = n = 0, this is simply the corresponding integral plus a firstorder endpoint correction.
Returns (s, e) where s is the EulerMaclaurin approximation and e is the estimated error (taken to be the magnitude of the first omitted term in the tail):
>>> from sympy.abc import k, a, b >>> from sympy import Sum >>> Sum(1/k, (k, 2, 5)).doit().evalf() 1.28333333333333 >>> s, e = Sum(1/k, (k, 2, 5)).euler_maclaurin() >>> s log(2) + 7/20 + log(5) >>> from sympy import sstr >>> print(sstr((s.evalf(), e.evalf()), full_prec=True)) (1.26629073187415, 0.0175000000000000)
The endpoints may be symbolic:
>>> s, e = Sum(1/k, (k, a, b)).euler_maclaurin() >>> s log(a) + log(b) + 1/(2*b) + 1/(2*a) >>> e Abs(1/(12*b**2)  1/(12*a**2))
If the function is a polynomial of degree at most 2n+1, the EulerMaclaurin formula becomes exact (and e = 0 is returned):
>>> Sum(k, (k, 2, b)).euler_maclaurin() (b**2/2 + b/2  1, 0) >>> Sum(k, (k, 2, b)).doit() b**2/2 + b/2  1
With a nonzero eps specified, the summation is ended as soon as the remainder term is less than the epsilon.

eval_zeta_function
(f, limits)[source]¶ Check whether the function matches with the zeta function. If it matches, then return a \(Piecewise\) expression because zeta function does not converge unless \(s > 1\) and \(q > 0\)

is_absolutely_convergent
()[source]¶ Checks for the absolute convergence of an infinite series.
Same as checking convergence of absolute value of sequence_term of an infinite series.
Examples
>>> from sympy import Sum, Symbol, sin, oo >>> n = Symbol('n', integer=True) >>> Sum((1)**n, (n, 1, oo)).is_absolutely_convergent() False >>> Sum((1)**n/n**2, (n, 1, oo)).is_absolutely_convergent() True
See also
References

is_convergent
()[source]¶ Checks for the convergence of a Sum.
We divide the study of convergence of infinite sums and products in two parts.
First Part: One part is the question whether all the terms are well defined, i.e., they are finite in a sum and also nonzero in a product. Zero is the analogy of (minus) infinity in products as \(e^{\infty} = 0\).
Second Part: The second part is the question of convergence after infinities, and zeros in products, have been omitted assuming that their number is finite. This means that we only consider the tail of the sum or product, starting from some point after which all terms are well defined.
For example, in a sum of the form:
\[\sum_{1 \leq i < \infty} \frac{1}{n^2 + an + b}\]where a and b are numbers. The routine will return true, even if there are infinities in the term sequence (at most two). An analogous product would be:
\[\prod_{1 \leq i < \infty} e^{\frac{1}{n^2 + an + b}}\]This is how convergence is interpreted. It is concerned with what happens at the limit. Finding the bad terms is another independent matter.
Note: It is responsibility of user to see that the sum or product is well defined.
There are various tests employed to check the convergence like divergence test, root test, integral test, alternating series test, comparison tests, Dirichlet tests. It returns true if Sum is convergent and false if divergent and NotImplementedError if it can not be checked.
Examples
>>> from sympy import factorial, S, Sum, Symbol, oo >>> n = Symbol('n', integer=True) >>> Sum(n/(n  1), (n, 4, 7)).is_convergent() True >>> Sum(n/(2*n + 1), (n, 1, oo)).is_convergent() False >>> Sum(factorial(n)/5**n, (n, 1, oo)).is_convergent() False >>> Sum(1/n**(S(6)/5), (n, 1, oo)).is_convergent() True
See also
Sum.is_absolutely_convergent
,Product.is_convergent
References

reverse_order
(*indices)[source]¶ Reverse the order of a limit in a Sum.
Usage
reverse_order(self, *indices)
reverses some limits in the expressionself
which can be either aSum
or aProduct
. The selectors in the argumentindices
specify some indices whose limits get reversed. These selectors are either variable names or numerical indices counted starting from the innermost limit tuple.Examples
>>> from sympy import Sum >>> from sympy.abc import x, y, a, b, c, d
>>> Sum(x, (x, 0, 3)).reverse_order(x) Sum(x, (x, 4, 1)) >>> Sum(x*y, (x, 1, 5), (y, 0, 6)).reverse_order(x, y) Sum(x*y, (x, 6, 0), (y, 7, 1)) >>> Sum(x, (x, a, b)).reverse_order(x) Sum(x, (x, b + 1, a  1)) >>> Sum(x, (x, a, b)).reverse_order(0) Sum(x, (x, b + 1, a  1))
While one should prefer variable names when specifying which limits to reverse, the index counting notation comes in handy in case there are several symbols with the same name.
>>> S = Sum(x**2, (x, a, b), (x, c, d)) >>> S Sum(x**2, (x, a, b), (x, c, d)) >>> S0 = S.reverse_order(0) >>> S0 Sum(x**2, (x, b + 1, a  1), (x, c, d)) >>> S1 = S0.reverse_order(1) >>> S1 Sum(x**2, (x, b + 1, a  1), (x, d + 1, c  1))
Of course we can mix both notations:
>>> Sum(x*y, (x, a, b), (y, 2, 5)).reverse_order(x, 1) Sum(x*y, (x, b + 1, a  1), (y, 6, 1)) >>> Sum(x*y, (x, a, b), (y, 2, 5)).reverse_order(y, x) Sum(x*y, (x, b + 1, a  1), (y, 6, 1))
See also
index
,reorder_limit
,reorder
References
 R72
Michael Karr, “Summation in Finite Terms”, Journal of the ACM, Volume 28 Issue 2, April 1981, Pages 305350 http://dl.acm.org/citation.cfm?doid=322248.322255

class
sympy.concrete.products.
Product
[source]¶ Represents unevaluated products.
Product
represents a finite or infinite product, with the first argument being the general form of terms in the series, and the second argument being(dummy_variable, start, end)
, withdummy_variable
taking all integer values fromstart
throughend
. In accordance with longstanding mathematical convention, the end term is included in the product.Finite Products
For finite products (and products with symbolic limits assumed to be finite) we follow the analogue of the summation convention described by Karr [1], especially definition 3 of section 1.4. The product:
\[\prod_{m \leq i < n} f(i)\]has the obvious meaning for \(m < n\), namely:
\[\prod_{m \leq i < n} f(i) = f(m) f(m+1) \cdot \ldots \cdot f(n2) f(n1)\]with the upper limit value \(f(n)\) excluded. The product over an empty set is one if and only if \(m = n\):
\[\prod_{m \leq i < n} f(i) = 1 \quad \mathrm{for} \quad m = n\]Finally, for all other products over empty sets we assume the following definition:
\[\prod_{m \leq i < n} f(i) = \frac{1}{\prod_{n \leq i < m} f(i)} \quad \mathrm{for} \quad m > n\]It is important to note that above we define all products with the upper limit being exclusive. This is in contrast to the usual mathematical notation, but does not affect the product convention. Indeed we have:
\[\prod_{m \leq i < n} f(i) = \prod_{i = m}^{n  1} f(i)\]where the difference in notation is intentional to emphasize the meaning, with limits typeset on the top being inclusive.
Examples
>>> from sympy.abc import a, b, i, k, m, n, x >>> from sympy import Product, factorial, oo >>> Product(k, (k, 1, m)) Product(k, (k, 1, m)) >>> Product(k, (k, 1, m)).doit() factorial(m) >>> Product(k**2,(k, 1, m)) Product(k**2, (k, 1, m)) >>> Product(k**2,(k, 1, m)).doit() factorial(m)**2
Wallis’ product for pi:
>>> W = Product(2*i/(2*i1) * 2*i/(2*i+1), (i, 1, oo)) >>> W Product(4*i**2/((2*i  1)*(2*i + 1)), (i, 1, oo))
Direct computation currently fails:
>>> W.doit() Product(4*i**2/((2*i  1)*(2*i + 1)), (i, 1, oo))
But we can approach the infinite product by a limit of finite products:
>>> from sympy import limit >>> W2 = Product(2*i/(2*i1)*2*i/(2*i+1), (i, 1, n)) >>> W2 Product(4*i**2/((2*i  1)*(2*i + 1)), (i, 1, n)) >>> W2e = W2.doit() >>> W2e 2**(2*n)*4**n*factorial(n)**2/(RisingFactorial(1/2, n)*RisingFactorial(3/2, n)) >>> limit(W2e, n, oo) pi/2
By the same formula we can compute sin(pi/2):
>>> from sympy import pi, gamma, simplify >>> P = pi * x * Product(1  x**2/k**2, (k, 1, n)) >>> P = P.subs(x, pi/2) >>> P pi**2*Product(1  pi**2/(4*k**2), (k, 1, n))/2 >>> Pe = P.doit() >>> Pe pi**2*RisingFactorial(1  pi/2, n)*RisingFactorial(1 + pi/2, n)/(2*factorial(n)**2) >>> Pe = Pe.rewrite(gamma) >>> Pe pi**2*gamma(n + 1 + pi/2)*gamma(n  pi/2 + 1)/(2*gamma(1  pi/2)*gamma(1 + pi/2)*gamma(n + 1)**2) >>> Pe = simplify(Pe) >>> Pe sin(pi**2/2)*gamma(n + 1 + pi/2)*gamma(n  pi/2 + 1)/gamma(n + 1)**2 >>> limit(Pe, n, oo) sin(pi**2/2)
Products with the lower limit being larger than the upper one:
>>> Product(1/i, (i, 6, 1)).doit() 120 >>> Product(i, (i, 2, 5)).doit() 120
The empty product:
>>> Product(i, (i, n, n1)).doit() 1
An example showing that the symbolic result of a product is still valid for seemingly nonsensical values of the limits. Then the Karr convention allows us to give a perfectly valid interpretation to those products by interchanging the limits according to the above rules:
>>> P = Product(2, (i, 10, n)).doit() >>> P 2**(n  9) >>> P.subs(n, 5) 1/16 >>> Product(2, (i, 10, 5)).doit() 1/16 >>> 1/Product(2, (i, 6, 9)).doit() 1/16
An explicit example of the Karr summation convention applied to products:
>>> P1 = Product(x, (i, a, b)).doit() >>> P1 x**(a + b + 1) >>> P2 = Product(x, (i, b+1, a1)).doit() >>> P2 x**(a  b  1) >>> simplify(P1 * P2) 1
And another one:
>>> P1 = Product(i, (i, b, a)).doit() >>> P1 RisingFactorial(b, a  b + 1) >>> P2 = Product(i, (i, a+1, b1)).doit() >>> P2 RisingFactorial(a + 1, a + b  1) >>> P1 * P2 RisingFactorial(b, a  b + 1)*RisingFactorial(a + 1, a + b  1) >>> simplify(P1 * P2) 1
See also
Sum
,summation
,product
References
 R73
Michael Karr, “Summation in Finite Terms”, Journal of the ACM, Volume 28 Issue 2, April 1981, Pages 305350 http://dl.acm.org/citation.cfm?doid=322248.322255
 R74
https://en.wikipedia.org/wiki/Multiplication#Capital_Pi_notation
 R75

is_convergent
()[source]¶ See docs of Sum.is_convergent() for explanation of convergence in SymPy.
The infinite product:
\[\prod_{1 \leq i < \infty} f(i)\]is defined by the sequence of partial products:
\[\prod_{i=1}^{n} f(i) = f(1) f(2) \cdots f(n)\]as n increases without bound. The product converges to a nonzero value if and only if the sum:
\[\sum_{1 \leq i < \infty} \log{f(n)}\]converges.
Examples
>>> from sympy import Interval, S, Product, Symbol, cos, pi, exp, oo >>> n = Symbol('n', integer=True) >>> Product(n/(n + 1), (n, 1, oo)).is_convergent() False >>> Product(1/n**2, (n, 1, oo)).is_convergent() False >>> Product(cos(pi/n), (n, 1, oo)).is_convergent() True >>> Product(exp(n**2), (n, 1, oo)).is_convergent() False
References

reverse_order
(*indices)[source]¶ Reverse the order of a limit in a Product.
Usage
reverse_order(expr, *indices)
reverses some limits in the expressionexpr
which can be either aSum
or aProduct
. The selectors in the argumentindices
specify some indices whose limits get reversed. These selectors are either variable names or numerical indices counted starting from the innermost limit tuple.Examples
>>> from sympy import Product, simplify, RisingFactorial, gamma, Sum >>> from sympy.abc import x, y, a, b, c, d >>> P = Product(x, (x, a, b)) >>> Pr = P.reverse_order(x) >>> Pr Product(1/x, (x, b + 1, a  1)) >>> Pr = Pr.doit() >>> Pr 1/RisingFactorial(b + 1, a  b  1) >>> simplify(Pr) gamma(b + 1)/gamma(a) >>> P = P.doit() >>> P RisingFactorial(a, a + b + 1) >>> simplify(P) gamma(b + 1)/gamma(a)
While one should prefer variable names when specifying which limits to reverse, the index counting notation comes in handy in case there are several symbols with the same name.
>>> S = Sum(x*y, (x, a, b), (y, c, d)) >>> S Sum(x*y, (x, a, b), (y, c, d)) >>> S0 = S.reverse_order(0) >>> S0 Sum(x*y, (x, b + 1, a  1), (y, c, d)) >>> S1 = S0.reverse_order(1) >>> S1 Sum(x*y, (x, b + 1, a  1), (y, d + 1, c  1))
Of course we can mix both notations:
>>> Sum(x*y, (x, a, b), (y, 2, 5)).reverse_order(x, 1) Sum(x*y, (x, b + 1, a  1), (y, 6, 1)) >>> Sum(x*y, (x, a, b), (y, 2, 5)).reverse_order(y, x) Sum(x*y, (x, b + 1, a  1), (y, 6, 1))
See also
index
,reorder_limit
,reorder
References
 R77
Michael Karr, “Summation in Finite Terms”, Journal of the ACM, Volume 28 Issue 2, April 1981, Pages 305350 http://dl.acm.org/citation.cfm?doid=322248.322255
Concrete Functions Reference¶

sympy.concrete.summations.
summation
(f, *symbols, **kwargs)[source]¶ Compute the summation of f with respect to symbols.
The notation for symbols is similar to the notation used in Integral. summation(f, (i, a, b)) computes the sum of f with respect to i from a to b, i.e.,
b ____ \ ` summation(f, (i, a, b)) = ) f /___, i = a
If it cannot compute the sum, it returns an unevaluated Sum object. Repeated sums can be computed by introducing additional symbols tuples:
>>> from sympy import summation, oo, symbols, log >>> i, n, m = symbols('i n m', integer=True)
>>> summation(2*i  1, (i, 1, n)) n**2 >>> summation(1/2**i, (i, 0, oo)) 2 >>> summation(1/log(n)**n, (n, 2, oo)) Sum(log(n)**(n), (n, 2, oo)) >>> summation(i, (i, 0, n), (n, 0, m)) m**3/6 + m**2/2 + m/3
>>> from sympy.abc import x >>> from sympy import factorial >>> summation(x**n/factorial(n), (n, 0, oo)) exp(x)
See also
Sum
,Product
,product

sympy.concrete.products.
product
(*args, **kwargs)[source]¶ Compute the product.
The notation for symbols is similar to the notation used in Sum or Integral. product(f, (i, a, b)) computes the product of f with respect to i from a to b, i.e.,
b _____ product(f(n), (i, a, b)) =   f(n)   i = a
If it cannot compute the product, it returns an unevaluated Product object. Repeated products can be computed by introducing additional symbols tuples:
>>> from sympy import product, symbols >>> i, n, m, k = symbols('i n m k', integer=True)
>>> product(i, (i, 1, k)) factorial(k) >>> product(m, (i, 1, k)) m**k >>> product(i, (i, 1, k), (k, 1, n)) Product(factorial(k), (k, 1, n))

sympy.concrete.gosper.
gosper_normal
(f, g, n, polys=True)[source]¶ Compute the Gosper’s normal form of
f
andg
.Given relatively prime univariate polynomials
f
andg
, rewrite their quotient to a normal form defined as follows:\[\frac{f(n)}{g(n)} = Z \cdot \frac{A(n) C(n+1)}{B(n) C(n)}\]where
Z
is an arbitrary constant andA
,B
,C
are monic polynomials inn
with the following properties:\(\gcd(A(n), B(n+h)) = 1 \forall h \in \mathbb{N}\)
\(\gcd(B(n), C(n+1)) = 1\)
\(\gcd(A(n), C(n)) = 1\)
This normal form, or rational factorization in other words, is a crucial step in Gosper’s algorithm and in solving of difference equations. It can be also used to decide if two hypergeometric terms are similar or not.
This procedure will return a tuple containing elements of this factorization in the form
(Z*A, B, C)
.Examples
>>> from sympy.concrete.gosper import gosper_normal >>> from sympy.abc import n
>>> gosper_normal(4*n+5, 2*(4*n+1)*(2*n+3), n, polys=False) (1/4, n + 3/2, n + 1/4)

sympy.concrete.gosper.
gosper_term
(f, n)[source]¶ Compute Gosper’s hypergeometric term for
f
.Suppose
f
is a hypergeometric term such that:\[s_n = \sum_{k=0}^{n1} f_k\]and \(f_k\) doesn’t depend on \(n\). Returns a hypergeometric term \(g_n\) such that \(g_{n+1}  g_n = f_n\).
Examples
>>> from sympy.concrete.gosper import gosper_term >>> from sympy.functions import factorial >>> from sympy.abc import n
>>> gosper_term((4*n + 1)*factorial(n)/factorial(2*n + 1), n) (n  1/2)/(n + 1/4)

sympy.concrete.gosper.
gosper_sum
(f, k)[source]¶ Gosper’s hypergeometric summation algorithm.
Given a hypergeometric term
f
such that:\[s_n = \sum_{k=0}^{n1} f_k\]and \(f(n)\) doesn’t depend on \(n\), returns \(g_{n}  g(0)\) where \(g_{n+1}  g_n = f_n\), or
None
if \(s_n\) can not be expressed in closed form as a sum of hypergeometric terms.Examples
>>> from sympy.concrete.gosper import gosper_sum >>> from sympy.functions import factorial >>> from sympy.abc import i, n, k
>>> f = (4*k + 1)*factorial(k)/factorial(2*k + 1) >>> gosper_sum(f, (k, 0, n)) (factorial(n) + 2*factorial(2*n + 1))/factorial(2*n + 1) >>> _.subs(n, 2) == sum(f.subs(k, i) for i in [0, 1, 2]) True >>> gosper_sum(f, (k, 3, n)) (60*factorial(n) + factorial(2*n + 1))/(60*factorial(2*n + 1)) >>> _.subs(n, 5) == sum(f.subs(k, i) for i in [3, 4, 5]) True
References
 R78
Marko Petkovsek, Herbert S. Wilf, Doron Zeilberger, A = B, AK Peters, Ltd., Wellesley, MA, USA, 1997, pp. 73–100