Find non-unique Hermitian solution to pair of matrix equations

This is a continuation of my previous question here, but I previously did not fully specify some of the restrictions on the matrices and so the solutions posted there were not that helpful.

I'm looking to find a Hermitian solution $$X$$ to the pair of equations,

I'm trying to solve a pair of simultaneous matrix equations of the form,

$$AX + XA^\dagger + BJB^\dagger = 0,\ \ X C^\dagger + B J D^\dagger = 0,$$ where, $$J = \text{diag}\left(\begin{bmatrix}1&0\\0&-1\end{bmatrix},\dots,\begin{bmatrix}1&0\\0&-1\end{bmatrix}\right).$$

and $$J \in \mathbb{R}^{m \times m}$$.

Here $$A \in \mathbb{C}^{n\times n}$$, $$B \in \mathbb{C}^{n\times m}$$, $$C \in \mathbb{C}^{m\times n}$$, $$D \in \mathbb{C}^{m\times m}$$ where $$n$$ and $$m$$ are even

The former is the Lyapunov equation, but the solution $$X$$ may not be unique since $$A$$ may be unstable. Furthermore, the matrices $$A, B, C$$ may all be singular. So there are many solutions, but I need just one which is Hermitian. Does anyone know if there is a general way to find this (or the non-existence of this) for matrices with symbolic elements?

I tried something like this:

X=Array[x,{n,n}];
X=Solve[Simplify[c.a.X==-c.b.j.ConjugateTranspose[b]+d.j.ConjugateTranspose[b].ConjugateTranspose[a]]
&&Simplify[X==ConjugateTranspose[X]],Flatten[X]];


but it takes far too long.

Thanks!

EDIT: here is a concrete example as requested, here $$\gamma > 0$$

$$A=\left( \begin{array}{cc} -\frac{\gamma }{2 \gamma ^2+1} & \frac{\gamma \sqrt{4 \gamma ^2+1}}{2 \gamma ^2+1} \\ \frac{\gamma \sqrt{4 \gamma ^2+1}}{2 \gamma ^2+1} & -\frac{4 \gamma ^3+\gamma }{2 \gamma ^2+1} \\ \end{array} \right),$$

$$B=\left( \begin{array}{cc} 0 & \frac{1}{\sqrt{2 \gamma ^2+1}} \\ \sqrt{\frac{2 \gamma ^2+1}{4 \gamma ^2+1}} & -\frac{2 \gamma ^2}{\sqrt{8 \gamma ^4+6 \gamma ^2+1}} \\ \end{array} \right),$$

$$C=\left( \begin{array}{cc} -\frac{4 \gamma ^2}{\sqrt{\frac{1}{\gamma ^2}+2}} & -\frac{2 \left(4 \gamma ^3+\gamma \right)}{\sqrt{8 \gamma ^4+6 \gamma ^2+1}} \\ -2 \gamma \sqrt{2 \gamma ^2+1} & 0 \\ \end{array} \right),$$

$$D=I_{2\times 2}.$$

Then a Hermitian solution is

$$X=\left( \begin{array}{cc} -\frac{1}{4 \gamma ^3+2 \gamma } & \frac{\gamma }{\left(2 \gamma ^2+1\right) \sqrt{4 \gamma ^2+1}} \\ \frac{\gamma }{\left(2 \gamma ^2+1\right) \sqrt{4 \gamma ^2+1}} & \frac{1}{4 \gamma ^3+2 \gamma } \\ \end{array} \right).$$

The above equations in mathematica code

