I'd like to take the modular inverse of a matrix like [[1,2],[3,4]] mod 7 in Python. I've looked at numpy (which does matrix inversion but not modular matrix inversion) and I saw a few number theory packages online, but nothing that seems to do this relatively common procedure (at least, it seems relatively common to me).
By the way, the inverse of the above matrix is [[5,1],[5,3]] (mod 7). I'd like Python to do it for me though.
It can be calculated using Sage (www.sagemath.org) as
Although Sage is huge to install and you have to use the version of python that comes with it which is a pain.
A hackish trick which works when rounding errors aren't an issue:
A less hackish way is to actually implement gaussian elimination. Here's my code using Gaussian elimination, which I wrote for my own purposes (rounding errors were an issue for me). q is the modulus, which is not necessarily prime.
EDIT: as commenters point out, there are some cases this algorithm fails. It's slightly nontrivial to fix, and I don't have time nowadays. Back then it worked for random matrices in my case (the moduli were products of large primes). Basically, the first non-zero entry might not be relatively prime to the modulus. The prime case is easy since you can search for a different row and swap. In the non-prime case, I think it could be that all leading entries aren't relatively prime so you have to combine them
This little piece of code seems to do it: link
Note the comment below for a little improvement. Seems to do the correct linear algebra as far as I can see. I have never found any option in regular packages so probably taking a code snippet from the web (there are a lot more available) is the easiest approach.
Unfortunately numpy does not have modular arithmetic implementations. You can always code up the proposed algorithm using row reduction or determinants as demonstrated here. A modular inverse seems to be quite useful for cryptography.
Okay...for those who care, I solved my own problem. It took me a while, but I think this works. It's probably not the most elegant, and should include some more error handling, but it works:
'sympy' package Matrix class function 'sqMatrix.inv_mod(mod)' computes modulo matrix inverse for small and arbitrarily large modulus. By combining sympy with numpy, it becomes easy to compute modulo inverse of 2-D numpy arrays (see the code snippet below):
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) print vminv vmtestnp = vm.dot(vminv)%p # test mtrx inversion print vmtestnp