# How to generate a matrix group?

I have three $7\times 7$ matrices (with real entries, lots of zeros) and I'd like to check if they generate a finite group (or, more precisely, if the group they generate is of precise order). Would it be possible to do with Mathematica?

These are the matrices. All of order 3. I expect this group to be of order 6048 or 12096.

gMatrix =
SparseArray[{{1, 3} -> 1, {2, 1} -> 1, {3, 2} ->
1, {4, 4} -> -(1/2), {4, 5} -> -(1/2), {4, 6} -> -(1/2), {4,
7} -> -(1/2), {5, 4} -> 1/2, {5, 5} -> -(1/2), {5, 6} ->
1/2, {5, 7} -> -(1/2), {6, 4} ->
1/2, {6, 5} -> -(1/2), {6, 6} -> -(1/2), {6, 7} -> 1/2, {7, 4} ->
1/2, {7, 5} ->
1/2, {7, 6} -> -(1/2), {7, 7} -> -(1/2), {_, _} -> 0}, 7];

hMatrix =
SparseArray[{{1, 5} ->
1, {2, 2} -> -(1/2), {2, 3} -> -(1/2), {2, 6} -> 1/2, {2, 7} ->
1/2, {3, 2} -> 1/2, {3, 3} -> -(1/2), {3, 6} -> -(1/2), {3, 7} ->
1/2, {4, 1} -> 1, {5, 4} -> 1, {6, 2} -> -(1/2), {6, 3} ->
1/2, {6, 6} -> -(1/2), {6, 7} ->
1/2, {7, 2} -> -(1/2), {7, 3} -> -(1/2), {7, 6} -> -(1/2), {7,
7} -> -(1/2), {_, _} -> 0}, 7];

kMatrix =
SparseArray[{{1, 4} -> -1, {2, 2} -> -(1/2), {2, 3} -> -(1/2), {2,
6} -> -(1/2), {2, 7} -> 1/2, {3, 2} ->
1/2, {3, 3} -> -(1/2), {3, 6} -> -(1/2), {3, 7} -> -(1/2), {4,
5} -> -1, {5, 1} -> 1, {6, 2} -> 1/2, {6, 3} ->
1/2, {6, 6} -> -(1/2), {6, 7} -> 1/2, {7, 2} -> -(1/2), {7, 3} ->
1/2, {7, 6} -> -(1/2), {7, 7} -> -(1/2), {_, _} -> 0}, 7];

-
Related: What is the command to find function invariant? where I construct a multiplication table. You need to generate such a table for your matrices to verify the closure of the set, or to identify missing elements. You're going to have to show the matrices to get more help. Perhaps best to list them here after applying ArrayRules[...] to each of them and posting the output. –  Jens Mar 22 '13 at 16:36
Will it work with a group of four-digits order? –  mathdonk Mar 22 '13 at 16:38
Some closely related techniques are described at mathematica.stackexchange.com/questions/13835/…. Initially identifying some simple relations among the matrices could be exploited to (substantially) reduce the billions of matrix multiplications and comparisons that might otherwise be needed. In particular, you should at a minimum precompute the order of each matrix and establish some kind of nontrivial relation between each matrix and at least one other matrix. –  whuber Mar 22 '13 at 16:53
OK, so these are three orthogonal matrices with determinant 1. That means they can be generated by $n(n-1)/2 = 21$ parameters. What remains to do is to find the subgroup to which they belong. –  Jens Mar 22 '13 at 22:58
@Jens That's a great observation. It's unclear how to exploit it, though. After all, orthogonal groups have finite subgroups of arbitrarily large order. –  whuber Mar 23 '13 at 2:31

The group has 6048 elements. (Could it be isomorphic to $U_3(3)$?--see below.)

