SymPy 离散具体数学（Concrete）：超几何项判定、Gosper 求和与序列猜想工具实战指南

SymPy 离散具体数学（Concrete）：超几何项判定、Gosper 求和与序列猜想工具实战指南

【免费下载链接】sympyA computer algebra system written in pure Python项目地址: https://gitcode.com/GitHub_Trending/sy/sympy

SymPy 的concrete子模块（源码位于 sympy/concrete）是纯 Python 计算机代数系统中面向“具体数学”（Concrete Mathematics）的离散符号计算工具集：它以超几何项（hypergeometric term）为理论核心，向上提供Sum/Product求和求积、Gosper 超几何求和算法，向下提供序列识别（guess）与生成函数猜想工具。读完本文，你将掌握如何判定一个序列是否为超几何项、如何用hypersimp得到相邻项之比的最简多项式商、如何用gosper_sum求封闭形式和、以及如何从一串有理数序列反推其递推关系与生成函数。

本文内容以 doc/src/modules/concrete.rst 为骨架，结合仓库内源码与测试（如 test_gosper.py、test_guess.py）展开，所有示例均可直接在 Python 中复制运行。

一、超几何项：离散求和与递推的中心舞台

在递推求解（recurrence solving）与求和（summation）中，超几何项（hypergeometric term）占据核心地位。形式化地说，超几何项是被一阶线性递推算子零化的序列：给定序列a(n)，若其相邻项之比a(n+1)/a(n)是n的有理函数，则称a(n)为超几何项。

直观理解：多项式、阶乘、组合数、上升/下降阶乘、伽马函数、指数函数等"具体数学"中的常见物种，其相邻项比值都能化为n的多项式商，因此都是超几何项。

1.1 用is_hypergeometric快速判定

SymPy 在Basic基类上提供了is_hypergeometric(k)方法（实现见 sympy/core/basic.py），它内部调用hypersimp(self, k)并判断结果是否为None：

from sympy import * n, k = symbols('n,k') # 多项式当然是超几何项 (n**2 + 1).is_hypergeometric(n) # True # 具体数学中的常见物种 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

需要注意，二项式系数以及上升、下降阶乘对两个参数都是超几何的（在另一个参数上同样成立）：

binomial(n, k).is_hypergeometric(k) # True rf(n, k).is_hypergeometric(k) # True ff(n, k).is_hypergeometric(k) # True

1.2 整数线性参数依然成立

上述所有示例对n的整数线性参数依然成立，这是求和算法能够处理形如factorial(2n)、binomial(3n+1, k)这类项的关键：

factorial(2*n).is_hypergeometric(n) # True binomial(3*n+1, k).is_hypergeometric(n) # True rf(n+1, k-1).is_hypergeometric(n) # True ff(n-1, k+1).is_hypergeometric(n) # True gamma(5*n).is_hypergeometric(n) # True (2**(n-7)).is_hypergeometric(n) # True

1.3 非线性参数会使序列失去超几何性

一旦参数变为非线性（如n**2、n**3），相邻项之比不再是有理函数，判定结果即变为False：

factorial(n**2).is_hypergeometric(n) # False (2**(n**3 + 1)).is_hypergeometric(n) # False

实现细节：从源码看，is_hypergeometric对Piecewise表达式直接返回None（见 basic.py），即分段定义序列不在超几何判定范围之内。

二、hypersimp：把相邻项之比化简为最小次数多项式商

如果不仅想知道"是否是超几何项"，还想得到相邻项之比化简后的最简形式，应使用hypersimp(f, k)函数（定义于 sympy/simplify/simplify.py）。

其工作流程在源码 docstring 中清晰给出，共三步：

1. 尽可能把所有函数改写为伽马函数（gamma）形式；
2. 把所有 gamma 改写为 gamma 与整数（或绝对常数）指数的上升阶乘的乘积；
3. 化简嵌套分式与幂；若结果恰为多项式商，则约分降低其总次数。

