Message Boards Message Boards

0
|
4610 Views
|
13 Replies
|
2 Total Likes
View groups...
Share
Share this post:

Finding the similarity transformation between two matrices

Hongyi Zhao
Posted 2 years ago

I've two matrices A and B, and I try to find the similarity transformation between them by the following code snippet adopted from the method described here:

In[74]:= (*Data*)

A = {{-I, 0, 0, 0}, {0, I, 0, 0}, {0, 0, 0, -I}, {0, 0, -I, 0}};
B = {{0, 1, 0, 0}, {-1, 0, 0, 0}, {0, 0, 0, I}, {0, 0, I, 0}};

(*Search for x and y based on characteristic polynomial*)
n = Length@A;
Id = IdentityMatrix@n;
solxy = SolveAlways[Det[A - l*Id] == Det[B - l*Id], l]

(*Update data*)
A = A /. solxy[[1]];
B = B /. solxy[[1]];

(*Solve for general P*)
P = Array[p, {n, n}];
solP = Solve[P . B == A . P, Flatten@P];
P = P /. solP[[1]]

(*Check*)
B == Inverse@P . A . P // Simplify


Out[78]= {{}}

During evaluation of In[74]:= Solve::svars: Equations may not give solutions for all "solve" variables.

Out[83]= {{p[1, 1], I p[1, 1], 
  p[1, 3], -p[1, 3]}, {p[2, 1], -I p[2, 1], p[2, 3], 
  p[2, 3]}, {p[3, 1], p[3, 2], p[3, 3], p[3, 4]}, {-I p[3, 2], 
  I p[3, 1], -p[3, 4], -p[3, 3]}}

Out[84]= True

In[88]:= P

Out[88]= {{p[1, 1], I p[1, 1], 
  p[1, 3], -p[1, 3]}, {p[2, 1], -I p[2, 1], p[2, 3], 
  p[2, 3]}, {p[3, 1], p[3, 2], p[3, 3], p[3, 4]}, {-I p[3, 2], 
  I p[3, 1], -p[3, 4], -p[3, 3]}}

As you can see, it gives a general form symbol matrix P as the result, and with all the symbols setting to 1, the following form will be obtained:

In[95]:= P /. {p[1, 1] -> 1, p[1, 3] -> 1, p[2, 1] -> 1, p[2, 3] -> 1,
   p[3, 1] -> 1, p[3, 2] -> 1, p[3, 3] -> 1, p[3, 4] -> 1}

Out[95]= {{1, I, 1, -1}, {1, -I, 1, 1}, {1, 1, 1, 1}, {-I, I, -1, -1}}

But this method is a little complicated and not intuitive. Any enhancements/suggestions/comments about this method will be highly appreciated. I also attached the code snippet shown above.

Regards, HZ

Attachments:
finding-the-simi...nb
POSTED BY: Hongyi Zhao
13 Replies
Sort By:
Hongyi Zhao
Hongyi Zhao
Posted 2 years ago

Do you mean that it may be limited by approximate accuracy, say, an exactly non-zero value been evaluated to zero internally, during the computation process?
POSTED BY: Hongyi Zhao

I cannot rule out that possibility. Then again, such an example would require careful construction.
POSTED BY: Daniel Lichtblau

I assume that JordanDecomposition in Mathematica works with arbitrary complex numbers, exact or approximate, not just rationals or algebraics:

JordanDecomposition[{{E, Pi}, {ArcTan[2], Sin[I]}}]
POSTED BY: Gianluca Gorni

Correct, it works with arbitrary complex numbers and/or symbolic expressions. It is however going to depend on functions such as PossibleZeroQ.
POSTED BY: Daniel Lichtblau
Hongyi Zhao
Hongyi Zhao
Posted 2 years ago

I've noticed the following technical explanation of the Jordan form based method implemented in sage math:

Warning

When the two matrices are similar, this routine may fail to find the similarity transformation. A technical explanation follows.

The similarity check is accomplished with rational form, which will be successful for any pair of matrices over the same field. However, the computation of rational form does not provide a transformation. So we instead compute Jordan form, which does provide a transformation. But Jordan form will require that the eigenvalues of the matrix can be represented within Sage, requiring the existence of the appropriate extension field. When this is not possible, a RuntimeError is raised, as demonstrated in an example below.

