深入 SymPy 级数展开模块Gruntz 极限算法、series 展开、Order 阶项与留数计算【免费下载链接】sympyA computer algebra system written in pure Python项目地址: https://gitcode.com/GitHub_Trending/sy/sympySymPy 的sympy.series模块是纯 Python 计算机代数系统sympy/中负责极限计算与级数展开的核心模块。本指南以官方文档 doc/src/modules/series/series.rst 为主线完整讲解该模块的六大功能limit()/Limit极限计算、Gruntz 算法及其理论、series()直观级数展开、Order大 O 阶项自动追踪、Richardson/Shanks 级数加速以及residue()留数计算。读完本文你既能熟练使用这些 API 求解具体问题也能理解其底层的增长阶比较MRV 集与渐近展开原理并知道如何在源码 sympy/series/ 中定位实现与调试。模块总览以极限计算为核心正如文档开头所述本模块的主要目的是极限的计算。模块由sympy.series包内多个文件分工实现sympy/series/文件职责limits.py用户入口limit()与未求值极限类Limitgruntz.pyGruntz 算法完整实现gruntz、mrv、rewrite等series.pyseries()函数式级数展开包装order.pyOrder大 O 阶项类acceleration.pyrichardson()、shanks()级数加速residues.pyresidue()留数计算这一架构体现了直观接口 底层严格算法的设计日常求极限走limit()它会先用启发式方法快速处理简单情形复杂情形才落到gruntz()。极限计算limit()与Limit类函数入口limit(e, z, z0, dir)文档给出的核心接口是 limits.py 中的limit()from sympy import limit, sin, oo from sympy.abc import x limit(sin(x)/x, x, 0) # 1 limit(1/x, x, 0) # oo默认 dir从右侧 limit(1/x, x, 0, dir-) # -oo从左侧 limit(1/x, x, 0, dir-) # zoo双侧极限发散为复无穷 limit(1/x, x, oo) # 0参数说明与 limit 源码 docstring 一致e要求极限的表达式z极限变量其他符号视为常量多元极限不支持z0趋于的点可以是任意表达式包括oo与-oodir方向表示右极限z→z0-表示左极限z→z0--表示双侧。当z0为无穷时dir由无穷的方向决定即对oo取dir-。limit()的实现limits.py只是Limit(e, z, z0, dir).doit(deepFalse)的语法糖返回的是Limit对象求值后的结果。未求值极限类LimitLimit表示一个尚未求值的极限表达式可以直接构造from sympy import Limit, sin from sympy.abc import x Limit(sin(x)/x, x, 0) # Limit(sin(x)/x, x, 0, dir) Limit(1/x, x, 0, dir-) # Limit(1/x, x, 0, dir-)从 构造逻辑 可以看出几个重要约束当z0本身含极限变量时抛出NotImplementedError不支持趋于动点dir只能是、-、-三者之一。其doit()方法limits.py内部的分支顺序本身就揭示了求解路径双侧极限分别求左右极限并比较趋于无穷时做代换、规范化符号对Float浮点做nsimplify精确化避免舍入误差导致条件分支误判若表达式在z0处亚纯meromorphic直接用leadterm提取首项判定兜底调用gruntz()见下节若 Gruntz 也失败PoleError/ValueError再退回heuristics()逐项求值。这里的关键事实是文档明确指出的极限计算的真正主力是gruntz()函数它实现了 Gruntz 算法。limit()只是先启发式、后 Gruntz的调度层。相关的启发式逻辑可在 heuristics() 中查看——它逐项计算各参数的极限再合并结果只对简单情形快而有效。Gruntz 算法极限计算背后的理论函数增长速度的序关系Gruntz 算法的第一步是建立函数之间的增长序。设f(x)、g(x)是两个实值函数且当x→∞时都趋于∞。文档定义了支配关系若对任意a, b ∈ ℝ₀都有lim_{x→∞} f(x)^a / g(x)^b ∞则称f(x)支配g(x)记作f(x) ≻ g(x)。若f(x)与g(x)互不支配则称它们属于同一个可比较类comparability class记作f(x) ≍ g(x)。该定义可通过f ~ 1/f、f ~ -f自然地推广到趋于 0 或±∞的函数。源码 gruntz.py 的模块 docstring 给出了等价的工程化判据计算L lim log|f(x)| / log|g(x)|x→∞L±∞则f ≻ gL0则f ≺ gL为有限非零值则同属一类。文档给出的一组经典例子在源码 docstring 中也重复出现e^x ≻ x^me^{x²} ≻ e^{mx}e^{e^x} ≻ e^{x^m}x^m ≍ x^n任何正幂同属一类e^{x1/x} ≍ e^{xlog x} ≍ e^x也就是说2 x exp(x) exp(x**2) exp(exp(x))而2 ~ 3 ~ -5、x ~ x**2 ~ x**3 ~ 1/x、exp(x) ~ exp(-x) ~ exp(2x)。最低的可比较类被设定为2~3~-5即常数类log(x)低于该类不参与比较。支撑极限判定的定理文档紧接着给出两条核心结论定理一主项支配求和设ω, g₁, g₂, …是x的函数lim ω 0且ω ≻ gᵢ对所有i成立设c₁ c₂ …为严格递增的实数。则lim_{x→∞} Σᵢ gᵢ ω^{cᵢ} lim_{x→∞} g₁ ω^{c₁}即求和式的极限由指数最小的首项决定其余项都更快地趋于零。推论对g·ω^c的三种情形c 0时lim g·ω^c 0c 0时lim g·ω^c ±∞符号由g的最终符号决定c 0时lim g·ω^0 lim g。这两条正是 limitinf() 的实现依据算出展开式的首项(c0, e0)后依据e0的符号直接返回0、±∞或递归地求lim c0。四步求解策略据此文档给出了计算lim_{x→∞} f(x)的标准流程找出f(x)的 MRV 集most rapidly varying subexpressions最快变化子表达式集从f的所有子表达式中选出在≻关系下极大的那些元素选取ω选择一个与 MRV 集元素同属一个可比较类、且lim ω 0的函数例如 MRV 集为{eˣ, e²ˣ}时可取ω e⁻ˣ按ω展开级数把f写成f c₀ω^{e₀} c₁ω^{e₁} … O(ω^{eₙ})e₀e₁…eₙc₀≠0使之满足定理一的前提应用定理并递归由首项判定极限必要时对g₁(x)递归执行以上步骤。源码 gruntz.py 的注释精确对应这四步求 MRV 集用mrv()选ω并改写用rewrite()求首项用mrv_leadterm()/leadterm()最终收敛在limitinf()中完成。实现中的关键函数文档 Reference 节逐一列出了 gruntz.py 的公开接口对应实现如下gruntz(e, z, z0, dir)算法入口。z0可为任意表达式含oo/-oodir求右极限、dir-求左极限z0无穷时dir无意义。内部先把一切极限通过代换z0 1/z、z0 - 1/z、-z归约为z→∞再交给limitinf()最后用rewrite(intractable)把结果美化回熟悉的形式。注意gruntz()不接收非符号作为第二个参数第 679 行。compare(a, b, x)计算L lim log(a)/log(b)返回、或。对exp类表达式会直接取指数部分log(exp(…))必须在此化简以保证算法终止。mrv(e, x)返回一个SubsSet包含e的所有 MRV 子表达式及其重写信息。对Add/Mul递归拆分比较mrv_max1/mrv_max3对log、exp、一般函数、导数分别处理未定义函数会抛ValueError多变量函数如BesselJ(x, x)与导数暂未实现NotImplementedError。mrv_leadterm(e, x)返回e的首项(c₀, e₀)。内部调用rewrite()后对w求leadterm失败时逐步增大_eval_nseries的阶数重试还用ilcm对分母取最小公倍数消除非整数幂展开的困难。limitinf(e, x)计算x→∞的极限。先做rewrite(tractable)、powdenest等规范化再依e₀符号返回0/±∞/递归求lim c₀。sign(e, x)返回e在x→∞时的最终符号1/0/-1对反复变号的函数如sin(x)结果未定义。实现先查假设is_positive等再对Mul、exp、Pow、log结构递归最后兜底走mrv_leadterm。rewrite(e, Omega, x, wsym)算法中最复杂的一步。docstring 明确说当算法失败时bug 通常在级数展开即 SymPy 本身或 rewrite 中。它借助 build_expression_tree() 按依赖高度排序 MRV 表达式把e改写为关于w趋零变量的级数形式并同步计算出log(w)。SubsSetMRV 表达式→哑变量的字典附带rewrites映射。类 docstring 解释了为何不能用普通的subs()它太聪明会错误化简文档中的-1/w exp(exp(p)*exp(-exp(-p))/(1 - 1/p))例子说明普通subs无法正确处理嵌套 MRV 表达式。mrv_max1 / mrv_max3两组 MRV 集合取更大者同可比较类时取并集gruntz.py。文档未展开的四个细节与调试手段文档Notes节诚实指出前面的叙述略过了四个关键问题如何判定f ≻ g、g ≻ f还是f ≍ g→ 由compare()实现如何求 MRV 集→ 由mrv()实现如何计算级数展开→ 依赖 SymPy 的leadterm/_eval_nseries算法为何必然终止→ 由 Gruntz 博士论文证明gruntz.py 的模块 docstring 还提供了调试方法Gruntz 算法高度递归难以在调试器里追踪可通过设置环境变量SYMPY_DEBUGTrue开启漂亮的调试打印例如SYMPY_DEBUGTrue python -m isympy然后执行limit(sin(x)/x, x, 0)控制台会输出完整的递归树limitinf(_x*sin(1/_x), _x) 1 -mrv_leadterm(_x*sin(1/_x), _x) (1, 0) | -mrv(_x*sin(1/_x), _x) set([_x]) | | -mrv(_x, _x) set([_x]) | | -mrv(sin(1/_x), _x) set([_x]) ... -sign(0, _x) 0 -limitinf(1, _x) 1逐行核对哪一步出错再回到对应源码函数定位问题。这得益于sympy.utilities.misc.debug_decorator见 gruntz.py对mrv、sign、limitinf、mrv_leadterm、rewrite等函数的装饰。此外算法实现几乎是 Gruntz 博士论文中 Maple 代码的直接改写gruntz.py 注释论文包含完整描述与大量示例是深入理解算法细节的首选资料。更直观的级数展开series()包装文档指出x*cos(x)这种写法比(x*cos(x)).series(x)更直观因此模块提供了一个围绕Basic.series()的包装函数 series()from sympy import Symbol, cos, series x Symbol(x) series(cos(x), x) # 1 - x**2/2 x**4/24 O(x**6)签名与参数series.pyexpr要展开的表达式x展开变量x0展开点取值范围-oo到oon展开到第n项默认 6dir方向为x→x0−为x→x0−当x0为无穷时方向由无穷本身决定对oo即dir-。例如在x2处对tan(x)展开from sympy import series, tan from sympy.abc import x series(tan(x), x, 2, 6, ) # tan(2) (1 tan(2)**2)*(x - 2) ... O((x - 2)**6, (x, 2)) series(tan(x), x, 2, 3, -) # tan(2) (2 - x)*(-tan(2)**2 - 1) ... O((x - 2)**3, (x, 2))注意noo会抛出TypeError源码 docstring 中即有此示例因为项数必须是整数。实现上 series() 只是sympify(expr)后调用expr.series(x, x0, n, dir)——完整语义请参见 Expr.series 的 docstring。阶项追踪Order类级数展开不可避免要记录截断误差本模块通过 Order别名O自动完成阶项的追踪。文档示例from sympy import Symbol, Order x Symbol(x) Order(x) x**2 # O(x) Order(x) 1 # 1 O(x)Order用大 O 记号表征函数的极限行为g(x) O(f(x))x→a当且仅当存在δ0与M0使得|x-a|δ时|g(x)| ≤ M|f(x)|等价于limsup |g(x)/f(x)| ∞。以sin(x)在 0 处展开为例sin(x) x - x³/3! O(x⁵)此时O(x⁵) x⁵/5! - x⁷/7! …且lim |O(x⁵)/x⁵| 1/5!order.py。Order 的构造逻辑 还包含若干自动规范化规则不传变量时默认使用表达式全部自由符号展开点默认为 0O(expr*f(x), x)规约为O(f(x), x)O(expr, x)规约为O(1)O(0, x)规约为0O(f(x), x)会自动取f(x).as_leading_term(x)的首项支持多元O(f(x,y), x, y)假定各变量取极限可交换支持在oo、-oo等展开点内部做1/Dummy()代换归化到 0提供contains()方法判定阶项包含关系如O(x) in O(1, x)为 Trueorder.pyO(1) in O(1, x)→ TrueO(1, x) in O(1)→ FalseO(x) in O(1, x)→ TrueO(x**2) in O(x)→ True算术规则O(x)*x→O(x**2)O(x) - O(x)→O(x)O(cos(x))→O(1)O(cos(x), (x, pi/2))→O(x - pi/2, (x, pi/2))order.py 示例。Order的removeO()返回 0、getO()返回自身是提取/丢弃阶项的常用工具_eval_nseries返回自身保证它不参与再次展开。级数加速Richardson 外推与 Shanks 变换文档中 Series Acceleration 一节标注为 TODO但 acceleration.py 中的两个函数实际上已完整实现参考书目为 Bender Orszag《Advanced Mathematical Methods for Scientists and Engineers》Shanks 变换见 pp. 368-375Richardson 外推见 pp. 375-377。richardson(A, k, n, N)用项A(n), A(n1), …, A(nN1)做 Richardson 外推近似lim_{k→∞} A(k)。经验上取N ≈ 2n效果较好。经典例子是用极限定义算ee (11/n)**n收敛极慢n100时只有两位有效数字2.7048138294而用richardson(e, n, 10, 20)一次即可得到2.7182818285。对 ζ(2) 级数1/k²的部分和取前 100 项仅得1.6349839002Richardson 外推则给出1.6449340668与精确值π²/6一致。shanks(A, k, n, m1)n 项 Shanks 变换S(A)(n)m1时做 m 重递归变换S(S(…S(A)…))(n)。它对在极点/奇点附近收敛慢的 Taylor 级数特别有效例如交错调和级数求log(2)前 100 项只有0.6881721793shanks(A, n, 25)得到0.6931396564shanks(A, n, 25, 5)达到0.6931471806与log(2)的精确值0.6931471805599453…高度吻合。两个函数的实现均基于文献中的递推/插值公式Richardson 使用(nj)^N与二项式系数的加权和acceleration.pyShanks 使用三对角(z*x - y²)/(z x - 2y)递推表acceleration.py。这两个函数的正确性在 test_series.py 中有断言覆盖。留数计算residue()文档 Residues 一节同样标注 TODO但 residues.py 中的 residue(expr, x, x0) 已实现。留数定义为表达式在xx0处幂级数展开中1/(x-x0)项的系数from sympy import Symbol, residue, sin x Symbol(x) residue(1/x, x, 0) # 1 residue(1/x**2, x, 0) # 0 residue(2/sin(x), x, 0) # 2实现原理源码注释当前实现基于级数展开先平移x→xx0用nseries在n ∈ (0,1,2,4,8,16,32)中逐步取到足够的展开阶collect后逐项扫描凡是1/x项就把系数累加进结果若出现无法识别的项则抛NotImplementedError。对于一般函数可参考 Bronstein《Symbolic Integration I》5.6 节纯有理函数则有基于结式resultant的简单算法见 residues.py 的注释。test_residues.py 提供了丰富的验证用例覆盖了各种情形residue(1/(x**2 1), x, I) # -I/2 residue(1/(x**2 1), x, -I) # I/2 residue(1/sin(x)**5, x, 0) # Rational(3, 8) residue(1/(x**4 1), x, exp(I*pi/4)).equals(-(Rational(1, 4) I/4)/sqrt(2))以及展开点不是极点则留数为 0的性质residue(1/x, x, 1) 0。留数定理residues.py 引用使其成为复分析积分与级数求和的核心工具。从文档到实战使用建议综合文档与源码可以总结出几条实用建议日常求极限直接使用limit()。它先走启发式快速路径失败才进入 Gruntz 算法绝大多数情形开箱即用。需要研究算法行为时把dir方向与Limit类配合使用Limit让你保留未求值形式的极限便于符号操作.doit()触发求值dir-会校验左右极限是否一致不一致时抛出带左右值的ValueError。级数展开优先用series()函数式写法并善用Order做截断管理需要丢弃阶项时用.removeO()。gruntz()直接调用的适用前提第二个参数必须是Symbolgruntz.py文档也提示它只在更快的limit()失败时才常被用到。当 Gruntz 结果异常时用SYMPY_DEBUGTrue开启递归打印从输出树中定位出错的mrv/rewrite环节再回到 gruntz.py 对应函数排查。数值上收敛很慢的极限/级数可尝试richardson/shanks加速符号留数计算则用residue()注意其实现基于nseries展开复杂奇点如本性奇点可能需要更高阶展开。如需验证行为或深入学习可查阅 sympy/series/tests/ 下的测试文件Gruntz 算法与极限的断言集中在 test_gruntz.py 与 test_limits.py级数展开见 test_series.py 与 test_nseries.py阶项语义见 test_order.py留数见 test_residues.py。【免费下载链接】sympyA computer algebra system written in pure Python项目地址: https://gitcode.com/GitHub_Trending/sy/sympy创作声明：本文部分内容由AI辅助生成（AIGC），仅供参考