$Assumptions = {\[Gamma] > 0}; {a,b,c,d}={{{-(\[Gamma]/(1 + 2 \[Gamma]^2)), (\[Gamma] Sqrt[1 + 4 \[Gamma]^2])/( 1 + 2 \[Gamma]^2)}, {(\[Gamma] Sqrt[1 + 4 \[Gamma]^2])/( 1 + 2 \[Gamma]^2), -((\[Gamma] + 4 \[Gamma]^3)/( 1 + 2 \[Gamma]^2))}}, {{0, 1/Sqrt[1 + 2 \[Gamma]^2]}, {Sqrt[( 1 + 2 \[Gamma]^2)/( 1 + 4 \[Gamma]^2)], -((2 \[Gamma]^2)/Sqrt[ 1 + 6 \[Gamma]^2 + 8 \[Gamma]^4])}}, {{-((4 \[Gamma]^2)/Sqrt[ 2 + 1/\[Gamma]^2]), -((2 (\[Gamma] + 4 \[Gamma]^3))/Sqrt[ 1 + 6 \[Gamma]^2 + 8 \[Gamma]^4])}, {-2 \[Gamma] Sqrt[ 1 + 2 \[Gamma]^2], 0}}, {{1, 0}, {0, 1}}}; X={{-(1/(2 \[Gamma]+4 \[Gamma]^3)),\[Gamma]/((1+2 \[Gamma]^2) Sqrt[1+4 \[Gamma]^2])},{\[Gamma]/((1+2 \[Gamma]^2) Sqrt[1+4 \[Gamma]^2]),1/(2 \[Gamma]+4 \[Gamma]^3)}}  EDIT 2: Here is a slightly more complex example which takes long to compute $Assumptions = {\[Gamma] > 0, s0 > 0};
{a, b, c, d} =
{{{-((\[Gamma] (1 - s0 \[Gamma] + \[Gamma]^2))/(
1 + s0 \[Gamma] + \[Gamma]^2)), Sqrt[
s0 \[Gamma] (s0^2 \[Gamma]^2 -
2 s0 \[Gamma] (-1 + \[Gamma]^2) + (1 + \[Gamma]^2)^2)]/(
1 + s0 \[Gamma] + \[Gamma]^2)}, {(
Sqrt[(s0^5 \[Gamma]^3)/(1 + s0 \[Gamma] + \[Gamma]^2)] +
2 Sqrt[(s0^7 \[Gamma]^5)/(1 + s0 \[Gamma] + \[Gamma]^2)] +
2 Sqrt[(s0^5 \[Gamma]^7)/(1 + s0 \[Gamma] + \[Gamma]^2)] + Sqrt[(
s0^9 \[Gamma]^7)/(1 + s0 \[Gamma] + \[Gamma]^2)] -
2 Sqrt[(s0^7 \[Gamma]^9)/(1 + s0 \[Gamma] + \[Gamma]^2)] + Sqrt[(
s0^5 \[Gamma]^11)/(1 + s0 \[Gamma] + \[Gamma]^2)])/(
s0^2 \[Gamma] Sqrt[(1 +
s0 \[Gamma] + \[Gamma]^2) (s0^2 \[Gamma]^2 -
2 s0 \[Gamma] (-1 + \[Gamma]^2) + (1 + \[Gamma]^2)^2)]), \
\[Gamma] (-1 - (2 s0 \[Gamma])/(1 + s0 \[Gamma] + \[Gamma]^2))}}, {{0,
1/Sqrt[1 + s0 \[Gamma] + \[Gamma]^2]}, {Sqrt[(
1 + s0 \[Gamma] + \[Gamma]^2)/(
s0^2 \[Gamma]^2 -
2 s0 \[Gamma] (-1 + \[Gamma]^2) + (1 + \[Gamma]^2)^2)], -2 \
\[Gamma] Sqrt[(
s0 \[Gamma])/((1 + s0 \[Gamma] + \[Gamma]^2) (s0^2 \[Gamma]^2 -
2 s0 \[Gamma] (-1 + \[Gamma]^2) + (1 + \[Gamma]^2)^2))]}}, \
{{-4 Sqrt[(s0 \[Gamma]^5)/(
1 + s0 \[Gamma] + \[Gamma]^2)], -2 \[Gamma] Sqrt[(
s0^2 \[Gamma]^2 -
2 s0 \[Gamma] (-1 + \[Gamma]^2) + (1 + \[Gamma]^2)^2)/(
1 + s0 \[Gamma] + \[Gamma]^2)]}, {-2 \[Gamma] Sqrt[
1 + s0 \[Gamma] + \[Gamma]^2], 0}}, {{1, 0}, {0, 1}}};


