1

I am trying to generate a random binary matrix and its inverse mod q where q is a power of 2. Sometimes when the determinant of my matrix is invertible modulo q (so the matrix over Z_q is invertible), I am getting the error "InvMod:inverse undefined Aborted (core dumped)" and other times the inverse is computed. What am I doing incorrectly?

#include <iostream>
//NTL files
#include <NTL/ZZ_p.h>
#include <NTL/vec_vec_ZZ_p.h>
#include <NTL/LLL.h>
#include <NTL/matrix.h>
#include <NTL/vector.h>
#include <NTL/tools.h>
#include <NTL/ZZ.h>
#include <NTL/vec_vec_ZZ.h>
using namespace std;
using namespace NTL;

int main(){//task generate a random matrix S with 0/1 entries stored as a ZZ_p matrix, then generate a random, invertible S 
int nn = 8;
ZZ n = ZZ(nn);
ZZ N = ZZ(0);
ZZ q; power2(q, 4);
ZZ_p::init(q);

mat_ZZ S; S.SetDims(nn,nn); 
for(int i = 0; i<nn; i++){
    for(int j = 0; j<nn; j++){
        S[i][j] = RandomBits_ZZ(1);     
    }
}
mat_ZZ_p S1; S1.SetDims(nn,nn);//copy to ZZ_P
mat_ZZ_p R; R.SetDims(nn,nn);//will set to inverse if 

cout<<"The random matrix is S = "<<endl; //print S
for(int i = 0; i<nn; i++){
    for(int j=0; j<n;j++){
        cout<<S[i][j]<<", ";
    } cout<<endl;
}

ZZ d; determinant(d,S); ZZ_p d1; conv(d1, d % q);
if(GCD(q,d) == 1){//convert to mod q datatype
    for(int i = 0; i<nn; i++){
        for(int j = 0; j<nn; j++){
            conv(S1[i][j], S[i][j]);        
        }
    }
    //let's invert the matrix and print it!
    cout<<"The random matrix is R = "<<endl; //print R
    R = inv(S1); //mul(R,R,S1);
    for(int i = 0; i<nn; i++){
        for(int j=0; j<n;j++){
            cout<<R[i][j]<<", ";
        } cout<<endl;
    }
}

cout<<endl<<"det of S is "<<d<<" and this mod q is "<<d1<<endl;
cout<<"Our modulus is "<< q <<endl;

return 0;
}
user11235
  • 11
  • 1

1 Answers1

1

If the determinant is invertable mod q this only means that there exists an inverse matrix. But the algorithm that computes this matrix can still come to a point where it would need to calculate the inverse of an element that don't has one.

You don't have this problem if q is prime.

By the way, here is a simplified version of your code.

#include <iostream>
//NTL files
#include <NTL/mat_ZZ_p.h>

using namespace std;
using namespace NTL;

int main()
{//task generate a random matrix S with 0/1 entries stored as a ZZ_p matrix, then generate a random, invertible S
    int nn = 8;
    ZZ q;
    power2(q, 4);
    ZZ_p::init(q);

    mat_ZZ_p S;
    S.SetDims(nn, nn);

    for(int i = 0; i < nn; i++)
    {
        for(int j = 0; j < nn; j++)
        {
            S[i][j] = conv<ZZ_p>(RandomBits_ZZ(1));
        }
    }

    mat_ZZ_p R;
    R.SetDims(nn, nn);//will set to inverse if

    cout << "The random matrix is S = " << endl << S;

    ZZ_p d;
    determinant(d, S);

    cout << endl << "det(S) = " << d << endl;
    cout << "q = " << q << endl;

    if(GCD(conv<ZZ>(d), q) == 1)
    {
        // let's invert the matrix and print it!
        R = inv(S);
        cout << "The random matrix is R = " << R << endl;
    }

    return 0;
}
AbcAeffchen
  • 14,400
  • 15
  • 47
  • 66
  • So the error goes away when we switch to a prime q. Is there an NTL algorithm that works for all q? I find it odd that say Maple has a function that computes the inverse mod n but NTL does not. – user11235 Jan 11 '18 at 22:38
  • @user11235 I'm not sure if NTL implements an algorithm that avoids this. I think the idea is that this function is only interesting if you want to work within a ring of matrices over a finite field. So maybe Shoup has just ignored the case for q not prime. – AbcAeffchen Jan 12 '18 at 00:43
  • Matrix inversion does work modulo powers of primes with zz_p's instead of ZZ_p's, I recently found out, so my problem is fixed. For arbitrary q, I think one would have to go through adjoints (the matrix with all the determinants of the minors). This wouldn't be too bad to implement and would have a complexity of about O(n^5). – user11235 Jan 14 '18 at 01:00
  • @user11235 This is interesting because the documentation only says that `zz_p` is the same as `ZZ_p` but p is single-precision. – AbcAeffchen Jan 14 '18 at 09:12
  • Here is the documentation I am viewing. It has a "relaxed" version of inversion and determinant which allow for the prime powers. [mat_lzz_p.cpp](http://www.shoup.net/ntl/doc/mat_lzz_p.cpp.html) PS-thank you for simplifying the code. all the includes was pretty messy. – user11235 Jan 14 '18 at 17:58