python解超越方程_初试在Python中使用PARI/GP

绝地灬酷狼 2022-11-01 00:52 65阅读 0赞

BoJone很喜欢Python，也很喜欢数论，所以就喜欢利用Python玩数论了。平时也喜欢自己动手写一些数论函数，毕竟Python支持大整数高精度运算，这点是非常好的；但是，在很多实际应用中，还是希望能有一个现成的数论函数库来调用。之前尝试过数学研发网的HugeCalc库，但是由于各种不熟悉不了了之。后来论坛上的无心老兄推荐了PARI/GP，小试一下，居然在Python上成功调用了。以后再也不用担心Python上的数论计算问题了，呵呵~

先来一段官方的关于PARI/GP的介绍：PARI/GP 是一种针对数论中的快速计算(大数分解，代数数论，椭圆曲线...) 而设计的广泛应用的计算机代数系统，同样具备大量实用的函数来对于数学实体的计算， 诸如矩阵，多项式，幂级数，代数数，以及相当多的超越方程等等。 PARI也可以作为快速计算的C语言库。

这个系统最初是由Henri Cohen和他的合作者们(来自法国的波尔多第一大学)发 展起来的。 PARI现在在公共授权条款(GPL)之下并在许多志愿者的帮助下由Karim Belabas进行维护。

\-PARI 是一个可以快速计算的C语言库。

\-gp是一个容易上手的交互式外壳，它赋予用户调用PARI函数的权力。

\-GP 是gp脚本语言的名字。

\-gp2c, GP到C的编译器，通过编译GP脚本至C语言以及透明载入结果函数至gp来结合两者的最优点。 (由gp2c编译的脚本一般情况下运行速度将快3至4倍。) 目前gp2c只能运行GP语言的一部分。

PARI/GP提供了源码包，可以编译使用，也提供了Windows下的可执行安装包，直接安装使用也可。本文只讨论在Windows平台下用Python调用PARI，Linux平台估计可以类似地调用，但写这篇文章之前，我还没亲自测试过。Windows下的调用只需要PARI/GP安装目录下的libpari.dll文件，读者可以从官网下载安装包安装后提取它，也可以在本文末尾的附件中下载它。

首先，需要注意的是，调用PARI的Python版本必须是32位的，Python 2和Python 3均可，但64位会出错(可以在64位的操作系统中，安装32位的Python)，当然，也可以用PyPy来调用。PyPy基于32位的Python，最新版本已经支持Python 3，读者可以选择自己喜欢的。先看在Python 2.x中的方法，3.x在细节上有些差别，在后面会说明。

DLL是动态链接库，加载DLL之后，就可以使用DLL里边的函数。这跟加载模块有点类似。不同的是，DLL一般是用C语言写的，运行效率比较高，而在Python中加载C/C++的函数，可以发挥C/C++的运行效率和Python的开发效率。在Python中加载DLL，需用导入ctypes模块，导入之后，利用ctypes.cdll.LoadLibrary()函数载入libpari.dll(linux下也有类似的链接库，后缀为.so)。代码如下

在调用PARI的功能之前，首先得用它的pari\_init()函数初始化，比如pari.pari\_init(4000000,2)，其中第一个参数表示提供给PARI使用的字节数，也就是PARI最大内存占用，一般来说不要小于500000，第二个参数是预先生成的素数表，预先生成部分素数表可以在之后处理素数相关问题时，降低一点计算量，但是这个不是必须的，你也可以设它为0，它同样还是可以完成素数相关运算。

下面以使用PARI的primes()函数为例。prime(n)函数的作用是输出前n个素数。代码如下

由于只是刚刚尝试，我也不是很懂restype的作用，望文生义，应该就是把函数的类型跟指针相对应。这代码是参考网络上的一些代码而来的，我测试也发现，就算没有restype的这两句，程序也能够正常工作，不知道有什么不同。程序的过程是先用primes()函数生成结果，结果是GEN类型的，在Python中没法识别，用GENtostr()函数转换为字符串，但这函数返回的是地址，要用string\_at()读取输出。在Python 2.7中，运行结果为(输出的是字符串)：\[2, 3, 5, 7, 11, 13, 17, 19, 23, 29\]

以上代码在Python 3同样可以正常运行，但是输出结果有小许不同。以下会说明。

第二个例子是isprime()函数，这是素性测试函数。这种函数有些特殊，它不是用C/C++的int类型，因此就不能简单输入Python的整数，估计PARI自定义了自己的大整数类型，或者用其他的方式储存大整数。我们要从字符串中读入数字，代码如下：

核心代码就是倒数第二句，首先将127转化为字符串，然后用gp\_read\_str()函数从字符串读入，最后用isprime()函数判断，本例的输出结果为1。

但是，如果把这个例子放到Python 3中运行，就会出现错误。原因在于，Python 3的字符串str，不等于Python 2.x的str了，在3中定义了一种新的类型bytes，事实上，2.x的str在这里更像3中的bytes，所以，在Python中，代码要改成

或者(用str的encode()函数转成bytes类型)

同样，在Python 3中，ctypes.string\_at(y)的返回对象也不是str，而是bytes，读者可以发现，输出结果多带了一个b，正是这个原因。

暂时就介绍到这里了，如果有新的发现，会继续跟大家分享~~希望高手路过多多指教。

更详细的转载事宜请参考：《科学空间FAQ》

如果您还有什么疑惑或建议，欢迎在下方评论区继续讨论。

如果您觉得本文还不错，欢迎分享/打赏本文。打赏并非要从中获得收益，而是希望知道科学空间获得了多少读者的真心关注。当然，如果你无视它，也不会影响你的阅读。再次表示欢迎和感谢！

打赏

微信打赏