Módulo de eliminación gaussiano p

Question

Estoy intentando reescribir un código en Java resolviendo un conjunto de ecuaciones lineales haciendo eliminación gaussiana en flotadores, para trabajar en ecuaciones módulo primo.El problema es que no funciona y no puedo entender qué pasa.Parece funcionar con conjuntos pequeños de ecuaciones, pero no con conjuntos grandes, lo que dificulta la depuración.

Mi algoritmo toma la primera fila, la normaliza encontrando el inverso del primer elemento y multiplica cada elemento de la fila por este inverso.Luego resta esta fila de las otras filas suficientes veces para que su primer elemento sea cero.En la siguiente iteración, pasa a la siguiente fila y realiza el mismo procedimiento, hasta que el elemento pivotante de la fila i esté en la columna i.Por último, resta cada fila de las filas anteriores para crear solo un elemento distinto de cero en cada columna (excepto la última).(A partir de ahora uso dobles, lo cual no es necesario, pero esto no debería ser un problema).Aquí está mi código:

// Transforms A so that the leftmost square matrix has at most one 1 per row,
    // and no other nonzero elements.
    // O(n^3)
 public static void gauss(int[][] A, int num_columns) {
        int n = A.length;
        int m = A[0].length;

        for (int i = 0; i < num_columns; i++) {
            // Finding row with nonzero element at column i, swap this to row i
            for(int k = i; k < num_columns; k++){
                if(A[k][i] != 0){
                    int t[] = A[i];
                    A[i] = A[k];
                    A[k] = t;
                }
            }
            // Normalize the i-th row.
            int inverse = (int)inverse((long)A[i][i], prime);
            for (int k = i ; k < m; k++) A[i][k] = (A[i][k]*inverse) % prime;

            // Combine the i-th row with the following rows.
            for (int j = 0; j < n; j++) {
                if(j == i) continue;
                int c = A[j][i];
                A[j][i] = 0;
                for (int k = i + 1; k < m; k++){
                    A[j][k] = (A[j][k] - c * A[i][k] + c * prime) % prime;
                }
            }
        }
    }

    public static void gauss(int[][] A) {
        gauss(A, Math.min(A.length, A[0].length));
    }
    public static long gcd(long a, long b){
        if(a < b){
            long temp = a;
            a = b;
            b = temp;
        }
        if(b == 0) return a;
        return gcd(b, a % b);
     }
    public static Pair ext_euclid(long a, long b){
        if(a < b){
            Pair r = ext_euclid(b,a);
            return new Pair(r.second, r.first);
        }
        if(b == 0) return new Pair(1, 0);
        long q = a / b;
        long rem = a - b * q;
        Pair r = ext_euclid(b, rem);
        Pair ret = new Pair(r.second, r.first - q * r.second);
        return ret;
    }

    public static  long inverse(long num, long modulo){
        num = num%modulo;
        Pair p = ext_euclid(num, modulo);
        long ret = p.first;
        if(ret < 0) return (modulo + ret) % modulo;
        return ret % modulo;
    }

    static class Pair{
        public long first;
        public long second;
        public Pair(long frst, long scnd){
            first = frst;
            second = scnd;
        }
    }

Esto funciona en pequeños ejemplos (mod 29):

matrix = {{1.0, 1.0, 1.0, 1.0}, {1.0, 2.0, 1.0, 2.0},{1.0, 0.0, 0.0‚ 3.0}};
answer= {{1.0, 0.0, 0.0, 0.0},{0.0, 1.0, 0.0, 1.0}, {0.0, 0.0, 1.0, 0.0}};

Lo cual es correcto (primera variable = 0, segunda = 1.0, tercera = 0), como lo puede verificar WolframAlpha para 0*k^0 + 1*k^1 + 0*k^2 para k = 1..3.

Para este ejemplo, teniendo 10 variables y las ecuaciones a*k^0 + b*k^1 + c*k^2...(mod 29) para k = 1..11, tengo esta matriz:

1 1 1 1 1 1 1 1 1 1 1 8 
1 2 4 8 16 3 6 12 24 19 9 5 
1 3 9 27 23 11 4 12 7 21 5 12 
1 4 16 6 24 9 7 28 25 13 23 12 
1 5 25 9 16 22 23 28 24 4 20 15 
1 6 7 13 20 4 24 28 23 22 16 0 
1 7 20 24 23 16 25 1 7 20 24 5 
1 8 6 19 7 27 13 17 20 15 4 1 
1 9 23 4 7 5 16 28 20 6 7 18 
1 10 13 14 24 8 22 17 25 18 7 20 
1 11 5 26 25 14 9 12 16 7 7 8  

Usando mi algoritmo obtengo la respuesta:

1 0 0 0 0 0 0 0 0 0 0 18 
0 1 0 0 0 0 0 0 0 0 0 8 
0 0 1 0 0 0 0 0 0 0 0 15 
0 0 0 1 0 0 0 0 0 0 0 11 
0 0 0 0 1 0 0 0 0 0 0 28 
0 0 0 0 0 1 0 0 0 0 0 27 
0 0 0 0 0 0 1 0 0 0 0 7 
0 0 0 0 0 0 0 1 0 0 0 21 
0 0 0 0 0 0 0 0 1 0 0 9 
0 0 0 0 0 0 0 0 0 1 0 24 
0 0 0 0 0 0 0 0 0 0 1 14  

¡Pero esto está mal!(se puede comprobar con WolframAlpha).La respuesta correcta debería ser (a b c...) = (8 13 9 13 4 27 18 10 12 24 15).

¿Alguien puede detectar mi error?¿O he entendido mal cómo hacer Gauss mod p?

Answer 1

No parece encontrar una fila con un valor distinto de cero en la iésima posición durante la primera etapa.Dado eso, y el hecho de que no sé cómo tu inverse función funciona, es posible que tenga o no un problema allí.

No entiendo por qué buscas "pivotes" en la segunda etapa.Sabes dónde están desde la primera etapa.En realidad, no veo por qué tienes una segunda etapa.En su lugar, haz toda la eliminación en la primera etapa.Esto aclarará enormemente su código.

No entiendo por qué estás usando doubles en su matriz.Tampoco entiendo por qué estás usando Math.abs por todo el lugar;Las comparaciones de igualdad son muy apropiadas aquí.

Finalmente, estás resolviendo un sistema Vandermonde.Si esa es su aplicación, no solo un caso de prueba, probablemente debería usar la interpolación de Lagrange.