# How to find a unitary transform between two matrices?

Given two matrices,

G1={{0, 0, 1, 0}, {0, 0, 0, 1}, {1, 0, 0, 0}, {0, 1, 0, 0}}


and

G2={{1, 0, 0, 0}, {0, 1, 0, 0}, {0, 0, -1, 0}, {0, 0, 0, -1}}


how do you find a unitary transformation matrix $$U$$ so that $$U.U^\dagger = I = U^\dagger.U$$ and $$G_2 = U.G_1.U^\dagger$$? Here $$^\dagger$$=ConjugateTranspose. I tried to solve it brute-force by setting up 16 unknown variables, but Mathematica is stuck trying to solve it. Is there a systematic/clever way of finding it?

• Welcome to Mathematica.SE! I hope you will become a regular contributor. To get started, 1) take the introductory tour now, 2) when you see good questions and answers, vote them up by clicking the gray triangles, because the credibility of the system is based on the reputation gained by users sharing their knowledge, 3) remember to accept the answer, if any, that solves your problem, by clicking the checkmark sign, and 4) give help too, by answering questions in your areas of expertise. Commented Jul 30, 2021 at 4:26
• JordanDecomposition is very similar to your goal. Commented Jul 30, 2021 at 6:20
• Can you, please, show the attempt in code that you have referenced? Thanks! Commented Jul 30, 2021 at 6:38

You can achieve this through the eigenvectors of the matrices:

G1 = {{0, 0, 1, 0}, {0, 0, 0, 1}, {1, 0, 0, 0}, {0, 1, 0, 0}};
G2 = {{1, 0, 0, 0}, {0, 1, 0, 0}, {0, 0, -1, 0}, {0, 0, 0, -1}};

v1 = Normalize /@ Eigenvectors[G1];
v2 = Normalize /@ Eigenvectors[G2];


Now we have

v1 . G1 . Transpose[v1] == v2 . G2 . Transpose[v2]
(*    True    *)


which means that what you want to achieve is possible. We define

U = Transpose[v2] . v1


$$\left( \begin{array}{cccc} \frac{1}{\sqrt{2}} & 0 & \frac{1}{\sqrt{2}} & 0 \\ 0 & \frac{1}{\sqrt{2}} & 0 & \frac{1}{\sqrt{2}} \\ -\frac{1}{\sqrt{2}} & 0 & \frac{1}{\sqrt{2}} & 0 \\ 0 & -\frac{1}{\sqrt{2}} & 0 & \frac{1}{\sqrt{2}} \\ \end{array} \right)$$

and verify that

U . G1 . Transpose[U] == G2
(*    True    *)

• This is the most elegant. Commented Jul 31, 2021 at 18:05

Given:

ClearAll["Global*"]
g1 = {{0, 0, 1, 0}, {0, 0, 0, 1}, {1, 0, 0, 0}, {0, 1, 0, 0}};
g2 = {{1, 0, 0, 0}, {0, 1, 0, 0}, {0, 0, -1, 0}, {0, 0, 0, -1}};


Assume matrix $$U$$ has the following Hermitian form:

upper = {{0, x12 + I y12, x13 + I y13, x14 + I y14},
{0, 0, x23 + I y23, x14 + I y24},
{0, 0, 0, x34 + I y34},
{0, 0, 0, 0}};
u = upper + ConjugateTranspose[upper] +
DiagonalMatrix[{x11, x22, x33, x44}] // ComplexExpand;
u // MatrixForm


$$\left( \begin{array}{cccc} \text{x11} & \text{x12}+i \text{y12} & \text{x13}+i \text{y13} & \text{x14}+i \text{y14} \\ \text{x12}-i \text{y12} & \text{x22} & \text{x23}+i \text{y23} & \text{x14}+i \text{y24} \\ \text{x13}-i \text{y13} & \text{x23}-i \text{y23} & \text{x33} & \text{x34}+i \text{y34} \\ \text{x14}-i \text{y14} & \text{x14}-i \text{y24} & \text{x34}-i \text{y34} & \text{x44} \\ \end{array} \right)$$

Use FindInstance to find a $$U$$ that satisfies $$G_2=U\cdot G_1\cdot U^\dagger$$, but is not necessarily unitary:

vars = Variables[Det[u]];
product = u.g1.ConjugateTranspose[u] // ComplexExpand;
eqns = Flatten@Table[g2 == product, {irow, 4}, {jcol, 4}];
soln = FindInstance[eqns, vars];

usol = u /. First[soln];
usol // MatrixForm


$$\left( \begin{array}{cccc} 0 & -\frac{i}{\sqrt{2}} & 0 & -\frac{i}{\sqrt{2}} \\ \frac{i}{\sqrt{2}} & 0 & \frac{i}{\sqrt{2}} & 0 \\ 0 & -\frac{i}{\sqrt{2}} & 0 & \frac{i}{\sqrt{2}} \\ \frac{i}{\sqrt{2}} & 0 & -\frac{i}{\sqrt{2}} & 0 \\ \end{array} \right)$$

Verify that $$U$$ is really a solution and happens to be unitary.

g2 == usol.g1.ConjugateTranspose[usol] // ComplexExpand // Simplify    (*  True  *)

UnitaryMatrixQ[usol]    (*  True  *)

• I tried something like this. The upper triangular is a nice idea in reducing the number of variables. This solution is different from the one by Roman, but is acceptable since the solution is not unique. Commented Aug 1, 2021 at 17:24
G1 = {{0, 0, 1, 0}, {0, 0, 0, 1}, {1, 0, 0, 0}, {0, 1, 0, 0}};
G2 = {{1, 0, 0, 0}, {0, 1, 0, 0}, {0, 0, -1, 0}, {0, 0, 0, -1}};

U = ToExpression[Table["a" <> ToString[i] <> ToString[j], {i, 4}, {j, 4}]];


To obtain values for the entries you can minimize a sum of squares of the equality constraints:

var = Variables[U];
varC = ToExpression[Flatten[{# <> "R", # <> "I"} &[ToString[#]] & /@ var]];

err = Total[ComplexExpand[Abs[{G2.U - U.G1, U.ConjugateTranspose[U] - IdentityMatrix[4]}]^2,
var, TargetFunctions -> {Re, Im}], ∞] /. conv;

SeedRandom[123]
cur = RandomReal[{-1, 1}, Length[varC]];
Do[cur -= PseudoInverse[D[err, {varC, 2}] /. #].(D[err, {varC, 1}] /. #) &[Thread[varC -> cur]], {15}]
err /. Thread[varC -> cur] // Abs
(* 3.3306691•(10^-16) *)


There are multiple solutions so PseudoInverse is more stable to use. The errors are small enough for Chop to discard them:

U = U /. Thread[var -> Partition[cur, 2].{1, I}];
G2 - U.G1.ConjugateTranspose[U] // Chop
U.ConjugateTranspose[U] // Chop
(* {{0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0}}
{{1., 0, 0, 0}, {0, 1., 0, 0}, {0, 0, 1., 0}, {0, 0, 0, 1.}} *)


In this particular case Solve can give a nice set of necessary equalities that could be substituted in the matrix to begin with to have fewer variables to optimize:

Solve[G2.U - U.G1 == 0]
(* {{a13 -> a11, a14 -> a12, a23 -> a21, a24 -> a22, a33 -> -a31, a34 -> -a32, a43 -> -a41, a44 -> -a42}} *)
`