¿La manera más fácil de realizar la inversión de matriz modular con Python?

Question

Me gustaría tomar la inversa modular de una matriz como [[1,2], [3,4]] mod 7 en Python. He observado numpy (que hace inversión de matriz pero no inversión de matriz modular) y vi algunos paquetes de teoría de números en línea, pero nada que parezca hacer este procedimiento relativamente común (al menos, me parece relativamente común).

Por cierto, la inversa de la matriz anterior es [[5,1], [5,3]] (mod 7). Aunque me gustaría que Python lo haga por mí.

Answer 1

Esta pequeña pieza de código parece hacerlo: link

Nota del comentario más abajo para una pequeña mejora. Parece que hace el álgebra lineal correcta hasta donde puedo ver. Nunca he encontrado ninguna opción en los paquetes regulares, por lo que probablemente tomar un fragmento de código de la web (hay muchos más disponibles) es el enfoque más fácil.

Answer 2

Desafortunadamente numpy no tiene implementaciones aritméticas modulares. Siempre puede codificar el algoritmo propuesto usando la reducción de fila o los determinantes como se demuestra here. Un inverso modular parece ser bastante útil para la criptografía.

Answer 3

Bien ... para los que les importa, resolví mi problema. Me tomó un tiempo, pero creo que esto funciona. Es probablemente no es el más elegante, y debe incluir un poco más de control de errores, pero funciona:

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 

Answer 4

Un truco hacker que trabaja cuando los errores de redondeo no son un problema:

• encontrar la inversa regulares (pueden tener entradas no enteros), y el determinante (un número entero), ambos implementado en numpy
• multiplicar la inversa por el determinante, y la ronda de números enteros (hacky)
• ahora multiplica todo por inverso multiplicativo del determinante (modulo tu modulu s, código de abajo)
• qué entrywise mod por su módulo de

Una manera menos hacker es para aplicar en la práctica la eliminación de Gauss. Aquí está mi código usando la eliminación de Gauss, que escribí para mis propios fines (los errores de redondeo fueron un problema para mí). q es el módulo, que no es necesariamente primo.

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 

EDITAR: como comentan, en algunos casos este algoritmo falla. No es trivial repararlo, y no tengo tiempo hoy en día. En aquel entonces funcionaba para matrices aleatorias en mi caso (los módulos eran productos de primos grandes). Básicamente, la primera entrada distinta de cero podría no ser relativamente primo para el módulo. El caso principal es fácil ya que puede buscar una fila diferente y cambiar. En el caso no primo, creo que podría ser que todo entradas principales no son primos entre sí lo que hay que combinarlos

Answer 5

Se puede calcular usando Sage (www.sagemath.org) como

Matriz (IntegerModRing (7), [[1, 2], [3,4]]).inversa()

Aunque Sage es enorme para instalar y usted tiene que utilizar la versión de Python que viene con él, que es un dolor.