/* This is a Magma program to generate all maximal graphs (up to isomorphism) of rank n. 
As input, it requires the data set of all connected graphs of order n in graph6 format which is available at 
Brendan McKay's page: http://users.cecs.anu.edu.au/~bdm/data/. */

F := OpenGraphFile("graph6c.g6", 0, 0);
n:=6;
count := 1;
Maximals:={ };
M:=[* *];
more, g := NextGraph(F);
while more do
     B:=AdjacencyMatrix(g);
     rk:=Rank(B);
     if rk eq n then
        A:=ChangeRing(B,RationalField());
        A:=A^-1;
        Vg:=Vertices(g);
        B0:=Subsets({Vg[i] : i in [1..n]}) diff { Neighbors(i) : i in Vg} diff {{}};
        XX:={ };
        for a in B0 do
          x:={};
          for i in [1..n] do 
            if i in a then x:=x join {i}; end if;
          end for;
          mx:=&+[A[i,j] : i in x, j in x];
          if mx eq 0 then XX:=XX join {x}; end if;
        end for;   
        nXX:=#XX;
        XX:=SetToIndexedSet(XX); 
        if nXX eq 0 then Maximals:=Maximals join {g}; end if;
        if nXX gt 0 then 
          edge:={ }; edge1:={};
          for i in [1..n] do for j in [i+1..n] do
            if B[i,j] eq 1 then edge1:= edge1 join {{i,j}}; end if; 
          end for; end for;
          for i in [1..nXX] do
            edge1:=edge1 join {{n+i,j} : j in XX[i]};             
            for j in [i+1..nXX] do
               mij:=&+[A[r,s] : r in XX[i], s in XX[j]];
               if mij eq 0 or mij eq 1 then edge:=edge join { {n+i,n+j} }; end if;
               if mij eq 1 then edge1:=edge1 join { {n+i,n+j} }; end if;  
             end for;
           end for;
           Compat:= Graph < {n+1..n+nXX} | edge>; 
           Compat1:= Graph < n+nXX | edge1>;
           V:=Vertices(Compat); V1:=Vertices(Compat1);
           AC:=AllCliques(Compat);
           cc:=0; 
           for C in AC do 
                h:=sub<Compat1 | {V1[i] : i in [1..n]} join {V1[n+Position(V,i)] : i in C} >;
                cc+:=1; M[cc]:=CanonicalGraph(h);
           end for;
           Maximals:=Maximals join {M[i] : i in [1..cc]}; 
       end if;
     end if;
   count +:= 1;
   if count mod 500 eq 0 then print count,#Maximals; end if;
   more, g := NextGraph(F);
end while;
print count,#Maximals; 
X:={* Order(h) : h in Maximals *};
print(X);
for h in Maximals do
 mm:=Order(h);
 M:=AdjacencyMatrix(h);
 for i in [1..mm-1] do
   for j in [i+1..mm] do
      printf "%o", M[i,j];
   end for;
   printf "\n";
 end for;
end for;