restart:
with(combinat):
structure:=proc(M)
#option trace;
local i,j,A,lambda,l,k,n,W,m,s,B,T,L,flag,m1;
    n:=nops(M);
    lambda:=Array(1..n,1..n,(i,j)->gcd(M[i],M[j]));
    W:=Array(1..n,1..n,1..n,(i,j,k)->[op({op(M[k]/gcd(M[k],lambda[i,j])*[seq(t,t=0..M[k]-1)] mod M[k])})]);
    T := cartprod([seq(seq(seq(W[i, j, k], i = 1 .. n), j = 1 .. n), k = 1 .. n)]):
    A:=Array();
    m:=1:
    while not T[finished] do
        A(m):= T[nextvalue](); 
        m := m+1;
    od:
    B:=Array();
    m:=1:
    L:=[seq(seq(seq([i, j, k], i = 1 .. n), j = 1 .. n), k = 1 .. n)]:
    for s from 1 to ArrayNumElems(A) do
        B(m):=Array(1..n,1..n,1..n);
        for l from 1 to nops(L) do 
            B(m)[op(L[l])]:=A[s][l];
        od:
        m:=m+1:
    od:
    A:=Array();
    m:=1:
    for s from 1 to ArrayNumElems(B) do 
        flag:=true:
        for i from 1 to n while flag do
            for j from 1 to n while flag do
                for k from 1 to n while flag do
                    for l from 1 to n while flag do 
                        if add(B[s][i,j,alpha]*B[s][alpha,l,k]-B[s][j,l,alpha]*B[s][i,alpha,k],alpha=1..n) mod M[k]<>0 then
                            flag:=false:  
                        fi;
                    od:
                od;
            od:
        od:
        if flag then
            A(m):=B[s];
            m:=m+1:
        fi:
    od:
RETURN(A);
end:
isomorphism:=proc(M,W,Z,B)
#option trace;
local i,j,P,PP,p,lambda,l,k,n,m,s,T,L,flag,A,alpha,beta,b,m1;
    n:=nops(M):
    for p from 1 to ArrayNumElems(B)  do
        flag:=false:
        for k from 1 to n while not flag do
            for i from 1 to n while not flag do
                for j from 1 to n while not flag do 
                    if (add(add(B[p][i,alpha]*B[p][j,beta]*W[alpha,beta,k],alpha=1..n),beta=1..n)-add(Z[i,j,l]*B[p][l,k],l=1..n)) mod M[k]<>0 then
                        flag:=true:
                    fi:
                od:
            od:
        od:
        if not flag then #print(B[p]);
            RETURN(true);
        fi:
    od:      
RETURN(false);
end:
classification:=proc(M)
#option trace;
local BB,AA,B,i,firsttime,firstbytes,secondtime,secondbytes,j,P,PP,p,lambda,l,k,n,m,s,T,L,flag,A,alpha,beta,b,m1,flagg,FLAG;
    firsttime,firstbytes:=kernelopts(cputime,bytesused);
    A:=structure(M):
    B:=Array();
    B(1):=A[2];print(seq(seq(seq(A[2][ii,j,k], ii = 1 .. nops(M)),j = 1 .. nops(M)),k = 1 .. nops(M)));    
    n:=nops(M);
    m:=max(op(M)):
    m1:=mul(M[i],i=1..n);m1:=m:
    lambda:=Array(1..n,1..n,(i,j)->gcd(M[i],M[j]));
    AA:=Array(1..n,1..n,(j,k)->[op({op(M[k]/gcd(M[k],M[j])*[seq(t,t=0..M[k]-1)] mod M[k])})]);
    T := cartprod([seq(seq(AA[j, k],j = 1 .. n), k = 1 .. n)]):
    P:=Array();
    m:=1:
    while not T[finished] do
        P(m):= T[nextvalue](); 
        m := m+1;
    od:
    L:=[seq(seq([j, k], j = 1 .. n), k = 1 .. n)]:
    BB:=Array();
    m:=1:
    for s from 1 to ArrayNumElems(P) do
        BB(m):=Array(1..n,1..n);
        for l from 1 to nops(L) do 
            BB(m)[op(L[l])]:=P[s][l];
        od:
        m:=m+1:
    od:
    T := cartprod([seq([seq(i,i=0..m1-1)],j=1..n)]);
    P:=Array(1..m1^n,i->T[nextvalue]());
    m:=0:
    PP:=Array();
    for b from 1 to ArrayNumElems(BB) do 
        FLAG:=true:
        for p from 1 to ArrayNumElems(P) while FLAG do
            flag:=true: 
            for k from 1 to n while flag do 
                if not (add(P[p][j]*BB[b][j,k],j=1..n) mod M[k]=0) then
                    flag:=false:
                fi:
            od: 
            flagg:=true:
            for i from 1 to n while flagg do
                if not (P[p][i] mod M[i]=0) then
                    flagg:=false;
                fi;
            od;
            if not ((not flag) or flagg) then
                FLAG:=false:
            fi:  
        od:  
        if FLAG then
            PP(m+1):=BB[b]:
            m:=m+1:
        fi:
    od:
    BB:=copy(PP):
    for i from 3 to ArrayNumElems(A) do
        flag:=false:
        for j from 1 to ArrayNumElems(B) while not flag do   
            flag:=isomorphism(M,B[j],A[i],BB); 
        od:
        if not flag then
            B(ArrayNumElems(B)+1):=A[i];print(seq(seq(seq(A[i][ii,j,k], k = 1 .. nops(M)),j = 1 .. nops(M)),ii = 1 .. nops(M)));
        fi;
    od:
    B(ArrayNumElems(B)+1):=A[1];print(seq(seq(seq(A[1][ii,j,k], k = 1 .. nops(M)),j = 1 .. nops(M)),ii = 1 .. nops(M)));
    secondtime,secondbytes:=kernelopts(cputime,bytesused);
    printf("%-1s %1s %1s %1s        : %3a %3a\n",The, cpu, time, is,secondtime-firsttime,(sec)):
    printf("%-1s %1s %1s        : %3a %3a\n",The,used,memory,secondbytes-firstbytes,(bytes)):
    printf("%-1s %1s %-1s             : %3g\n",Num,of,Rings,ArrayNumElems(B)):
RETURN();
end:
M:=[5^2]:
n:=nops(M):
A:=classification(M):

M:=[5,5]:
n:=nops(M):
A:=classification(M):

#for s from 1 to ArrayNumElems(A) do 
#    seq(seq(seq(A[s][i,j,k], i = 1 .. n),j = 1 .. n),k = 1 .. n); 
#od;