I wonder if the same limitation applies for the method discussed here.

Regards, HZ
POSTED BY: Hongyi Zhao
Hongyi Zhao
Hongyi Zhao
Posted 2 years ago

Thank you, Hans Dolhaine,

What are the differences and pros and cons between this approach and the approach suggested by Gianluca Gorni, that is, one based on JordanDecomposition?

On the surface, your method is cumbersome. So what other advantages does this method have?

Regards, HZ
POSTED BY: Hongyi Zhao

The differences? Gianlucas method use 2 Jordandecomposition, my method translates your problem to a sort of Gauss-elimination, which after all could be described as similar. Pros and cons? I don't know, I just wanted to see if this method works.

You could run each method (my could be a bit shorter: don't calculate the determinats of mats) a 1000 (or so) times and look which is faster.
POSTED BY: Hans Dolhaine

Hello Zhao, here is another approach to your problem :

(*Your matrices  *)
A = {{-I, 0, 0, 0}, {0, I, 0, 0}, {0, 0, 0, -I}, {0, 0, -I, 0}};
B = {{0, 1, 0, 0}, {-1, 0, 0, 0}, {0, 0, 0, I}, {0, 0, I, 0}};
(* Are they conjugates , do they have the same Eigenvalues (Thanks to Daniel)  ?*)
Eigenvalues[A]
Eigenvalues[B]
 (* The potential Matrix for the transformation *)
mm = Table[m[i, k], {i, 1, 4}, {k, 1, 4}];
mm // MatrixForm
(* Make an equation for your problem *)
e1 = # == 0 & /@ Flatten[A.mm - mm.B]
(* and transform it to an algebraic problem *)
e2 = CoefficientArrays[e1, Flatten[mm]]
(* get the according matrix *)
e3 = Normal[e2][[2]];
e3 // MatrixForm
(*check whether a solution is possible ( that means det = 0 ) *)
Det[e3]
(* There are 16 - rank of matrix = 8 solutions of the problem *)
MatrixRank[e3]
(*Find the solutions to the enhanced problem *)
vecs = NullSpace[e3]
(* Transform the solution vectors back to matrices *)
mats = Partition[#, 4] & /@ vecs
(* each fulfills the original problem, as it must be *)
Simplify[A.# - #.B] & /@ mats
(* but none is invertible  *)
Det[#] & /@ mats
(* so sum them up *)
mat = Total[mats]
Det[mat]
(* and *)
A.mat - mat.B
POSTED BY: Hans Dolhaine
Hongyi Zhao
Hongyi Zhao
Posted 2 years ago

Hi Daniel Lichtblau,

Exactly, the problem's background is here:

Suppose I already have a unitary irreducible complex representation for a group. In this case, can I find the conjugate/similar transformation matrices connecting the two unitary irreducible complex representations, i.e., the one used by GAP and the one I already have, of the same group?
POSTED BY: Hongyi Zhao

Do you expect the matrices to have the same eigenvalue?
POSTED BY: Daniel Lichtblau

You may try with JordanDecomposition:

a = {{-I, 0, 0, 0}, {0, I, 0, 0}, {0, 0, 0, -I}, {0, 0, -I, 0}};
b = {{0, 1, 0, 0}, {-1, 0, 0, 0}, {0, 0, 0, I}, {0, 0, I, 0}};
{a0, a1} = JordanDecomposition[a];
{b0, b1} = JordanDecomposition[b];
a1 == b1
b == b0 . b1 . Inverse[b0] == b0 . Inverse[a0] . a . a0 . Inverse[b0]
POSTED BY: Gianluca Gorni
Hongyi Zhao
Hongyi Zhao
Posted 2 years ago

Thank you, Gianluca Gorni,

In this method, if a1 != b1, what should I do?
POSTED BY: Hongyi Zhao
Hongyi Zhao
Hongyi Zhao
Posted 2 years ago

It seems that the is_similar function provided by sagemath is exactly for this purpose, as described here. But I'm not sure if there is an equivalent implementation in Wolfram.

Regards, HZ
POSTED BY: Hongyi Zhao
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Tag limit exceeded
Note: Only the first five people you tag will receive an email notification; the other tagged names will appear as links to their profiles.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract