SymPy 中全纯函数Holonomic Functions的表示从微分算子到零化子Annihilator【免费下载链接】sympyA computer algebra system written in pure Python项目地址: https://gitcode.com/GitHub_Trending/sy/sympy全纯函数是 SymPyholonomic模块的核心研究对象——这类函数是带多项式系数的线性齐次常微分方程的解涵盖exp、sin、cos、log、Bessel 函数乃至超几何函数与 Meijer G 函数等大量特殊函数。本指南基于仓库文档 represent.rst讲解如何用「微分算子Differential Operator」构造零化子annihilator并用HolonomicFunction类以「零化子 初值条件」的方式精确定义一个全纯函数读完你将能手动搭建 Weyl 代数、书写任意全纯函数的零化子并理解其背后的源码实现。为什么需要「表示」从微分方程到零化子在 about.rst 中定义了全纯函数的本质存在多项式p_0, p_1, ..., p_r ∈ K[x]K为特征 0 的域如QQ、RR使得p_0·f(x) p_1·f(x) ... p_r·f^(r)(x) 0这个方程可紧凑地写成算子形式L·f 0其中L p_0 p_1·D p_2·D² ... p_r·D^r这里D是微分算子作用于函数即求导满足D^n·y(x) y^(n)(x)而L就是该函数的零化子annihilator。由于常微分方程的通解是一族函数而非单个函数要唯一确定一个全纯函数还必须搭配一组初值条件initial conditions——这正是represent.rst所阐述的表示方案一个全纯函数 一个零化子 一组初值条件例如exp(x)对应L D - 1, f(0) 1sin(x)对应L D² 1, f(0) 0, f(0) 1。用sin(x)走通完整表示流程原文档以sin(x)为例给出了最典型的构造代码。sin(x)满足的微分方程是y(x) y(x) 0通解为C₁·sin(x) C₂·cos(x)通过初值y(0) 0, y(0) 1才唯一锁定为sin(x)。方程改写为(D² 1)·y(x) 0因此零化子就是D² 1。在 SymPy 中表示如下 from sympy.holonomic import DifferentialOperators, HolonomicFunction from sympy.abc import x from sympy import ZZ R, D DifferentialOperators(ZZ.old_poly_ring(x), D) HolonomicFunction(D**2 1, x, 0, [0, 1]) HolonomicFunction((1) (1)*D**2, x, 0, [0, 1])关键点解读DifferentialOperators(ZZ.old_poly_ring(x), D)返回两个对象R是微分算子代数Weyl 代数D是其中的求导算子。多项式系数属于环ZZ[x]整数系数多项式环。HolonomicFunction(D**2 1, x, 0, [0, 1])的四个参数依次是零化子D²1、自变量x、初值点x0 0、初值向量[y(0), y(0)] [0, 1]。输出HolonomicFunction((1) (1)*D**2, x, 0, [0, 1])按D的升幂打印常系数项(1)加上D²项中间的D项系数为 0 被省略。源码层面DifferentialOperators定义于 holonomic.py它内部构造一个DifferentialOperatorAlgebra并返回(ring, ring.derivative_operator)。其中derivative_operator是用多项式列表[base.zero, base.one]构造的即0 1·D对应D算子本身。三个核心类Algebra、Operator、Function原文档通过autoclass指令挂载了HolonomicFunction、DifferentialOperator、DifferentialOperatorAlgebra与函数DifferentialOperators它们的实现集中在 holonomic.py1.DifferentialOperatorAlgebra——Weyl 代数的载体这是微分算子的父环。它本质上是 Ore 代数以Dx为中间变量、以基础多项式环A中的元素为系数的非交换多项式集合遵循交换规则Dx·a σ(a)·Dx δ(a) a ∈ A当σ取恒等映射、δ取标准求导时该代数退化为微分算子代数Weyl 代数即Dx·a a·Dx a。构造时见 holonomic.py生成元generator可以是字符串或非交换Symbol如Dx、D若传None默认使用Symbol(Dx, commutativeFalse)——注意生成元被显式标记为非交换这是 Weyl 代数正确性的关键。base决定系数环。文档指出当前实现为优先级机制使用 SymPy 中较老的 ring 实现old_poly_ring因此示例中采用ZZ.old_poly_ring(x)或QQ.old_poly_ring(x)。2.DifferentialOperator——可像普通表达式一样运算的算子DifferentialOperator由一个多项式列表每个Dx幂次的系数和父代数实例定义列表长度减一即为算子阶数order见 holonomic.py。它的强大之处在于完整实现了运算符重载让零化子可以像 SymPy 表达式一样直接书写 from sympy.holonomic.holonomic import DifferentialOperator, DifferentialOperators from sympy import ZZ from sympy import symbols x symbols(x) R, Dx DifferentialOperators(ZZ.old_poly_ring(x), Dx) # 直接按列表构造1·Dx x²·Dx² DifferentialOperator([0, 1, x**2], R) (1)*Dx (x**2)*Dx**2 # 非交换运算Dx*x 会展开为 x·Dx 1 Dx*x (1) (x)*Dx底层支撑holonomic.py__mul__按交换规则Dx·a a·Dx a实现乘法_mul_Dxi_b迭代计算Dx^i·b并累加__add__、__sub__、__neg__、__truediv__、__pow__均已实现D**2这类写法会走__pow__的快速幂逻辑对纯D幂直接生成[0, 0, 1]这样的稀疏列表__eq__比较多项式列表与父代数is_singular(x0)通过求解最高阶系数多项式的根判断方程在x0处是否奇异holonomic.py。3.HolonomicFunction——零化子与初值条件的封装构造函数签名holonomic.pyHolonomicFunction(annihilator, x, x00, y0None)annihilatorDifferentialOperator对象满足L·f 0x函数自变量x0初值条件所在点通常为整数默认0y0初值条件为保证函数唯一向量长度应大于等于微分方程阶数。初值条件的两种格式普通点初值向量形式[y(x₀), y(x₀), y(x₀), ...]例如sin(x)的[0, 1]、exp(x)的[1]。正则奇异点字典格式{s0: [C₀, C₁, ...], s1: [C¹₀, C¹₁, ...], ...}其中s0, s1, ...是**指标方程indicial equation**的根各向量是对应幂级数的首项系数。源码中is_singularics()方法即据此区分dict为奇异初值、list为普通初值holonomic.py。 from sympy.holonomic.holonomic import HolonomicFunction, DifferentialOperators from sympy import QQ from sympy import symbols, S x symbols(x) R, Dx DifferentialOperators(QQ.old_poly_ring(x),Dx) p HolonomicFunction(Dx - 1, x, 0, [1]) # e^x q HolonomicFunction(Dx**2 1, x, 0, [0, 1]) # sin(x) # 加法自动求出 e^x sin(x) 的零化子与初值 p q HolonomicFunction((-1) (1)*Dx (-1)*Dx**2 (1)*Dx**3, x, 0, [1, 2, 1]) # 乘法自动求出 e^x * sin(x) 的零化子与初值 p * q HolonomicFunction((2) (-2)*Dx (1)*Dx**2, x, 0, [0, 1])加法与乘法的实现分别位于 holonomic.py 的__add__与 holonomic.py 的__mul__它们以两个零化子为行构造线性方程组用DomainMatrix上的高斯消元_find_nonzero_solution求非平凡解若解为零矩阵则逐步抬高阶数重试直至求出结果零化子初值则通过_extend_y0扩展后按对应规则加法对应相加、乘法对应二项式卷积合成。奇异初值的完整示例类文档中给出了指标方程仅有一个根1/2的例子 HolonomicFunction(-S(1)/2 x*Dx, x, 0, {S(1)/2: [1]}) HolonomicFunction((-1/2) (x)*Dx, x, 0, {1/2: [1]}) HolonomicFunction(-S(1)/2 x*Dx, x, 0, {S(1)/2: [1]}).to_expr() sqrt(x)零化子x·Dx - 1/2的解正是sqrt(x)而to_expr()会把全纯函数还原为初等函数表达式。to_expr的实现为hyperexpand(self.to_hyper()).simplify()holonomic.py即先转换为超几何级数再展开化简。基于零化子的常用操作一览一旦完成「零化子 初值」的表示模块便能在其上执行大量闭包运算全纯函数对加减乘、微分、积分封闭操作方法示例源码 docstring 提供微分diff()HolonomicFunction(Dx**2 1, x, 0, [0, 1]).diff().to_expr()→cos(x)不定积分integrate((x, 0, x))HolonomicFunction(Dx - 1, x, 0, [1]).integrate((x, 0, x))→e^x - 1定积分integrate((x, a, b))上限为数字时返回数值上限为x时返回新的HolonomicFunction复合composition(expr)与代数函数复合初值可另行提供初值平移change_ics(b)expr_to_holonomic(sin(x)).change_ics(1)→ 初值变为[sin(1), cos(1)]转初等函数to_expr()HolonomicFunction(x**2*Dx**2 x*Dx (x**2 - 1), x, 0, [0, S(1)/2]).to_expr()→besselj(1, x)数值计算evalf(points)默认 RK4 方法沿给定点路径做数值积分支持复数点可用于 matplotlib 绘图其中integrate的实现很有代表性holonomic.py不定积分只需将零化子右乘Dself.annihilator * D并把初值前插一个0得到新初值定积分则先求原函数再代入上下限NaN时回退为求极限无法符号化时改用evalf求值。evalf的用法示例来自类 docstring绘制sin(x)**2/ximport sympy.holonomic # 注册 expr_to_holonomic from sympy import var, sin import matplotlib.pyplot as plt import numpy as np var(x) r np.linspace(1, 5, 100) y sympy.holonomic.expr_to_holonomic(sin(x)**2/x, x01).evalf(r) plt.plot(r, y, labelholonomic function) plt.show()expr_to_holonomic是另一个入口holonomic.py给定任意可全纯化的表达式它会自动反解出零化子并计算相应初值无需手工推导微分方程。用测试用例印证表示的正确性仓库测试 test_holonomic.py 大量使用了与本文相同的构造模式R, Dx DifferentialOperators(ZZ.old_poly_ring(x), Dx)或QQ.old_poly_ring(x)分别测试整数环与有理数域系数随后断言零化子的相等性如assert q.annihilator p.annihilator并通过to_expr()把结果还原为标准函数做对拍。这意味着本文所有示例语法——包括D**2 1的构造、HolonomicFunction(..., 0, [0, 1])的参数顺序、{S(1)/2: [1]}的奇异初值字典格式——都被 CI 持续验证可直接复制运行。小结表示方案的完整心智模型零化子L p₀ p₁D ... pᵣD^r是「微分方程」的算子化封装由DifferentialOperator多项式列表 父代数承载系数位于基础环ZZ[x]或QQ[x]代数DifferentialOperatorAlgebra是带非交换生成元D的 Weyl 代数交换规则Dx·a a·Dx a是全部算子运算的基石函数HolonomicFunction把零化子、自变量、初值点x0与初值y0绑定普通点用列表、正则奇异点用{指标根: 首项系数}字典闭环使用构造 → 代数运算、*、diff、integrate→ 还原to_expr、evalf全程无需离开该表示体系。更深入的背景全纯函数定义与闭包性质可阅读 about.rst各类运算的完整说明见 operations.rst从普通表达式与超几何/Meijer G 函数反解全纯形式的入口见 convert.rst。【免费下载链接】sympyA computer algebra system written in pure Python项目地址: https://gitcode.com/GitHub_Trending/sy/sympy创作声明：本文部分内容由AI辅助生成（AIGC），仅供参考