Tutorial: Group theory with Susyno

0. Getting started

Susyno is a Mathematica package which can make various calculations related to Lie groups and the permutation group $S_n$ (even though the main aim of the program is a different one, a substantial part of the code is group theory related). This page shows some examples.

To get started, the program needs to be installed: it can be download here . Decompressing the downloaded zip file will generate a folder Susyno. This folder should be placed in a directory visible to Mathematica (typing $Path in Mathematica shows a list of the possible locations). A good choice is to place it (the whole Susyno folder and not just its contents!) in

(Mathematica base directory)/AddOns/Applications
Linux, Mac OS
(Mathematica base directory)\AddOns\Applications

Once this is done the package is installed. To load Susyno, type in the front end

The program's reference is

"Renato M. Fonseca, Calculating the renormalisation group equations of a
SUSY model with Susyno
, Computer Physics Communications 183 (2012) 2298"

For help, there is an easy-to-use built-in documentation for each of the functions described in the following text (and other ones as well); it becomes accessible from within Mathematica once the package is installed. For questions, comments or bug reports, please contact me at


1. Working with representations of a Lie group

To indicate a group, just use it's name:

U1, SU2, SU3, SU5, SO10, E6, E8, G2, ...
If a group is the product of $U(1)$'s and/or simple groups, a list of the factor groups should be given. For example, the groups $SU(3)\times SU(2) \times U(1)$, $SU(5)\times SU(5)$ and $SO(10)\times U(1)$ are written as follows:

In fact, even groups with a single factor might have the brackets {} around. For example, one might write $SU(2)$ as


To indicate a representation of each factor group one needs to write its Dynkin coefficients (unless it is a $U(1)$, in which case the charge [i.e., one number] is enough). These are a list of $n$ non-negative integers, where $n$ is the group's rank. In turn, the group's rank corresponds to the maximal number of elements of the group's algebra which can be made simultaneously diagonal, and in practice this can be obtained with Susyno by typing

such as


This means that $SU(5)$'s rank is four, and so each of its representations can be uniquely identified by a list {n1,n2,n3,n4} of four non-negative integers. The user can call the RepName function to help identify each representation from the Dynkin coefficients. For example:


With the help of the function RepsUpToDimN which lists all the representations of a group up to some size, it is then very easy to get a table of representations. For example, the list of all $SU(3)$ representations up to size 30 can be obtained as follows:

