我想在 Python 中采用像 [[1,2],[3,4]] mod 7 这样的矩阵的模逆。我看过 numpy （它做矩阵求逆但不是模矩阵求逆），我在网上看到了一些数论包，但似乎没有什么可以做这个相对常见的过程（至少，它对我来说似乎相对常见）。
顺便说一下，上面矩阵的逆矩阵是[[5,1],[5,3]]（mod 7）。不过，我希望 Python 为我做这件事。
【问题讨论】：
sympy 内置了模块化矩阵求逆（可能是新问题，因为有人提出了这个问题），并且模块化行减少也相当容易。见stackoverflow.com/a/37015283/2747370
标签： python matrix number-theory matrix-inverse
好的...对于那些关心的人，我解决了我自己的问题。我花了一段时间，但我认为这是可行的。它可能不是最优雅的，应该包含更多的错误处理，但它确实有效：
import numpy
import math
from numpy import matrix
from numpy import linalg
def modMatInv(A,p): # Finds the inverse of matrix A mod p
n=len(A)
A=matrix(A)
adj=numpy.zeros(shape=(n,n))
for i in range(0,n):
for j in range(0,n):
adj[i][j]=((-1)**(i+j)*int(round(linalg.det(minor(A,j,i)))))%p
return (modInv(int(round(linalg.det(A))),p)*adj)%p
def modInv(a,p): # Finds the inverse of a mod p, if it exists
for i in range(1,p):
if (i*a)%p==1:
return i
raise ValueError(str(a)+" has no inverse mod "+str(p))
def minor(A,i,j): # Return matrix A with the ith row and jth column deleted
A=numpy.array(A)
minor=numpy.zeros(shape=(len(A)-1,len(A)-1))
p=0
for s in range(0,len(minor)):
if p==i:
p=p+1
q=0
for t in range(0,len(minor)):
if q==j:
q=q+1
minor[s][t]=A[p][q]
q=q+1
p=p+1
return minor
【讨论】：
当舍入错误不是问题时有效的黑客技巧：
一种不那么骇人听闻的方法是实际实现高斯消除。这是我使用高斯消除的代码，我为自己的目的编写的（舍入错误对我来说是个问题）。 q 是模数，不一定是素数。
def generalizedEuclidianAlgorithm(a, b):
if b > a:
return generalizedEuclidianAlgorithm(b,a);
elif b == 0:
return (1, 0);
else:
(x, y) = generalizedEuclidianAlgorithm(b, a % b);
return (y, x - (a / b) * y)
def inversemodp(a, p):
a = a % p
if (a == 0):
print "a is 0 mod p"
return None
if a > 1 and p % a == 0:
return None
(x,y) = generalizedEuclidianAlgorithm(p, a % p);
inv = y % p
assert (inv * a) % p == 1
return inv
def identitymatrix(n):
return [[long(x == y) for x in range(0, n)] for y in range(0, n)]
def inversematrix(matrix, q):
n = len(matrix)
A = np.matrix([[ matrix[j, i] for i in range(0,n)] for j in range(0, n)], dtype = long)
Ainv = np.matrix(identitymatrix(n), dtype = long)
for i in range(0, n):
factor = inversemodp(A[i,i], q)
if factor is None:
raise ValueError("TODO: deal with this case")
A[i] = A[i] * factor % q
Ainv[i] = Ainv[i] * factor % q
for j in range(0, n):
if (i != j):
factor = A[j, i]
A[j] = (A[j] - factor * A[i]) % q
Ainv[j] = (Ainv[j] - factor * Ainv[i]) % q
return Ainv
编辑：正如评论者指出的那样，在某些情况下该算法会失败。修复它有点不重要，而且我现在没有时间。当时它适用于我的情况下的随机矩阵（模数是大素数的乘积）。基本上，第一个非零条目可能不是模数的相对质数。主要情况很容易，因为您可以搜索不同的行并交换。在非素数情况下，我认为可能是 所有 前导条目不是相对素数，因此您必须将它们组合起来
【讨论】：
long 替换为 int 并将 / 替换为 //（最后一个真的让我措手不及）。
inversemodp 函数可能会返回一个值即使 a 没有逆模 p。例如，inversemodp(2, 4) 返回 1，但 2*1 =/= 1 mod 4。
可以使用 Sage (www.sagemath.org) as 来计算
Matrix(IntegerModRing(7), [[1, 2], [3,4]]).inverse()
虽然 Sage 安装起来很庞大，但你必须使用它附带的 python 版本，这很痛苦。
【讨论】：
不幸的是 numpy 没有模运算实现。您始终可以使用行缩减或行列式对建议的算法进行编码，如here 所示。模逆似乎对密码学非常有用。
【讨论】：
'sympy' 包矩阵类函数 'sqMatrix.inv_mod(mod)' 计算小模和任意大模的模矩阵逆。通过将 sympy 与 numpy 相结合，计算二维 numpy 数组的模逆变得容易（参见下面的代码 sn-p）：
enter code here
import numpy
from sympy import Matrix
def matInvMod (vmnp, mod):
nr = vmnp.shape[0]
nc = vmnp.shape[1]
if (nr!= nc):
print "Error: Non square matrix! exiting"
exit()
vmsym = Matrix(vmnp)
vmsymInv = vmsym.inv_mod(mod)
vmnpInv = numpy.array(vmsymInv)
print "vmnpInv: ", vmnpInv, "\n"
k = nr
vmtest = [[1 for i in range(k)] for j in range(k)] # just a 2-d list
vmtestInv = vmsym*vmsymInv
for i in range(k):
for j in range(k):
#print i, j, vmtrx2[i,j] % mod
vmtest[i][j] = vmtestInv[i,j] % mod
print "test vmk*vkinv % mod \n:", vmtest
return vmnpInv
if __name__ == '__main__':
#p = 271
p =
115792089210356248762697446949407573530086143415290314195533631308867097853951
vm = numpy.array([[1,1,1,1], [1, 2, 4, 8], [1, 4, 16, 64], [1, 5, 25, 125]])
#vminv = modMatInv(vm, p)
vminv = matInvMod(vm, p)
打印 vminv
vmtestnp = vm.dot(vminv)%p # 测试 mtrx 反转
打印 vmtestnp
【讨论】：