Re: Find an invertible integer matrix that satisfies given conditions.

[Bill Allombert <Bill.Allombert@math.u-bordeaux.fr>]

On Wed, Apr 12, 2023 at 09:24:45AM +0200, Aurel Page wrote:
> Dear Zhao,
> 
> Your problem is equivalent to finding an isomorphism between two
> Z[t]-lattices. This is a nontrivial task, but there are algorithms for it,
> if you want an efficient solution. You can look at the work of Tommy
> Hofmann, he has several papers on this topic. If you don't care too much
> about efficiency, I think a solution that is reasonable to implement is to
> first parametrise the set of solutions of the linear equation (with mathnf
> or matsolvemod), and then backtrack the coefficients, for instance column by
> column, using the fact that for k columns to admit an extension to an nxn
> invertible matrix over Z, a necessary and sufficient condition is that their
> HNF must have k pivots equal to 1. There are also many quick necessary
> conditions, for instance M1 and M2 must be conjugate over Q and mod p for
> every p.

For 2x2 matrices, the program in attachment solves it.

M=[1,2;3,4];
N=[5,-2;-1,0];
P=isom(M,N)
M*P-P*N
matdet(P)

? P=isom(M,N)
%3 = [1,0;2,-1]
? M*P-P*N
%4 = [0,0;0,0]
? matdet(P)
%5 = -1

Over the rationals, the problem can be solved with matfrobenius(,1). Sometime
we are lucky and get an integral matrix.

M1 = [0, 1, 0, 0; -4, 0, 0, 0; 0, 0, 0, 1; 0, 0, -4, 0];
M2 = [0, 1, 4, 0; -4, 0, 0, -4; 0, 0, 0, 1; 0, 0, -4, 0];
[F1,P1]=matfrobenius(M1,2)
[F2,P2]=matfrobenius(M2,2)
P = P1^-1*P2
matdet(P)
M1*P-P*M2

? P = P1^-1*P2
%5 = [1,0,0,1;0,1,0,0;0,0,1,0;0,0,0,1]
? matdet(P)
%6 = 1
? M1*P-P*M2
%7 = [0,0,0,0;0,0,0,0;0,0,0,0;0,0,0,0]

So here P is the solution.

In general:
1) compute the matrix of the linear map P -> M1*P-P*M2 
2) compute a basis of its integral kernel (K_1,...,K_r) with matkerint
3) solve the Diophantine equation matdet(sum x_i*K_i) = ±1 (with unknown x_1,...,x_r).
4) Return P = sum x_i*K_i

The difficult part is obviously 3).
For n=2, the Diophantine equation is a binary quadratic form, so can be solved by PARI.
Hofmann shows that solving the equation can be reduced to solve a set of
norm equations in orders of number fields.

Cheers,
Bill

isom(M,N)=
{
  my(R=
  [M[1,1]-N[1,1],M[2,1],-N[1,2],0;
   M[1,2],M[2,2]-N[1,1],0,-N[1,2];
   -N[2,1], 0, M[1,1]-N[2,2],M[2,1];
   0,-N[2,1],M[1,2],M[2,2]-N[2,2]]);
  my(K=matkerint(R),Kx=K*[x,1]~);
  if(#K!=2,error(K));
  my(P=matdet(Mat([Kx[1..2],Kx[3..4]])));
  my(S);
  if(P==1 || P==-1,
     S = [0,1]~
   , poldegree(P)==1,
     my([u,v,d]=bezout(polcoeff(P,1),polcoeff(P,0)));
     if(d!=1,return(0));
     S=[u,v]~;
   , my(Q=Qfb(Vec(P)));
     my(s1=qfbsolve(Q,1),s2=qfbsolve(Q,-1));
     S = if(s1,s1~,s2,s2~,return(0)));
  S = K*S;
  Mat([S[1..2],S[3..4]])~;
}

test()=
{
  M=[1,2;3,4];
  N=[5,-2;-1,0];
  P=isom(M,N);
  M*P-P*N
}