若f(k)是超几何项，函数返回最小次数的多项式商f(k+1)/f(k)；否则返回None表示该序列不是超几何项：

>>> from sympy import hypersimp, factorial >>> from sympy.abc import n >>> hypersimp(factorial(2*n), n) 2*(n + 1)*(2*n + 1) >>> hypersimp(factorial(n**2), n) # 返回 None（空行）

第一行结果2*(n+1)*(2*n+1)正是factorial(2n+2)/factorial(2n)的约分结果；第二行因为相邻项之比不是有理函数，返回None。

源码补充：hypersimp先用g = f.subs(k, k+1) / f构造相邻项之比，再经rewrite(gamma)、expand_func、powsimp等化简，最后用is_rational_function(k)判断是否为有理函数，是则返回simplify约分结果，否则返回None。同一文件中的hypersimilar(f, g, k)还提供"超相似"判定——两个项的商为k的有理函数时返回True，在求解递推关系时很有用。

三、Gosper 算法：超几何求和的封闭形式

Gosper 算法是超几何求和的经典算法，用于计算形如s_n = sum(f(k), (k, 0, n-1))的和式封闭形式，其中f为不依赖于n的超几何项。其核心思想是寻找另一个超几何项g_n满足g_{n+1} - g_n = f_n（即f的"不定和分"），从而将求和问题转化为简单的代值计算。算法参考了 Marko Petkovsek、Herbert S. Wilf 与 Doron Zeilberger 的著作《A = B》（AK Peters, 1997, pp. 73–100）。

sympy.concrete.gosper模块提供三个逐层递进的函数（sympy/concrete/gosper.py）：

3.1gosper_normal：Gosper 正规形

gosper_normal(f, g, n, polys=True)把互素单变量多项式f(n)/g(n)改写为如下正规形：

f(n)/g(n) = Z · (A(n)·C(n+1)) / (B(n)·C(n))

其中Z为任意常数，A、B、C是n的首一多项式，且满足三条互素性条件（gcd(A(n), B(n+h)) = 1对所有自然数h成立、gcd(B(n), C(n+1)) = 1、gcd(A(n), C(n)) = 1）。这种"有理分解"是 Gosper 算法与差分方程求解的关键步骤，也可用于判定两个超几何项是否相似。

>>> 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)

返回三元组(Z*A, B, C)。测试 test_gosper.py 同时验证了polys=True（返回Poly对象）与polys=False（返回普通表达式）两种模式结果一致。

3.2gosper_term：寻找不定和分

gosper_term(f, n)对给定的超几何项f，返回满足g_{n+1} - g_n = f_n的超几何项g_n。其内部流程为：先用hypersimp求相邻项之比，再用gosper_normal做有理分解，之后求解一个关于待定系数的线性方程组（源码中通过构造H = A*x.shift(1) - B*x - C并solve系数实现）。若f不是超几何项、或不可 Gosper 求和，则返回None。

>>> from sympy.concrete.gosper import gosper_term >>> from sympy 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)

3.3gosper_sum：直接得到封闭和

gosper_sum(f, k)是对用户最友好的入口，给定超几何项f，计算g_n - g(0)（其中g_{n+1} - g_n = f_n）；若和式无法表示为超几何项的封闭形式，返回None。它接受两种调用形式：定和(k, a, b)与不定和（仅传符号k）。

>>> from sympy.concrete.gosper import gosper_sum >>> from sympy import factorial >>> from sympy.abc import 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

注意到 docstring 与测试都用subs(n, 2) == sum(...)做数值回验，这是验证封闭形式正确性的标准做法。测试文件 test_gosper.py 还覆盖了更多经典结果：

