通过 GF(2) 快速计算矩阵秩

Question

对于我当前的项目，我需要能够计算具有 GF(2) 条目的 64*64 矩阵的排名。我想知道是否有人有好的解决方案。

我一直在使用pyfinite来实现这一点，但它相当慢，因为它是一个纯 python 实现。我还尝试对我一直在使用的代码进行 cythonise，但由于依赖 pyfinite 而遇到了问题。

我的下一个想法是用 cython 编写我自己的类，但这对于我的需要来说似乎有点矫枉过正。

我需要以下功能

matrix = GF2Matrix(size=64) # creating a 64*64 matrix
matrix.setRow(i, [1,0,1....,1]) # set row using list
matrix += matrix2 # addition of matrices
rank(matrix) # then computing the rank

感谢您的任何想法。

Answer 1

在 GF(2) 上有效表示矩阵的一种方法是将行存储为整数，将每个整数解释为位串。例如，4×4 矩阵

\n\n

[0 1 1 0]\n[1 0 1 1]\n[0 0 1 0]\n[1 0 0 1]\n

\n\n

（排名为 3）可以表示为[6, 13, 4, 9]整数列表。这里我认为第一列对应于整数的最低有效位，最后一列对应于最高有效位，但相反的约定也可以工作。

\n\n

通过这种表示，可以使用 Python 的按位整数运算高效地执行行运算：^加法、&乘法。\n然后您可以使用标准高斯消元法计算秩。

\n\n

这是一些相当有效的代码。给定表示上述矩阵的非负整数集合rows，我们重复删除列表中的最后一行，然后使用该行消除1与其最低有效位对应的列中的所有条目。如果该行为零，则它没有最低有效位，并且对排名没有贡献，因此我们只需丢弃它并继续。

\n\n

def gf2_rank(rows):\n    """\n    Find rank of a matrix over GF2.\n\n    The rows of the matrix are given as nonnegative integers, thought\n    of as bit-strings.\n\n    This function modifies the input list. Use gf2_rank(rows.copy())\n    instead of gf2_rank(rows) to avoid modifying rows.\n    """\n    rank = 0\n    while rows:\n        pivot_row = rows.pop()\n        if pivot_row:\n            rank += 1\n            lsb = pivot_row & -pivot_row\n            for index, row in enumerate(rows):\n                if row & lsb:\n                    rows[index] = row ^ pivot_row\n    return rank\n

\n\n

让我们在 GF2 上运行一些随机 64×64 矩阵的计时。random_matrices是一个创建随机 64×64 矩阵集合的函数：

\n\n

import random\n\ndef random_matrix():\n    return [random.getrandbits(64) for row in range(64)]\n\ndef random_matrices(count):\n    return [random_matrix() for _ in range(count)]\n

\n\n

这是计时代码：

\n\n

import timeit\n\ncount = 1000\nnumber = 10\ntimer = timeit.Timer(\n    setup="ms = random_matrices({})".format(count),\n    stmt="[gf2_rank(m.copy()) for m in ms]",\n    globals=globals())\nprint(min(timer.repeat(number=number)) / count / number)\n

\n\n

在我的机器（2.7 GHz Intel Core i7、macOS 10.14.5、Python 3.7）上打印的结果是0.0001984686384，因此对于单列计算来说，这在 200\xc2\xb5s 以下。

\n\n

200\xc2\xb5s 对于纯 Python 排名计算来说是相当不错的，但如果这不够快，我们可以按照您的建议使用 Cython。这是一个 Cython 函数，它采用 dtype 的 1d NumPy 数组np.uint64，再次将数组的每个元素视为 GF2 上 64×64 矩阵的一行，并返回该矩阵的秩。

\n\n

# cython: language_level=3, boundscheck=False\n\nfrom libc.stdint cimport uint64_t, int64_t\n\ndef gf2_rank(uint64_t[:] rows):\n    """\n    Find rank of a matrix over GF2.\n\n    The matrix can have no more than 64 columns, and is represented\n    as a 1d NumPy array of dtype `np.uint64`. As before, each integer\n    in the array is thought of as a bit-string to give a row of the\n    matrix over GF2.\n\n    This function modifies the input array.\n    """\n    cdef size_t i, j, nrows, rank\n    cdef uint64_t pivot_row, row, lsb\n\n    nrows = rows.shape[0]\n\n    rank = 0\n    for i in range(nrows):\n        pivot_row = rows[i]\n        if pivot_row:\n            rank += 1\n            lsb = pivot_row & -pivot_row\n            for j in range(i + 1, nrows):\n                row = rows[j]\n                if row & lsb:\n                    rows[j] = row ^ pivot_row\n\n    return rank\n

\n\n

对 64×64 矩阵（现在表示为 dtypenp.uint64和 shape的 NumPy 数组）运行等效计时(64,)，我得到的平均排名计算时间为 7.56\xc2\xb5s，比纯 Python 版本快 25 倍以上。

\n