首页 > 解决方案 > Numpy 256 位计算

问题描述

我有一些用于椭圆曲线密码学的大数字。到目前为止,我一直在使用“bignum”整数类型在原生支持任意大小的数字的普通 python 中进行所有计算。但是,我想对其中一些操作进行矢量化(出于性能原因),并想使用一个包含 256 位无符号整数的数组。np.longlong但是 AFAIK 通过使用(或np.ulonglong对于无符号版本),numpy 中整数的最大大小为 64 位。以下帖子的一个想法是将每个 256 位数字表示为具有四个 64 位整数的数组:

large = 105951305240504794066266398962584786593081186897777398483830058739006966285013
arr = np.array([large.to_bytes(32,'little')]).view('P')
arr

输出:

array([18196013122530525909, 15462736877728584896, 12869567647602165677,
       16879016735257358861], dtype=uint64)

我可以通过执行以下操作将其转换回本机 python“bignum”:

int.from_bytes(arr.tobytes(), 'little')

输出:

105951305240504794066266398962584786593081186897777398483830058739006966285013

但是,我想在将其转换回“bignum”之前进行一些操作:

arr2 = arr + 1
int.from_bytes(arr2.tobytes(), 'little')

输出:

105951305240504794072543500697971467357257258687906003363414235534976478560982

这当然不等于large + 1


是否可以在 numpy 中处理 256 位数字?如果没有,是否存在任何合适的解决方法?

标签: pythonnumpy

解决方案


假设您有几个四个字节的数组:

a = np.array([0b0000_1000, 0b1000_1100, 0b0111_1111, 0b1111_0000], dtype=np.uint8)  # [  8, 140, 127, 240]
b = np.array([0b1111_0000, 0b1111_0000, 0b1111_0000, 0b1111_0000], dtype=np.uint8)  # [240, 240, 240, 240]

现在你想把它们加起来得到类似的东西

c = np.array([0b1111_1001, 0b0111_1101, 0b0111_0000, 0b1110_0000], dtype=np.uint8)  # [249, 125, 112, 224]

我故意避免需要在最高字节中结转的东西。在这种情况下,结果是用

ai = int.from_bytes(a.tobytes(), byteorder='big')
bi = int.from_bytes(b.tobytes(), byteorder='big')
c = np.frombuffer((ai + bi).to_bytes(4, 'big'), dtype=np.uint8)

那么如何在不转换为和 from 的情况下做到这一点int?你不能天真c = a + b

[1111_1000, 0111_1100, 0110_1111, 1110_0000]
          +1         +1         +1

请注意,标记位置的进位都被丢弃了。您可以使用检查进位而不触发溢出

0xFF - a < b

好消息是每个字节只能携带一位:0xFF + 0xFF + 1 == 0x1FF. 如果总和刚刚结转,您也不能通过加 1 结转:(0xFF + 0xFF) & 0xFF == 0xFE。这是正确添加的简单实现:

c = a.copy()
while b.sum():
    carry = np.r_[0xFF - c[1:] < b[1:], False].view(np.uint8)
    c += b
    b = carry

这个循环通常会运行 1-2 次迭代。在绝对最坏的情况下,它可能会运行多达 N-1 次迭代,N=256在你的情况下,或者N=4在这个例子中。使用更大的数据类型意味着您可以减少最坏情况下的迭代次数(并且可能也加快最好情况下的迭代次数)。唯一需要注意的是,您将无法使用该...view(np.uint8)技巧来转换布尔掩码。但这实际上是因祸得福,因为你可以写

c = a.copy()
carry = np.zeros_like(a)
cbuf = carry.view(bool)[::a.itemsize]
# On big-endian use this instead:
#cbuf = carry.view(bool)[a.itemsize - 1::a.itemsize]
cbuf = cbuf[:-1]
m = np.iinfo(a.dtype).max
while b.sum():
    c += b
    np.less(m - c[1:], b[1:], out=cbuf)
    b = carry

请注意,此版本不仅避免了重新分配进位缓冲区,而且与itemsizeofa和完全无关b

读者练习:

  • 正确处理 little-endian 顺序(我的示例中的单词按 big-endian 顺序排列,与问题不同)
  • 添加符号支持(如果您有固定宽度的二进制补码,则开箱即用)
  • 实际上跨 N 维向量化这个操作
  • 实现诸如乘法之类的操作

推荐阅读