>>> from sympy.concrete.gosper import gosper_sum >>> from sympy import factorial, binomial >>> from sympy.abc import k, n gosper_sum(1, (k, 0, n)) # n + 1 gosper_sum(k, (k, 0, n)) # n*(n + 1)/2 gosper_sum(k**2, (k, 0, n)) # n*(n + 1)*(2*n + 1)/6 gosper_sum(k**3, (k, 0, n)) # n**2*(n + 1)**2/4 gosper_sum(2**k, (k, 0, n)) # 2*2**n - 1 gosper_sum(factorial(k), (k, 0, n)) # None（不可超几何求和） gosper_sum(binomial(n, k), (k, 0, n)) # None

实践要点：factorial(k)与binomial(n, k)的定和返回None，说明它们不满足 Gosper 可和性（这正是为什么sum(binomial(n,k), k)需要借助其他机制，如二项式定理/超几何恒等式）。遇到None时应转而使用Sum对象或其数值求值能力，而不是强行期望封闭形式。

四、Sum、Product 与求和/求积入口函数

concrete模块的类参考（见 concrete.rst）包含三个核心类：

• sympy.concrete.summations.Sum（含ExprWithIntLimits基类）：表示未求值的符号和式，可通过.doit()求值；
• sympy.concrete.products.Product：表示未求值的符号连乘，可通过.doit()求值；
• sympy.concrete.expr_with_intlimits.ExprWithIntLimits：Sum/Product共用的"整数上下限"抽象基类，承载了换元、下限平移等公共逻辑。

对应的函数级入口为summation与product（定义于 sympy/concrete/summations.py 与 sympy/concrete/products.py），两者语法与Integral一致，都是f后跟(i, a, b)这样的元组，且支持多重求和/求积（重复传入多个符号元组）：

>>> from sympy import summation, product, symbols, oo, log >>> i, n, m = symbols('i n m', integer=True) # 求和：计算失败时返回未求值的 Sum 对象 >>> 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 import factorial >>> from sympy.abc import x >>> summation(x**n/factorial(n), (n, 0, oo)) exp(x) # 求积：与 Sum 对称 >>> i, k, m = symbols('i k m', integer=True) >>> product(i, (i, 1, k)) factorial(k) >>> product(m, (i, 1, k)) m**k

从源码看，summation的实现就是return Sum(f, *symbols, **kwargs).doit(deep=False)，product类似地构造Product并求值（若无法求值则返回未求值对象），因此掌握Sum/Product的doit机制即可理解这两个入口函数。Sum内部还集成了多种求值策略（telescopic裂项、eval_sum多项式求和等，见 summations.py 中的telescopic_direct/telescopic辅助函数），Gosper 算法也是其超几何项求和的候选策略之一。

五、序列猜想工具：从数列反推公式与生成函数

sympy.concrete.guess模块（sympy/concrete/guess.py）提供一组"从若干项猜公式"的工具，核心函数guess改编自 Christian Krattenthaler 的 Mathematica 软件包Rate.m。整套工具在测试文件 test_guess.py 中有完整验证。

5.1find_simple_recurrence：识别线性递推

find_simple_recurrence(v, A=Function('a'), N=Symbol('n'))从若干个整数（或有理数）项中检测并返回递推关系。返回表达式中函数名默认为a，主变量默认为n，且最小下标恒为n（不会是n-1、n-2等）。其底层函数find_simple_recurrence_vector(l)返回长度为n的系数向量（当发现n阶递推时）；若只返回[0]，则说明未找到关系——注意该函数对二次无理数等特殊实数需谨慎使用（源码 docstring 有明确警告）。

>>> from sympy.concrete.guess import find_simple_recurrence >>> from sympy import fibonacci >>> find_simple_recurrence([fibonacci(k) for k in range(12)]) -a(n) - a(n + 1) + a(n + 2) # 自定义函数名与主变量 >>> from sympy import Function, Symbol >>> a = [1, 1, 1] >>> for k in range(15): a.append(5*a[-1]-3*a[-2]+8*a[-3]) >>> find_simple_recurrence(a, A=Function('f'), N=Symbol('i')) -8*f(i) + 3*f(i + 1) - 5*f(i + 2) + f(i + 3) >>> from sympy.concrete.guess import find_simple_recurrence_vector >>> find_simple_recurrence_vector([fibonacci(k) for k in range(12)]) [1, -1, -1]

