Problem

The trace algebra of two matrices is the polynomial ring generated by the traces of products of matrices A and B, like Tr(A), Tr(B), Tr(A.B), Tr(A.A), Tr(A.A.B.A.B). It is known that this trace algebra is the complete set of polynomials that are invariant under simultaneous conjugation of A and B (ie, if one replaces A and B by T^(-1).A.T and T^(-1).B.T, respectively, then the value of the polynomial doesn’t change).

Also, in the case of 2-by-2 matrices, one knows that in fact, Tr(A), Tr(B), Tr(A.B), Tr(A.A) and Tr(B.B) are enough to generate the ring: every invariant polynomial can in theory be written as a polynomial in these five traces.

In particular, the polynomial f=Tr(A.B.A.B – A.A.B.B), being invariant, is a polynomial in these five traces. For my research, I wanted to find out how this polynomial looked like.

My general idea was quite simple. Looking at the degree of the polynomial f, it was clear that it would have to be a polynomial of degree 4 and less. Now, there are 126 different monomials of degree 4 or less in the 5 traces, so we can write:

f=a_0 + a_1 Tr(A) + a_2 Tr(B) + … + a_125 Tr(B^2)^4

Now, if we take two artibrary matrices A and B and fill this in, we get a linear equation in a_0,…,a_125. Doing this 126 times, we get a linear system that should, in principe, completely determine the 126 coefficients, giving us f as a polynomial in the traces.

(Of course, this is a well-known method, but the fact that I have polynomials in which to evaluate the points rather than just a 126-dimensional space makes it a bit less standard. Even though, quite possibly Mathematica has some easier way of doing it than how I proceeded…)

So, for starters, I generated a list of all the possible monomials of degree 4 or less in the traces:

l0 = {1, ta, tb, tab, taa, tbb}
{1, ta, tb, tab, taa, tbb}
l = Union[Map[Fold[Times, 1, #] &, Tuples[l0, 4]]]
{1, ta, ta^2, ..., tb tbb^3, (tbb^4)}

And looked how these look like polynomials:

A := {{a, b}, {c, d}}; B := {{e, f}, {g, h}}
Ltt = l /. {ta -> Tr[A], tb -> Tr[B], taa -> Tr[A.A],
      tbb -> Tr[B.B], tab -> Tr[A.B]}
{1, a + d, (a + d)^2, ..., (e^2 + 2 f g + h^2)^4}

Next, I generated 126 random tuples of 2-by-2 matrices and evaluated the monomials in all these tuples:

R = RandomReal[1, {126, 2, 2, 2}]
evalfun = (r1 = #;
   s = {a -> r1[[1, 1, 1]], b -> r1[[1, 1, 2]], c -> r1[[1, 2, 1]],
     d -> r1[[1, 2, 2]], e -> r1[[2, 1, 1]], f -> r1[[2, 1, 2]],
     g -> r1[[2, 2, 1]], h -> r1[[2, 2, 2]]}; Map[# /. s &, Lt]) &
mat = Map[evalfun, R]

This indeed gives a 126-by-126 matrix with all the monomials evaluated. Now, we also need to find the list of values of our polynomial f in the same way:

ourfun = Tr[A.B.A.B - A.A.B.B]
evalone = (r1 = #;
   s = {a -> r1[[1, 1, 1]], b -> r1[[1, 1, 2]], c -> r1[[1, 2, 1]],
     d -> r1[[1, 2, 2]], e -> r1[[2, 1, 1]], f -> r1[[2, 1, 2]],
     g -> r1[[2, 2, 1]], h -> r1[[2, 2, 2]]}; ourfun /. s) &
list = Map[evalone, R]

Finally, using LinearSolve we can try to find our coefficients:

sol = LinearSolve[mat, list]
{6.18939*10^-13, 4.08603*10^-13, -2.01183*10^-12,
 1.59512*10^-12, -6.31879*10^-13, ... }

As you can see, most of the coefficients are very close to being zero, while some other coefficients turned out to be nicely close to round values such as 1 or 1/2. To find our hypothesised polynomial relations, we can use:

For[i = 1, i <= 126, i++,
 If[Abs[li[[i]]] > 0.01,
  Print[li[[i]], l[[i]]]
  ]
 ]
1.tab^2
-1.ta tab tb
0.5taa tb^2
0.5ta^2 tbb
-1.taa tbb

And indeed, Mathematica can confirm this was the polynomial we were looking for:

Simplify[(Tr[A.B]^2 - Tr[A] Tr[A.B] Tr[B] + 1/2 Tr[A.A] Tr[B]^2 +
    1/2 Tr[A]^2 Tr[B.B] - Tr[A.A] Tr[B.B]) - Tr[A.B.A.B - A.A.B.B]]
0

Nice!

Discussion

One interesting point is that the 126-by-126 matrix we calculated seemes almost singular: the determinant calculated was something like 10^(-111). This would mean that in general we would not expect LinearSolve to work. And wouldn’t this suggest some linear dependence between the monomials (which shouldn’t be the case: they should be algebraically independent)?

Also, I’m quite curious how this scales: in the case of two 2-by-2 matrices and a degree<5 relation, which is almost the simplest case imaginable, one gets 126-by-126 matrices which are quite viable to calculate with. But how about higher degrees? Or larger matrices?

Last Modified : April 9th, 2009 Filed under : Mathematica Cookbook Navigate : Previous post / Next post