su3Reps = RepsUpToDimN[SU3, 30];
Grid[Prepend[{#, RepName[SU3, #]}&/@ su3Reps,{"Dynkin coefficients", "Name"}],
Frame -> All, FrameStyle -> LightGray]

2. Dimension, Casimir, Dynkin index and triangular anomaly

Given a representation $R$, its dimension $d(R)$, Casimir $C(R)$ and Dynkin index $T(R)$ can be calculated with the functions DimR, Casimir and DynkinIndex:

DimR[SO10, {0, 0, 0, 0, 1}]
Casimir[SO10, {0, 0, 0, 0, 1}]
DynkinIndex[SO10, {0, 0, 0, 0, 1}]

(Recall that $C(R)$ and $T(R)$ are defined by the relations $C(R) Id=\sum_{a}T^{a}T^{a}$ and $T(R)\delta^{ab}=\textrm{Tr}\left(T^{a}T^{b}\right)$ where $T^a$ are the representation matrices of the algebra generators.)
These functions actually accept symbolic Dynkin coefficients. For example, the $SU(2)$ representation $\left\{ d-1\right\}$ has dimension $d$ and its Casimir is $\frac{1}{4} (d-1) (d+1)$:

DimR[SU2, {d - 1}]
Casimir[SU2, {d - 1}]

Gauge triangular anomalies are known to be associated to the quantity $\textrm{Tr}\left(\left\{ T^{a},T^{b}\right\} T^{c}\right)\equiv\kappa(R) d^{abc}$ where the symmetric tensor $d^{abc}$ can be taken to be fixed (for a given group) while $\kappa(R)$ depends on the representation. It should be noted that the only groups for which there might be anomalies are those with $SU(n>2)$ and/or with $U(1)$ factors, so at least one of the indices $a$, $b$ and $c$ must refer to one of these groups.
The triangular gauge anomalies associated to a representation can be calculated with the TriangularAnomalyValue function:

TriangularAnomalyValue[{SU3}, {{1, 0}}]
TriangularAnomalyValue[{U1}, {y}]

More complex cases, where the gauge group has more than one factor, are also possible. Consider for example the representation of left-handed quarks $Q$ in the Standard Model:

Q = {1/6 , {1}, {1, 0}};
TriangularAnomalyValue[{U1, SU2, SU3}, Q]

There is more than one anomaly in the output, and one case use the Verbose -> True option to identify each of them:

Q = {1/6 , {1}, {1, 0}};
TriangularAnomalyValue[{U1, SU2, SU3}, Q, Verbose -> True]

One can use this function, for example, to quickly show that all anomalies cancel between each Standard Model generation of fermions:

uc = {-2/3 , {0}, {0, 1}};
dc = {1/3 , {0}, {0, 1}};
Q = {1/6 , {1}, {1, 0}};
ec = {1, {0}, {0, 0}};
L = {-1/2 , {1}, {0, 0}};

fields = {uc, dc, Q, ec, L};
TriangularAnomalyValue[{U1, SU2, SU3}, #] & /@ fields

3. Products of representations

The decomposition of products of representations can be achieved with ReduceRepProduct:

tripletSU3 = {1, 0};
ReduceRepProduct[SU3, {tripletSU3, tripletSU3, tripletSU3}]

In this case, the output is saying that $\mathbf{3}\times\mathbf{3}\times\mathbf{3}$ in $SU(3)$ (note that {1,0} is the $\mathbf{3}$) contains the representation {3,0} once, {1,1} twice, and {0,0} once. If desired, the option UseName -> True can be used to convert these Dynkin coefficients into the names of the representations:

ReduceRepProduct[SU3, {tripletSU3, tripletSU3, tripletSU3}, UseName -> True]

The group can have various factors, for example:

uc = {-2/3 , {0}, {0, 1}};
Q = {1/6 , {1}, {1, 0}};
H = {1/2 , {1}, {0, 0}};
ReduceRepProduct[{U1, SU2, SU3}, {uc, Q, H}, UseName -> True]

Also, there is no limit to the number of representations being multiplied:

rep10 = {1, 0, 0, 0, 0};
ReduceRepProduct[SO10, {rep10, rep10, rep10, rep10, rep10}, UseName -> True]

Sometimes, the product of representation being calculated contains repeated representations, such as the $\mathbf{3}\times\mathbf{3}\times\mathbf{3}$ example above. In those cases, there is a permutation symmetry involved: for example, it is well known that the singlet in $\mathbf{3}\times\mathbf{3}\times\mathbf{3}$ is completely anti-symmetric. To calculate these symmetries, see section 6.

4. Explicit representation matrices

Consider the generators $T^a$ of some representation $R$. For some invertible matrix $O$ one can make the transformation

$T^{a}\rightarrow{T'}^{a}\equiv O_{ab}T^{b}\,, $

and the ${T'}^{a}$ matrices, just like the original $T^a$, will form a basis of the algebra in the $R$ representation. The same thing happens with the following transformation, as long as $U$ is an invertible matrix:

$T^{a}\rightarrow{T'}^{a}\equiv U^{-1}T^{b}U \,. $

In physics, usually one requires that the group transformations represented by the matrices $\exp\left(i\varepsilon^{a}T^{a}\right)$ are unitary, for some real parameters $\varepsilon^{a}$, which implies that the $T^{a}$ need to be hermitian matrices (${T^{a}}^\dagger=T^{a}$). Furthermore, it is convenient to have $C(R) Id=\sum_{a}T^{a}T^{a}$ and $T(R)\delta^{ab}=\textrm{Tr}\left(T^{a}T^{b}\right)$ (see section 2 above). These requirements constrain somewhat the above transformation freedom: if $T^{a}$ obeys these conditions, then the $O$ and $U$ above must be orthogonal and unitary, respectively, in order to keep having $C(R) Id=\sum_{a}{T'}^{a}{T'}^{a}$ and $T(R)\delta^{ab}=\textrm{Tr}\left({T'}^{a}{T'}^{b}\right)$.

Up to now we have been discussing the computation of quantities which are invariant under these $O$/$U$ basis transformations. However, one often wants to obtain basis-depend quantities (such as the Lagrangian). This inevitably requires taking a particular choice of basis.

Keeping this in mind, Susyno can make various basis-dependent calculations. The most basic one is to compute the representation matrices $T^a$ with the RepMatrices function. For example, the representation matrices of doublets of $SU(2)$ are half the Pauli matrices:

MatrixForm /@ RepMatrices[SU2, {1}] (* {1}=doublet representation *)

Just as easily one can calculate, say, the quintuplet representation matrices:

MatrixForm /@ RepMatrices[SU2, {4}] (* {4}=quintuplet representation *)

This function works for any representation of any group:

MatrixForm /@ RepMatrices[SU3, {1,0}] (* {1,0}=fundamental representation *)

Note carefully that these matrices are similar, but not exactly equal, to half the Gell-Mann matrices. However, both these and $1/2$ the GellMan matrices make a perfectly fine basis for the algebra of $SU(3)$, in the fundamental representation and therefore, one can check that there is an $O$/$U$ basis transformations which converts one set of matrices into the other.

It turns out that the basis used by Susyno for the $T^a$'s is such that there is a maximal number of diagonal generators (which are always given last by RepMatrices). This is usually a very convenient basis choice, but there is an important consequence for real representations which we will now discuss. A group transformation $V=\exp\left(i\varepsilon^{a}T^{a}\right)$ which is real, $V=V^*$, implies that the $T^{a}$ algebra matrices must be anti-symmetric, $T^{a}=-{T^{a}}^T$, and therefore none of them can be diagonal. This means that Susyno's basis choice is never the real one. For example, even though the triplet representation of $SU(2)$ is real, the $T^a$'s given by RepMatrices do not lead to a real $V=\exp\left(i\varepsilon^{a}T^{a}\right)$ as can be seen from the fact that not all them are anti-symmetric:

MatrixForm /@ RepMatrices[SU2, {2}]

However, the dedicated user can figure out which is the unitary transformation rotation $T^{a}\rightarrow{T'}^{a}\equiv U^{\dagger}T^{b}U$ which makes the ${T'}^{a}$'s anti-symmetric. In the present case, it turns out that

$U=\left(\begin{array}{ccc} \frac{1}{\sqrt{2}} & 0 & -\frac{i}{\sqrt{2}}\\ 0 & -i & 0\\ \frac{1}{\sqrt{2}} & 0 & \frac{i}{\sqrt{2}} \end{array}\right)$

will do the trick:

repMs = RepMatrices[SU2, {2}];
U = {{1/Sqrt[2], 0, -(I/Sqrt[2])}, {0, -I, 0}, {1/Sqrt[2], 0, I/Sqrt[2]}};
MatrixForm[ConjugateTranspose[U].#.U] & /@ repMs 

5. Group invariant combinations of representations

Often one needs to know how to contract, in a group invariant way, the components of a product of representations. The simplest example would be two $SU(2)$ doublets --- let us call them $D=(D_1,D_2)^T$ and $D'=({D'}_1,{D'}_2)^T$: it is well known that the combination $D_1 {D'}_2-D_2 {D'}_1$ is left invariant under the action of the $SU(2)$ group. To calculate group invariant combinations with Susyno, use the Invariants function:

Invariants[SU2,{{1},{1}}] (*product of two SU(2) doublets*)

The syntax is Invariants[<group>,{<rep1>, <rep2>, ...}], with an arbitrary number of representations. In the output, the program considers that the components of <rep1> are named a[1], a[2], etc., that those of <rep2> are b[1], b[2], etc., and so forth.

Notice that the output above is not just a[2] b[1] - a[1] b[2]; this expression is surrounded by curly brackets. The reason is, in general, there might more than one independent way of contracting the representations in an invariant way. To illustrate this point, consider the product of four $SU(2)$ doublets, which is known to have two independent invariants (this statement can be confirmed with the ReduceRepProduct function, by counting the number of singlets in the product of four doublets):

Invariants[SU2, {{1}, {1}, {1}, {1}}] (*product of four SU(2) doublets*)

The user should keep in mind that if $I_{1}$, $I_{2}$, ... are a list of independent group invariants then any linear combination of them will also be invariant. Therefore instead of the $I_a$ one could take instead $I'_a\equiv P_{ab} I_b$ for some invertible matrix $P$. For example, in the example above of the product of 4 doublets of $SU(2)$, the two invariants can be combined and presented in a different way. In fact, even when there is just one invariant expression, it can be arbitrarily multiplied by some numerical factor.

Another important point is that Invariants assumes that the representations transform in the basis given by RepMatrices. We can test this of course: if we write an invariant as $I=c_{ijk\cdots}a\left[i\right]b\left[j\right]c\left[k\right]\cdots$, the reader can check that the condition to be verified is


for all $a$, where $T^{a}$, ${T'}^{a}$, ${T''}^{a}$, ... are the representation matrices of $\left(a[1],a[2],...\right)^{T}$, $\left(b[1],b[2],...\right)^{T}$, $\left(c[1],c[2],...\right)^{T}$, ... respectively. The two invariants which can be formed with the product of four $SU(2)$ above can be tested as follows:

invariants = Invariants[SU2, {{1}, {1}, {1}, {1}}] ;

subsRuleA = Table[MapThread[Rule, {{a[1], a[2]}, mat.{a[1], a[2]}}],
{mat, RepMatrices[SU2, {1}]}];
subsRuleB = subsRuleA /. a -> b;
subsRuleC = subsRuleA /. a -> c;
subsRuleD = subsRuleA /. a -> d;

(invariants /. subsRuleA) + (invariants /. 
   subsRuleB) + (invariants /. subsRuleC) + (invariants /. subsRuleD)

The 3 groups of {0,0} in the result show the two invariants do not change under any of the three $SU(2)$ generators ($a=1,2,3$).

6. Permutation symmetries of products of representations

Section 3 explains how one can break a product of representations $R_{1}\times R_{2}\times\cdots\times R_{n}$ in its irreducible parts. However, sometimes, that is not enough: if some of the $R_i$ are the same, then there is a permutation symmetry to consider as well. For example, in $SU(2)$ we have $\mathbf{2}\times\mathbf{2}=\mathbf{1}+\mathbf{3}$ and it is well known that the singlet ($\mathbf{1}$) is an antisymmetric contraction of the two $\mathbf{2}$'s, while the triplet ($\mathbf{3}$) is symmetric. Usually this is indicated as $\mathbf{2}\times\mathbf{2}=\mathbf{1}_A+\mathbf{3}_S$. The function ReduceRepProduct introduced in section 3 does not provide this information about what happens under permutations of equal representations in a product $R_{1}\times R_{2}\times\cdots\times R_{n}$, but PermutationSymmetryOfTensorProductParts does.

This function computes such information with full generality, however some users might be surprised to know that at stake is not simply a question of adding S's and A's to each irreducible representation in a product. Perhaps is is then instructive to consider first what happens in a product of four $SU(2)$ doublets: $\mathbf{2}\times\mathbf{2}\times\mathbf{2}\times\mathbf{2}=\mathbf{1}+\mathbf{1}+\mathbf{3}+\mathbf{3}+\mathbf{3}+\mathbf{5}$. Are these irreps symmetric (S) or antisymmetric (A)? It turns out that such question is not well formulated: in the case of $\mathbf{2}\times\mathbf{2}=\mathbf{1}+\mathbf{3}$ there are just two equal factors in the product, so the relevant permutation group is $S_2$, which has two 1-dimensional irreducible representations: the symmetric one (S) and the antisymmetric one (A). In the case of $\mathbf{2}\times\mathbf{2}\times\mathbf{2}\times\mathbf{2}$, one should identify the irreducible $SU(2)$ representations contained in it with $S_4$ irreps. In general then, if there is a product with $n$ repeated representations $R$, one should look to the irreps of the $S_n$ permutation group, which are given by the partitions of $n$ (for example {4}, {3,1}, {2,2}, {2,1,1}, {1,1,1,1} in the case of $n=4$).

Having said this, we are ready to take a look at the output of PermutationSymmetryOfTensorProductParts. Consider first $\mathbf{2}\times\mathbf{2}$ (please be careful with the brackets):

PermutationSymmetryOfTensorProductParts[{SU2}, {{{1}}, {{1}}}]

The first part of the output, {{1,2}}, informs us that representation #1 and #2 are the same (which is quite obvious to the user). The second part says that

  1. there is an $SU(2)$ triplet (={2}) in the product which is in the {2}=symmetric representation of the $S_2$ induced by permutations of the doublets (multiplicity is 1);
  2. there is an $SU(2)$ singlet (={0}) in the product which is in the {1,1}=anti-symmetric representation of the $S_2$ induced by permutations of the doublets (multiplicity is 1).

In the case of $\mathbf{2}\times\mathbf{2}\times\mathbf{2}\times\mathbf{2}$, the code would be the following:

PermutationSymmetryOfTensorProductParts[{SU2}, {{{1}}, {{1}}, {{1}},{{1}}}]

It says that, under the $SU(2)_{S_{4}}$ group, the product breaks as $\mathbf{5}_{S}+\mathbf{3}_{\mathbf{3}}+\mathbf{1}_{\mathbf{2}}$ (using the notation $\left\{4\right\}=S$, $\left\{ 3,1\right\}=\mathbf{3}$, $\left\{ 2,2\right\}=\mathbf{2}$ for the $S_4$ representations). So, for example, neither of the two singlets in this product is completely symmetric nor anti-symmetric: the pair of singlets is converted into one-another by permutations of the doublets (the permutation symmetry entangles them), thereby they form an 2-dimensional irreducible representation of $S_4$.

Finally, note that the PermutationSymmetryOfTensorProductParts function is capable of handling more general cases, such as when there are more than one type of repeated representations being multiplied or when the Lie group is not simple:

PermutationSymmetryOfTensorProductParts[{SU2}, {{{1}}, {{2}}, {{1}}, {{2}}}]

PermutationSymmetryOfTensorProductParts[{SU3, SU2}, {{{1, 0}, {1}},
{{1, 0}, {1}}, {{1, 0}, {1}}}]

7. Symmetry breaking: branching rules and more

Concerning symmetry breaking, one usually wants to know how a given representation of a group breaks into irreducible representations of some subgroup. This can be calculated with DecomposeRep, which requires 3 elements from the user:

  1. the group $G$;
  2. the subgroup $H\subset G$;
  3. information on how $H$ is embedded in $G$ — the so-called projection matrix, to be specific;
  4. the representation $R$ of $G$ to be decomposed.
(1), (2) and (4) are trivial to provide, while (3) might require a little bit of work, but not much (more on this later). In any case, the user needs to figure out (3) only once for a given symmetry breaking pattern $G\rightarrow H$.

Consider the simple case where $SU(3)$ breaks into $SU(2)\times U(1)$. The projection matrix in this case can be chosen to be

$\left(\begin{array}{cc} 1 & 0\\ 1 & 2 \end{array}\right)$

We may then go ahead and see, for example, how does the triplet of $SU(3)$ decomposes:

rep={{1,0}}; (* 3 of SU(3) *)
prjMat = {{1, 0}, {1, 2}};


Likewise, consider $SU(5)\rightarrow SU(3)\times SU(2)\times U(1)$. Once the group, subgroup and projection matrix have been given, it is easy to write a few lines of code which will generate a table with the decomposition of $SU(5)$ representations:

Grid[Prepend[data,Style[#,Bold]&/@{"SU(5)","SU(3) x SU(2) x U(1) content"}],

So how can one compute the projection matrix? The function RegularSubgroupProjectionMatrix exists for this purpose. It requires some knowledge of Dynkin diagrams. For example, the diagram of $SU(5)$ is

while the one of $SU(3)\times SU(2)\times U(1)$ is

The latter can be obtained from the first by dropping one of the dots (representing a simple root) and converting it to a $U(1)$,

This is the information needed by RegularSubgroupProjectionMatrix: one must describe how the dots/simple roots of the subgroup are related to those of the original group. In this case, the $SU(3)$ dots are formed from the 1st and 2nd dots of $SU(5)$, so we can write this as {1,{1,2}}, with the first 1 referring to $SU(5)$, the only group factor in the original group. Likewise, the $SU(2)$ dot matches the 4th dot of $SU(5)$, so we may write this as {1,{4}}. The program will then notice that one dot of $SU(5)$ was unused (the third one) and assume that it was dropped. There is a $U(1)$ group associated to each dropped dot which may or may not, depending on the user's wish, be part of the subgroup.

Consider first the case $SU(5)\rightarrow SU(3)\times SU(2)$:

group = {SU5};
subgroup = {SU3, SU2};
breakingInfo = {{1, {1, 2}}, {1, {4}}};
prjMat = RegularSubgroupProjectionMatrix[group, subgroup, breakingInfo];
prjMat // MatrixForm

Now take the case $SU(5)\rightarrow SU(3)\times SU(2)\times U(1)$. As mentioned before, the program is aware that a $U(1)$ is "available" since the third dot of is not being used, so the user can only needs to say which linear combination of these "available" $U(1)$'s are indeed in the subgroup. In this case, there is just one $U(1)$, so one can choose at most a normalization factor $Z$:

group = {SU5};
subgroup = {SU3, SU2, U1};
breakingInfo = {{1, {1, 2}}, {1, {4}}, {Z}};
prjMat = RegularSubgroupProjectionMatrix[group, subgroup, breakingInfo];
prjMat // MatrixForm

(From here, it is clear $Z=1$ was used earlier to get the table in Out[65].) For more examples on how to use this RegularSubgroupProjectionMatrix, see the build-in documentation files of the program.

8. Functions related to the $S_n$ permutation group

Up to now, we have been discussing the calculation of quantities related to Lie groups. However, as explained in section 6, there is a connection with the discrete permutation group $S_n$, and for this reason the Susyno program also contains some $S_n$ related functions which we will now go through.

Recall that the irreducible representations of $S_n$ are given (i.e., can be labeled) by the partitions of $n$:

(* Example: the representations of S5 *)

(IntegerPartitions is a Mathematica built-in function.) Likewise, it is well known that the conjugacy classes of $S_n$ can also be labeled by the partitions of $n$, since this can be used to describe the length of the cycles of the group elements. For example, the identity element of $S_5$ can be written as $(1)(2)(3)(4)(5)$ in the cycle notation, meaning that it belongs to the {1,1,1,1,1} conjugacy class, while the $(143)(52)$ group element belongs to the class {3,2}.

With this reminder, we move on to the function SnIrrepDim which computes the dimension of $S_n$ irreducible representations. For example, {2,1} is the 2-dimensional representation of $S_3$:

SnIrrepDim[{2, 1}]

As a second example, we may compute the dimensions of the $S_5$ representations:

S5reps = IntegerPartitions[5]
dimS5reps = SnIrrepDim /@ S5reps
dimS5reps.dimS5reps (*check that the sum of the dimensions squared is 120=5!*)

The function SnClassCharacter is more general, as it calculates the character of a given conjugacy class in some irreducible representation of $S_n$. With it we may for example build the character table of $S_5$:

n = 5;
characterTable = Reverse /@ Table[SnClassCharacter[i, j],
{i, IntegerPartitions[n]}, {j, IntegerPartitions[n]}];

Print["Character table of ", Subscript[S, n], ":"];
characterTable // MatrixForm

The order of each conjugacy class (i.e., the number of elements in it) can be computed with SnClassOrder:

SnClassOrder[{1, 1, 1}]
SnClassOrder[{2, 1}]

With this data, one can easily decompose the product of $S_n$, and indeed that is the purpose of DecomposeSnProduct:

DecomposeSnProduct[{{2, 1}, {2, 1}}]

The output indicates the multiplicity of each $S_n$ irrep in the given product, as ordered by IntegerPartitions[n], so in this case the output is saying that $\left\{ 2,1\right\} \times\left\{ 2,1\right\} =\left\{ 3\right\} +\left\{ 2,1\right\} +\left\{ 1,1,1\right\}$. Here is a slightly more elaborate example with $S_5$:

multiplicities = DecomposeSnProduct[{{3, 2}, {3, 2}, {2, 2, 1}}];
Print["{3,2} x {3,2} x {2,2,1} contains ..."];
Do[Print["... ", multiplicities[[i]], " time(s) the representation ", 
    IntegerPartitions[5][[i]]];, {i, Length[multiplicities]}];

Finally, there is the function SnIrrepGenerators which calculates explicitly the representation matrices of the group elements $(12)$ and $(12\cdots)$ of $S_n$. Note that these two elements are enough to generate the whole group, with size $n!$. The matrices returned by this function are always real and orthogonal/unitary. For example:

MatrixForm /@ SnIrrepGenerators[{2, 1}]
MatrixForm /@ SnIrrepGenerators[{3, 2}]

With successive matrices multiplications, one can recover all group elements:

gens = SnIrrepGenerators[{2, 1}];
MatrixForm /@ FixedPoint[Sort[DeleteDuplicates[Join[#, 
Dot @@@ Tuples[#, 2]]]] &, gens]

Renato Fonseca


Last updated
26 November 2016