5.2rationalize：从浮点数识别有理数

rationalize(x, maxcoeff=10000)通过连分数从浮点值（或mpmath.mpf）识别有理数。算法在检测到超过阈值（默认 10000）的大部分商（partial quotient）时停止。与Fraction.from_decimal、mpmath.identify、nsimplify等方法不同，它关注的是部分商的量级而非全局近似精度——如果该实数"已知是有理数"，即使分母很大也能在默认参数下正确识别。

>>> from sympy.concrete.guess import rationalize >>> from mpmath import cos, pi >>> rationalize(cos(pi/3)) 1/2 >>> from mpmath import mpf >>> rationalize(mpf("0.333333333333333")) 1/3 # 提高 maxcoeff 阈值可用作近似 >>> rationalize(pi, maxcoeff=250) 355/113

5.3guess_generating_function：猜生成函数

guess_generating_function(v, X=Symbol('x'), types=['all'], maxsqrtn=2)尝试为有理数序列v"猜"出生成函数，返回一个字典，键为生成函数类型名。目前实现了六种类型：

参数说明与使用要点：

• types默认为['all']；只关心部分类型时可传入类型列表以节省时间。注意：丢弃某类型只是不为其做额外计算，结果字典中仍可能包含该类型（因为可以从其他类型轻松转换而来）；
• lgdogf与lgdegf在序列首项为 0 时不会被计算，此时可先去掉前导零再重试；
• maxsqrtn（默认 2）指定要测试的"有理函数的 n 次方根"的最大阶数，用于识别生成函数为某有理函数的平方根等情形。

>>> from sympy.concrete.guess import guess_generating_function as ggf >>> ggf([k+1 for k in range(12)], types=['ogf', 'lgf', 'hlgf']) {'hlgf': 1/(1 - x), 'lgf': 1/(x + 1), 'ogf': 1/(x**2 - 2*x + 1)} >>> from sympy import sympify >>> l = sympify("[3/2, 11/2, 0, -121/2, -363/2, 121]") >>> ggf(l) {'ogf': (x + 3/2)/(11*x**2 - 3*x + 1)} >>> from sympy import fibonacci >>> ggf([fibonacci(k) for k in range(5, 15)], types=['ogf']) {'ogf': (3*x + 5)/(-x**2 - x + 1)} >>> from sympy import factorial >>> ggf([factorial(k) for k in range(12)], types=['ogf', 'egf', 'lgf']) {'egf': 1/(1 - x)} >>> ggf([k+1 for k in range(12)], types=['egf']) {'egf': (x + 1)*exp(x), 'lgdegf': (x + 2)/(x + 1)} # n 次方根检测（对应 OEIS A108626 序列） >>> ggf([1, 2, 5, 14, 41, 124, 383, 1200, 3799, 12122, 38919])['ogf'] sqrt(1/(x**4 + 2*x**2 - 4*x + 1))

其中guess_generating_function_rational(v, X=Symbol('x'))是只处理"有理生成函数"的低层版本：先用find_simple_recurrence_vector求分母q，再按卷积关系求分子p。其返回(3*x + 5)/(-x**2 - x + 1)之类的分式，当未找到时返回None。源码 docstring 同时提示它与 sympy/series/approximants.py 中的approximants（Padé 近似型序列逼近）功能相关，可互为补充。

5.4guess：从序列猜组合公式

guess(l, all=False, evaluate=True, niter=2, variables=None)从一串有理数序列猜出闭式公式，返回一个公式列表（可能是多个等价结果）。参数语义：

