一种计算机代数系统的设计与实现
2023-12-14汪明
汪 明
（孔智科技（徐州）有限公司 研发部，江苏 徐州 221118）
0 引言
计算机代数是一门主要研究符号计算算法的学科，它所运算的基本对象是符号，而不是数值。[1]这里的符号可以看作变量，用以表示数学领域中的众多对象，如实数、复数、多项式、方程、函数、微积分和矩阵等。而数值计算处理的基本对象是数值，由于计算机内部数值是采用浮点类型进行表示的，因此，在运算过程中会出现结果不精确的情况。由于计算机代数系统可以对数学公式进行自动推导，也可以求解某些方程的精确解，因此，在科研和工程计算中都发挥着至关重要的作用。
通过对国内外相关资料进行整理发现，国外对于计算机代数系统的研究起步较早，约始于20世纪60 年代，经过几十年的研究积累，创造了在业界享有盛誉的产品，如Maple、Mathematica 和Matlab。目前，国内对于计算机代数系统的研究还处于初级阶段，很多学者的研究都是基于计算机代数系统的应用研究，而对其内部技术原理的研究较为鲜见。[2]
从文献上来看，计算机代数系统的研究始于James Slagle 创建的Saint 程序，它能够对函数进行简单的积分[3]，在此基础上，Moses Joel 改进了Saint 程序的算法，并建立了一个可以进行符号积分的程序SIN[4]。1960 年，Anthony C.Hearn 基于Lisp 语言开发了一款名为Reduce 的计算机代数系统。[5-6]1968 年，麻省理工学院的Carl Engleman和Moses Joel 等人开发了Macsyma 系统。[7]1980年，滑铁卢大学研究人员基于C 语言开发了Maple系统。[4，8]1982 年，William Schelter 开始对Macsyma系统进行分支，并命名为Maxima 系统。1986 年，Stephen Wolfram 用C 语言实现了Mathematica 系统。[9]1989 年，MathWorks 公司发布了Matlab 中的Symbolic Math Toolbox 工具箱，可以执行微分、积分、化简、变换和方程求解。[10]2006 年，Ondrˇej Cˇertík 基于Python 语言开发了SymPy 计算机代数系统。2007 年，Waldek Hebisch 基于Axiom 系统分支出了FriCAS 系统，并持续进行开发和改进。
1 计算机代数系统设计
1.1 系统功能设计
开发一个功能完备、运行稳定的计算机代数系统是非常大的工程，除了需要具备计算机领域的专业知识外，还需要对数学领域的相关理论比较熟悉。
本文设计的计算机代数系统重点是探索符号计算相关的可行技术路径，并不完全覆盖某一个数学领域的所有推导规则。为了方便软件的部署，这里采用B/S 架构，用户可以通过浏览器与系统进行交互，客户端和服务器通过Ajax 技术进行数据通信。另外，系统基于F# SDK 进行开发。F#是微软推出的一种开源免费的函数式编程语言，语法类似于OCaml。它可以非常方便地实现特定领域语言（Domain Specific Language，DSL），并进行模式匹配，以达到数学表达式解析和执行的目的。计算机代数系统的总体功能设计示意图见图1。
图1 系统的总体功能设计示意图
1.1.1 解析器模块
实现DSL 语言的词法、语法和语义分析，并构建抽象语法树（Abstract Syntax Tree，AST）。
1.1.2 化简模块
可基于多种化简规则实现对表达式的化简操作。一般来说，在经过多次递归的符号计算过程后，计算结果会出现多层嵌套的现象，且不够简洁。而化简后，可让结果看起来更加直观。
1.1.3 展开模块
实现对表达式的展开操作，特别是对多个表达式的乘法进行展开得到多项式。
1.1.4 求导模块
实现对表达式的求导操作，它是一个递归方法，内部将常见的求导规则进行逐条归纳，并支持复合函数的求导处理。
1.1.5 极限模块
实现对表达式的求极限操作，与求导模块类似，内部需要归纳多个极限规则，并支持复合函数的极限处理。
1.1.6 泰勒级数模块
实现对表达式的泰勒级数展开操作，内部基于求导模块，可根据参数实现N 阶泰勒级数展开。
1.1.7 积分模块
实现对表达式的不定积分操作，不定积分比求导复杂得多，主要理由有两点：首先，不定积分的推导规则非常多（需维护一个符号积分表）；其次，很多不定积分推导规则涉及条件判断，有的不定积分甚至需要借助复杂的算法才能求解。
1.1.8 打印模块
实现对抽象语法树的文本打印操作，它会递归解析语法树，并根据规则输出文本格式，其中的难点是对括号的优化，它需要根据运算符的优先级和结合性，来尽量减少括号的嵌套。
1.1.9 绘图模块
实现对表达式的绘图操作，目前对2D 图形支持较好，内部可对表达式进行解析，并用循环生成离散的二维数据，借助XPlot 库实现2D 图形的绘制，它会生成HTML 格式的图形数据，之后返回到前端并显示。
1.1.10 求值模块
实现对表达式的求值操作，如计算sin（3）*cos（7）的值。下面给出主要功能模块对应的函数设计说明，具体如表1 所示。
表1 系统主要功能函数
1.2 领域特定语言设计
领域特定语言指的是专注于某个应用程序领域的计算机编程语言。对于计算机代数系统而言，设计与实现一种可处理数学表达式的编译器是一个关键和难点。[2]一般来说，数学领域的公式sin（x^2）-x 可以看作是一个表达式，而表达式可以被计算机代数算法处理，比如对两个表达式进行四则运算，或者对其进行求导和化简。表达式主要由以下要素组成：符号变量、常数和操作符。符号变量本质是一个字符串，类似于编程语言中的变量。常见的变量为：α、β、x、y。常数可以表示数值（如7、3.0、3.5），也可以表示数学常量（如e）。操作符包括中缀运算符和前缀运算符，中缀运算符如+、-、*、/；前缀运算符可以看作是一种数学函数的应用，如sin，cos，sqrt。下面给出DSL 语言的BNF语法描述（语法1）。
语法1：
|
| "csc" | "ln" | "..." | "exp"
:：=
| "5" | "6" | "7" | "8" | "9"
1.3 解析器设计
有了DSL 语法定义后，下一步的核心工作就是对解析器进行设计。这里基于FParsec 库，可对多个语法解析器进行组合，从而形成快速处理数学表达式的解析器。
FParsec 库中的OperatorPrecedenceParser 实例对象，具有AddOperator 方法，可以让用户自行设置操作符的解析规则。下面给出操作符的解析规则说明，具体如表2 所示。
表2 操作符的解析规则
表2 中的exp 操作符上一行的...表示省略号，而不是操作符。F#语言支持用type 关键字自定义复合数据类型，结合语法1 给出的BNF 语法描述，下面给出F#语言对数学表达式类型的定义。
代码1：
type Expr =
| Const of float | Var of string
| Add of Expr * Expr | Subt of Expr * Expr
| Prod of Expr * Expr | Div of Expr * Expr
| Pow of Expr * Expr | Sin of Expr
| Cos of Expr | Tan of Expr | Sec of Expr
| Cot of Expr | Csc of Expr
| ArcSin of Expr | ArcCos of Expr
| Exp of Expr | Ln of Expr | Neg of Expr
| Diff of Expr * Expr
|...
| E | Pi | I | PositiveInfinity | NegativeInfinity
利用设计的解析器可以将数学表达式sin（x^2）+log（3×y）-5 解 析 为：Subt（Add（Sin（Pow（Var "x"，Const 2.0）），Log（Prod（Const 3.0，Var "y"））），Const 5.0）。为了更直观地展示，下面给出解析过程示意图（见图2）。
图2 解析过程示意图
1.4 模式匹配算法
当用户输入文本命令后，解析器会对文本进行语义解析，并返回抽象语法树。为了让计算机能自动进行公式推理，还需要对抽象语法树进行模式匹配，最终根据事先定义的推理规则求出符号计算的结果。
F#本身提供了功能非常强大的模式匹配算法功能，可以对表达式进行模式匹配。下面给出F#模式匹配的语法，如语法2 所示。
语法2：
// Match Expression
match expression with
| rule1Left [ when condition1 ] -> rule1Right
| rule2Left [ when condition2 ] -> rule1Right2
|...
从语法2 可知，F#模式匹配中还提供了可选的when 条件判断支持，这样就可以更加灵活地对匹配规则进行细节处理。
计算机代数系统中所涉及的数学表达式自动推导，从本质上来说，都需要事先定义好相关的推理规则。当然，这个推理规则需要结合前面的DSL语义规则以及解析器来定制。下面给出符号计算模式匹配算法：
（1）模块M 从外部获取抽象语法树AST；
（2）遍历系统内基于抽象语法树构建的推理规则表T，利用F#的模式逐条进行规则匹配；
（3）依次取出当前推理规则R 中的左部RL和右部RR，如果RL 匹配AST，则解析出匹配各参数的值，并替换到RR 中；
（4）如果匹配的RR 中无待递归的函数，则返回当前结果。否则，返回第（1）步。
为了便于对该算法的理解，下面结合diff（sin（x^2），x）的求导过程来说明：
首先，给出3 条基于抽象语法树构建的求导规则：
（1）Diff（Sin（e），Var "x"） → Prod（Cos（e），Diff（e，Var "x"））
（2）Diff（Pow（e，Const m），Var "x"） → Prod（Prod（m，Pow（e，Const（m-1））），Diff（e，Var "x"））
（3）Diff（Var "x"，Var "x"）→Const 1
diff（sin（x^2），x）经解析器解析后，生成的抽象语法树为：
Diff（Sin（Pow（Var "x"，Const 2.0）），Var "x"）
其次，遍历求导规则表，匹配规则（1）的左部，其中的参数e 匹配值Pow（Var "x"，Const 2.0），返回规则（1）的右部，并替换参数e，即：
Prod（Cos（Pow（Var "x"，Const 2.0）），Diff（Pow（Var "x"，Const 2.0），Var "x"））
由于右部表达式中有Diff（Pow（Var "x"，Const 2.0），Var "x"），则继续执行求导过程，匹配规则（2）的左部，其中的参数m 值为2.0，参数e 为Var "x"，并返回规则（2）的右部，替换参数m 和e，即：
Prod（Cos（Pow（Var "x"，Const 2.0）），
Prod（Prod（Const 2.0，Pow（Var "x"，Const（2.0-1））），
Diff（Var "x"，Var "x"）
））
由于右部表达式中有Diff（Var "x"，Var "x"），因此继续执行求导过程，匹配规则（3）的左部，并返回规则（3）的右部Const 1，即：
Prod（Cos（Pow（Var "x"，Const 2.0）），
Prod（Prod（Const 2.0，Pow（Var "x"，Const（2.0-1））），
Const 1
））
此时的表达式无Diff 函数的调用，因此返回。上述表达式经打印处理后为：
cos（x^2）*2*x^（2-1）*1 →cos（x^2）*2*x
上述算法适用于符号计算中的求导、求积分、求极限、化简、打印以及其他基于规则的功能模块。
2 系统实验分析
2.1 系统运行环境
基于F#语言实现计算机代数系统，需要一定的软硬件环境。下面给出系统运行的软硬件平台配置，如表3 所示。
表3 实验软硬件平台参数
2.2 系统效果分析
为了验证系统能否对数学表达式进行正确的自动推导，下面分类别给出多个测试题，并将系统推理的结果和SymPy 求出的正确结果进行对比分析。
如果推理结果二者一致，则评价为正确，此时在系统测试评价表的评价一栏中打勾。具体实验对比分析结果如表4 所示。
表4 系统测试评价表
下面给出本计算机代数系统的交互Web 界面，其中的计算结果会转换成LaTeX 文本格式，Web 端利用KaTeX 库进行渲染呈现，交互界面如图3 所示。
图3 系统交互Web 界面图
通过实验对比分析可知，只要是系统中对符号计算规则进行了归纳，系统就可以自动且正确求出数学表达式的推导结果，由此可说明本文所探讨的技术路线是可行的，且系统运行效果良好。
3 结语
通过实验发现，基于F#强大的模式匹配能力和FParsec 库实现的解析器，可以对数学表达式进行符号计算。当前，本计算机代数系统还相对初级，功能并不完备，但其中的设计思想和技术路线是可行的。
未来，需要引入更多的推导规则和符号计算算法，特别是不定积分的规则。下一步研究的重点工作有以下两个方面：一方面，可以参考Mathematica 的积分库Rubi[11，12]，其规则大约有5000 条；另一方面，对于Rubi 规则库不能求出的积分，可引入Risch 算法[13]进行求解。