A Hermitian solution is

X = {{-(1/(2 \[Gamma] (1+s0 \[Gamma]+\[Gamma]^2))),Sqrt[(s0 \[Gamma])/(s0^2 \[Gamma]^2-2 s0 \[Gamma] (-1+\[Gamma]^2)+(1+\[Gamma]^2)^2)]/(1+s0 \[Gamma]+\[Gamma]^2)},{Sqrt[(s0 \[Gamma])/(s0^2 \[Gamma]^2-2 s0 \[Gamma] (-1+\[Gamma]^2)+(1+\[Gamma]^2)^2)]/(1+s0 \[Gamma]+\[Gamma]^2),1/(2 \[Gamma]+2 s0 \[Gamma]^2+2 \[Gamma]^3)}};

• Now I see what you mean. This is a nice problem. But you need to formulate it better focusing on the essential part, namely, on analytically solving linear equations with complicated coefficients. I think the problem is that in the course of solution MA does not simplify intermediate expressions and, therefore, cannot decide if certain terms are zero. Jul 8, 2020 at 19:36
• @yarchik Yes that would seem so, I worry that maybe there isn't a neat and fast solution! Which would be a shame as I need this for my quantum network synthesis algorithm. I think I'll have to solve it on a case-by-case basis for now Jul 9, 2020 at 8:05
• It seems so. MA can still be helpful in performing intermediate simplifications. Jul 9, 2020 at 11:35

Combining equations we obtain a conventional linear equation in fully unique way without any constraints on the matrices:

$$\left\{ \begin{array}{l}XA^\dagger C^\dagger=ABJD^\dagger-BJB^\dagger C^\dagger,\\ X C^\dagger + B J D^\dagger = 0;\end{array}\right.$$ or transposing

$$\left\{ \begin{array}{l}CAX^\dagger=D J B^\dagger A^\dagger-C B JB^\dagger,\\ C X^\dagger =- D J B^\dagger; \end{array}\right.$$

Just solve it using LinearSolve when matrices are not singular:

xd=LinearSolve[C.A, D.J.ConjugateTranspose[B.A]-C.B.J.ConjugateTranspose[B]]
x=ConjugateTranspose[xd]


It works also for symbolic matrices.

Notice that without any additional information on $$A,\ldots, D$$ it is impossible to tell whether $$X$$ is Hermitian.

In order to obtain the solution for the modified question we combine two equations and use Solve as in OP:

\$Assumptions = {g > 0};

j={{1,0},{0,-1}};
a={{-(g/(1 + 2 g^2)), (g Sqrt[1 + 4 g^2])/(1 + 2 g^2)}, {(g Sqrt[1 + 4 g^2])/(1 + 2 g^2), -((g + 4 g^3)/(1 + 2 g^2))}};
b={{0, 1/Sqrt[1 + 2 g^2]}, {Sqrt[(1 + 2 g^2)/( 1 + 4 g^2)], -((2 g^2)/Sqrt[1 + 6 g^2 + 8 g^4])}};
c={{-((4 g^2)/Sqrt[ 2 + 1/g^2]), -((2 (g + 4 g^3))/Sqrt[ 1 + 6 g^2 + 8 g^4])}, {-2 g Sqrt[1 + 2 g^2], 0}};
d={{1, 0}, {0, 1}};

u=Join[c,c.a//FullSimplify]
v=Join[-d.j.ConjugateTranspose[b]//FullSimplify,
d.j.ConjugateTranspose[a.b] - c.b.j.ConjugateTranspose[b]//FullSimplify]
X=Array[x,{2,2}]

Solve[u.X==v]//FullSimplify

• 1. It is not unique since CA is not invertible 2. This does not guarantee X to be Hermitian Jul 8, 2020 at 12:26
• I added a concrete example. It is not clear to me that they are superfluous, especially since it is possible to find a solution to the former equation which is not a solution to the latter using "LyapunovSolve" Jul 8, 2020 at 12:37
• Unfortunately this becomes really slow on slightly more complicated examples, see the second edit on my question. I left it running for about an hour and solve was still running Jul 8, 2020 at 15:07