• all=False（默认）：一旦某次迭代出结果即停止计算，加速流程；设为True则继续更多迭代，可能返回更多（可能与前序等价）的公式；
• evaluate=True（默认）：对结果中的连乘进行求值；设为False则保留未求值的Product对象（便于观察结构）；
• niter=2：迭代次数，最大可取len(l)-1。迭代阶数越高结果越复杂：
  • 第一次迭代返回多项式或有理函数；
  • 第二次迭代返回上升阶乘及其倒数的乘积；
  • 第三次迭代返回"上升阶乘乘积的乘积"；
  • 依此类推。
• variables=None：返回公式默认包含符号i0, i1, i2, ...，其中主变量是i0（辅助变量为i1, i2, ...）；也可传入自定义符号列表（长度应不小于niter，主变量取列表第一个符号）。

>>> from sympy.concrete.guess import guess >>> guess([1,2,6,24,120], evaluate=False) [Product(i1 + 1, (i1, 1, i0 - 1))] >>> from sympy import symbols >>> r = guess([1,2,7,42,429,7436,218348,10850216], niter=4) >>> i0 = symbols("i0") >>> [r[0].subs(i0,n).doit() for n in range(1,10)] [1, 2, 7, 42, 429, 7436, 218348, 10850216, 911835460]

注意事项（来自源码 docstring 与实现）：若序列除最后一项外含有 0，guess直接返回空列表[]（源码第 448 行检查any(a==0 for a in l[:-1])）；内部通过rational_interpolate（多项式有理插值，见 sympy/polys/polyfuncs.py）与相邻项比值变换逐层推进。

六、完整工作流：从"猜"到"证"再到"求"

将上述工具串联起来，可形成一条经典的"具体数学"研究闭环——例如从斐波那契数列出发：

from sympy import fibonacci, simplify, summation from sympy.concrete.guess import guess_generating_function, find_simple_recurrence fib = [fibonacci(k) for k in range(15)] # 第 1 步：识别递推关系 print(find_simple_recurrence(fib)) # -a(n) - a(n+1) + a(n+2) # 第 2 步：猜生成函数并解析验证 print(guess_generating_function(fib, types=['ogf'])) # 第 3 步：对超几何项使用 Gosper 求和得到封闭形式 from sympy.concrete.gosper import gosper_sum from sympy.abc import k, n print(gosper_sum(k**3, (k, 0, n))) # n**2*(n + 1)**2/4

先由find_simple_recurrence得到递推结构，再由guess_generating_function获得生成函数解析式，最后对超几何项用gosper_sum/summation求封闭和——这正是concrete模块设计上"以超几何项为中心，统一支撑递推、求和与序列识别"的体现。

七、常用资源与扩展阅读

• 模块文档：doc/src/modules/concrete.rst
• 求和与求积实现：sympy/concrete/summations.py、sympy/concrete/products.py、sympy/concrete/expr_with_intlimits.py
• Gosper 算法实现：sympy/concrete/gosper.py，配套测试 sympy/concrete/tests/test_gosper.py
• 序列猜想实现：sympy/concrete/guess.py，配套测试 sympy/concrete/tests/test_guess.py
• 超几何判定与化简：is_hypergeometric（sympy/core/basic.py）、hypersimp（sympy/simplify/simplify.py）
• 相关算法参考文献：W. Koepf《Algorithms for m-fold Hypergeometric Summation》（J. Symbolic Computation, 1995）为hypersimp的算法依据；Graham、Knuth、Patashnik《Concrete Mathematics》与 OEIS 生成函数词条为guess_generating_function的参考来源。

总结：concrete模块以"超几何项"为统一视角，将离散求和、递推求解与序列识别整合为一套可操作的符号计算工具链。掌握is_hypergeometric/hypersimp的判定与化简、gosper_sum的封闭求和、summation/product的符号求值，以及guess系列工具的逆向猜想能力，即可在日常研究与工程实践中完成从数列观察到封闭公式的完整闭环。

【免费下载链接】sympyA computer algebra system written in pure Python项目地址: https://gitcode.com/GitHub_Trending/sy/sympy