count = 0; (matrices = NestWhile[(Print[count++];
Union[#~Join~Flatten[Outer[Dot, {gMatrix, hMatrix, kMatrix}, #, 1], 1]]) &,
{IdentityMatrix[7]}, Length[#2] != Length[#1] &, 2, 99]) // Length // Timing


$\{2.2, 6048\}$

This code begins with the identity matrix (IdentityMatrix) and iteratively adjoins (Join) any new matrices (Union) created by left multiplication of all previous matrices (Outer, Dot) by the three generators defined in the question. (Flatten unravels the 2D array of matrices into a simple list of matrices.) It stops when no new matrices appear (the Length test). In the process, it prints out an iteration count to verify that the loop did not terminate because the calculation had simply gone on too long (limited to $99$ iterations here; only $13$ were needed). In addition to reporting the total time ($2.2$ seconds), the output consists of the number of unique matrices ($6048$) and a list of them all (matrices), useful for further processing. For instance, we can compute and visualize the orders of the group elements:

matrixOrder[a_, maxIter_: 9999] := Module[{n = 1, i = IdentityMatrix[Length[a]]},
NestWhile[(n++; a.#) &, a, ! And @@ (PossibleZeroQ /@ Flatten[# - i]) &, 1, maxIter]; n];
Histogram[matrixOrder /@ matrices, {1/2, 25/2, 1}]


$U_3(3)$, a sporadic simple group of Lie type modeled on $G_2$, has a seven-dimensional integral representation generated by matrices of orders $2$ and $6$,

a = {{-1, 0, 0, 0, 0, 0, 0}, {0, -1, 0, 0, 0, 0, 0}, {0, 1, 0, 0, 0, 1, 0}, {-1, 0, 0, 0, 0, 0, 1}, {0, 0, 0, 0, 1, 0, 0}, {0, 1, 1, 0, 0, 0, 0}, {-1, 0, 0, 1, 0, 0, 0}};
b = {{0, -1, 0, 0, 0, -1, 0}, {0, 1, 1, 0, 0, 0, 0}, {0, -1, 0, 0, 0, 0, 0}, {1, 0, 0, 0, 0, 0, 0}, {0, 0, 0, 0, 0, 0, -1}, {0, 0, 0, -1, 0, 0, 1}, {0, 0, 0, 0, 1, 0, 1}};


Using the preceding code, it takes a quarter second to generate all $6048$ elements from these two matrices. The histogram of orders of its elements is identical to the previous figure. An exhaustive search among our group for such matrices (c. $17$ seconds) finds at least one such pair of matrices that satisfy all the required relations, proving this group is isomorphic to $U_3(3)$.

(* By construction, the last element of matrices is the identity. *)
i7 = IdentityMatrix[7];
order2 = Select[Most@matrices, MatrixPower[#, 2] == i7 &];
order6 = Select[Most@matrices, MatrixPower[#, 6] == i7 &&
MatrixPower[#, 3] != i7 && MatrixPower[#, 2] != i7 &];
commutator[a_, b_] := a.b.Inverse[a].Inverse[b];
(* generate assumes a has order 2 and b has order 6. *)
generate[a_, b_] := MatrixPower[a.b, 7] == i7  &&
commutator[a, MatrixPower[a.b.b, 3]] == i7  &&
MatrixPower[b, 3].MatrixPower[commutator[b.b, a.MatrixPower[b, 3].a], 2] == i7;
(* Does any pair at all qualify as generators? *)
Or @@ Flatten@Outer[generate, order2, order6, 1]


$\text{True}$

(Specifically, among many possibilities, $a'=(ghk)^2$ (which is integral) and $b'=ghg$ generate $U_3(3)$.)

-
That's better than using Tuples - the scaling is dramatically improved. –  Jens Mar 23 '13 at 6:21

A brute-force approach would be to define the set of matrices and form all products of them:

FixedPoint[Union[#, Dot @@@ Tuples[#, 2]] &, set]


The Union sorts the entries so that duplicates are eliminated. The length of the resulting list is the order of the group.

Edit to explain the code

In the command FixedPoint, the first argument is a function that will be applied repeatedly to the second argument (set), until the result no longer changes.

So I need a function, and that's provided by the construct Union[#,...]& where the # stands for the Slot (or function parameter) and the & says that it's a Function.

Inside the function is the instruction to form the Union of two sets, namely the original one (from the previous iteration), and the outcome of forming the Dot products between all Tuples of matrices in the previous list. To get the latter, the Dot command is applied to all members of the list generated by Tuples; this is what the command @@@ does.

Edit to avoid crashes

Since the number of tuples explodes with the number of elements, brute force approaches like this can easily cause the computation to hang. Either try MemoryConstrained as nikie mentioned, or alternatively add a termination condition to FixedPoint that stops it before the explosion happens. This could be done by adding the option SameTest -> ((Length[#] == Length[#2]) || Length[#] > 50 &) to FixedPoint. Here, the length cutoff 50 is meant to limit the number of Tuples to a manageable size.

-
What do "#" and "@" stand for here? It's probably smart but I don't understand how it would work... –  mathdonk Mar 22 '13 at 17:11
That this works can be established by defining set[n_] := {RotationMatrix[2 \[Pi]/n], {{0, 1}, {1, 0}}} and applying it to, say, set[3] or even set[12] (a group of 24 2D matrices). But its fundamental flaw will become apparent when you change this to the equally innocent set[25] (a 50-element group): your machine will crash due to the lack of exabytes of RAM it needs to complete the task. –  whuber Mar 22 '13 at 17:12
@whuber The question doesn't say that this is going to be a large group - the order obviously can't be deduced from the matrix dimension. So this is all you can do in the absence of more information. If there's reason to believe that the group order is in the four digits, the structure of the matrices would have to be provided. I'm not telepathic today. –  Jens Mar 22 '13 at 17:21
The main reason is to avoid having to reboot one's system :-). –  whuber Mar 22 '13 at 17:35
@whuber Brute force methods without crashes aren't true brute force... But I'll add some safety remarks, I think. –  Jens Mar 22 '13 at 18:11