首页 > 解决方案 > 如何使用 Tuple/Array/Vector 从 Python (ctypes) 调用 PARI/GP?

问题描述

我想从 Python调用PARI/GP 。我需要使用ellisdivisible(E; P; n;{&Q})PARI 的功能 (请参阅此链接中第 441 页上的第 3.15.35 号功能:),所以我必须传递 2 个向量或数组(例如, E = ellinit([0,-1,1,0,0], K);P = [0,0];),我该怎么做?

要从 Python (由 Thomas Baruchel 给出)调用单个参数/变量的 PARI 函数(在 C 中),我们有如下内容 -

import ctypes

# load the library
pari=ctypes.cdll.LoadLibrary("libpari.so")

# set the right return type of the functions
pari.stoi.restype = ctypes.POINTER(ctypes.c_long)
pari.nextprime.restype = ctypes.POINTER(ctypes.c_long)

# initialize the library 
pari.pari_init(2**19,0)

def nextprime(v):
  g = pari.nextprime(pari.stoi(ctypes.c_long(v))) # nextprime(argument) is a PARI function
  return pari.itos(g)



print( nextprime(456) )

例如我试过 -

h=(0,0,0, 4,6)
pari.stoi.restype = ctypes.POINTER(ctypes.c_long*5)
pari.ellinit.restype = ctypes.POINTER(ctypes.c_long)
def ellinit(v):
  g = pari.ellinit(pari.stoi(ctypes.c_long(v)*5))
  return pari.itos(g)


print(ellinit(h))

我得到以下错误 -

   File "C:\Users\miron\Desktop\trash5\x\f.py", line 68, in <module>
    print( ellinit(h) )
  File "C:\Users\miron\Desktop\trash5\x\f.py", line 62, in ellinit
    g = pari.ellinit(pari.stoi(ctypes.c_long(v)*5))
TypeError: an integer is required (got type tuple)

如何传递元组/数组/向量?谢谢。

编辑: 尝试获取失败ellisdivisible(E; P; n;{&Q})-

from ctypes import *

pari = cdll.LoadLibrary("C:\\Program Files\\Python37\\Pari64-2-11-3\\libpari.dll")

pari.stoi.restype = POINTER(c_long)
pari.cgetg.restype = POINTER(POINTER(c_long))
pari.ellinit.restype = POINTER(POINTER(c_long))
#-------------------------CHANGE 1
pari.ellisdivisible.restype = c_long
Q = pari.stoi(c_long(0))
#-------------------------
(t_VEC, t_COL, t_MAT) = (17, 18, 19)  # incomplete
precision = c_long(38)
pari.pari_init(2 ** 19, 0)
def t_vec(numbers):
    l = len(numbers) + 1
    p1 = pari.cgetg(c_long(l), c_long(t_VEC))
    for i in range(1, l):
        p1[i] = pari.stoi(c_long(numbers[i - 1]))
    return p1
def main():
    h = (0, 0, 0, 0, 1)
    P=(0,0)
    res = pari.ellinit(t_vec(h), pari.stoi(c_long(1)), precision)
#---------------CHANGE 2
   # res = pari.ellinit(t_vec(h), pari.stoi(c_long(1)), precision).disc
    y = pari.ellisdivisible(res, t_vec(P), pari.stoi(c_long(5)), byref(Q))
    print(pari.itos(y))
#---------------
    for i in range(1, 13):
        print(pari.itos(res[i]))

if __name__ == '__main__':
    main()

错误是 -

Traceback (most recent call last):
  File "C:\Users\miron\Desktop\trash5\x\ex - Copy (2).py", line 34, in <module>
    main()
  File "C:\Users\miron\Desktop\trash5\x\ex - Copy (2).py", line 28, in main
    print(pari.itos(y))
OSError: exception: access violation reading 0x0000000000000009

标签: pythoncarraysctypespari

解决方案


Python 元组或 C 数组不能直接使用,因为 PARI 正在使用 PARI/GP 特定向量,其中类型/长度编码在第一个元素中。

在第 4.4.1 节Creation of PARI objects中它说:

创建 PARI 对象的基本函数是 GEN cgetg(long l, long t) l 指定要分配给对象的长字数,t 是对象的类型,以符号形式(参见第 4.5 节的列表这些)。这个函数的确切作用如下:它首先在 PARI 栈上创建一个长度为 longwords 的内存块,并保存它最终返回的块的地址。如果堆栈已用完,则打印“PARI 堆栈溢出”的消息,并引发错误。否则,它设置 PARI 对象的类型和长度。实际上,它填充了它的第一个代码字 (z[0])。

https://pari.math.u-bordeaux.fr/pub/pari/manuals/2.7.6/libpari.pdf

在本文档的示例中,您可以看到要创建一个包含两个元素的向量,调用它的大小为l=3以获得合适的向量。实际数字向量的第一个元素不是从索引 0 开始,而是从索引 1 开始(参见本 PDF 文档中的第 4.5.15 节)。

git clone http://pari.math.u-bordeaux.fr/git/pari.git   

可以获取 PARI 的源代码。

您可以在 src/headers/parigen.h 的末尾看到不同的类型。它是一个枚举,我们需要的类型是 t_VEC。对应的整数是 17。

因此,我们现在可以定义一个将元组转换为 GP 向量的小函数,如下所示:

def t_vec(numbers):
    l = len(numbers) + 1
    p1 = pari.cgetg(c_long(l), c_long(t_VEC))
    for i in range(1, l):
        p1[i] = pari.stoi(c_long(numbers[i - 1]))
    return p1

然后我们可以这样称呼ellinit

h = (0, 0, 0, 4, 6)
res = pari.ellinit(t_vec(h), pari.stoi(c_long(1)), precision)

为了使用您的 [0, 0, 0, 4, 6] 参数对其进行测试,您可以从命令行调用 GP:

? ellinit([0, 0, 0, 4, 6], 1)
%1 = [0, 0, 0, 4, 6, 0, 8, 24, -16, -192, -5184, -19648, 110592/307, Vecsmall([1]), [Vecsmall([128, -1])], [0, 0, 0, 0, 0, 0, 0, 0]]

引用的 PDF 文档第 441 页上示例的小型、独立的 Python 程序可能如下所示:

from ctypes import *

pari = cdll.LoadLibrary("libpari.so")

pari.stoi.restype = POINTER(c_long)
pari.cgetg.restype = POINTER(POINTER(c_long))
pari.ellinit.restype = POINTER(POINTER(c_long))
pari.ellisdivisible.restype = c_long
pari.nfinit0.restype = POINTER(c_long)
pari.polcyclo_eval.restype = POINTER(c_long)
pari.fetch_user_var.restype = c_long
pari.pol_x.restype = POINTER(c_long)

(t_VEC, t_COL, t_MAT) = (17, 18, 19)  # incomplete
precision = c_long(38)

pari.pari_init(2 ** 19, 0)


def t_vec(numbers):
    l = len(numbers) + 1
    p1 = pari.cgetg(c_long(l), c_long(t_VEC))
    for i in range(1, l):
        p1[i] = pari.stoi(c_long(numbers[i - 1]))
    return p1


def main():
    t = pari.pol_x(pari.fetch_user_var(bytes("t", "utf8")))
    Q = pari.pol_x(pari.fetch_user_var(bytes("Q", "utf8")))
    K = pari.nfinit0(pari.polcyclo_eval(11, t), c_long(0), precision)
    h = (0, -1, 1, 0, 0)
    res = pari.ellinit(t_vec(h), K, precision)
    P = (0, 0)
    y = pari.ellisdivisible(res, t_vec(P), pari.stoi(c_long(5)), byref(Q))

    pari.pari_printf(bytes("Q: %Ps\n", "utf8"), Q)

    print("ellisdivisible =", y)


if __name__ == '__main__':
    main()

测试

现在我们可以调用 Python 程序,并将其与交互式 GP 程序的输出进行比较,它实际上给出了相同的结果:

Q: [Mod(-t^7 - t^6 - t^5 - t^4 + 1, t^10 + t^9 + t^8 + t^7 + t^6 + t^5 + t^4 + t^3 + t^2 + t + 1), Mod(-t^9 - 2*t^8 - 2*t^7 - 3*t^6 - 3*t^5 - 2*t^4 - 2*t^3 - t^2 - 1, t^10 + t^9 + t^8 + t^7 + t^6 + t^5 + t^4 + t^3 + t^2 + t + 1)]
ellisdivisible = 1

推荐阅读