# Solve matrix equation $X A Y + Y A^T X = B$ with diagonal $X$ and $Y$

Given square matrices $$A$$ and $$B$$, whereas $$B$$ is symmetric.

How to solve matrix equation (unknowns are $$X$$ and $$Y$$)

$$X A Y + Y A^T X = B$$ with diagonal $$X$$ and $$Y$$ (if a solution exists)? And for which $$A$$ and $$B$$ is there a solution?

From inspecting the diagonal component, we get $$B_{ii} = X_{ii} A_{ii} Y_{ii} + Y_{ii} A_{ii} X_{ii}$$ such that $$X_{ii} Y_{ii} = \frac{1}{2} \frac{B_{ii}}{A_{ii}}.$$

• So the unknowns are $X$ and $Y$?
– md5
May 8, 2020 at 22:37
• @md5 Yes, I will clarify this
– Jiro
May 8, 2020 at 22:38
• Multiplication with diagonal matrices is cummutative and leaves $XY(A+A^T)=B$ which is fearly easy to solve. Or am I wrong here? May 8, 2020 at 22:44
• @Laray I made the same mistake at first, the terms on the diagonal do not have to be the same.
– J P
May 8, 2020 at 22:45
• I have no idea if this would converge, but you could try setting $X$ or $Y$ equal to the identity matrix and solve for the other matrix by solving the resulting Lyapunov equation. Then set all the non-diagonal elements to zero and solve the new Lyapunov equation for the other matrix. Repeat until it hopefully converges. May 9, 2020 at 17:47

This isn't a complete answer but I think it might be useful. First observe that for any matrix, $$e_i^T M e_j = M_{ij}$$ and for a diagonal matrix $$e_i^T D = D_{ii} e_i^T$$ and $$D e_j = D_{jj} e_j$$ where $$e_i$$ is the vector with $$i$$th element is $$1$$ and all other elements are $$0$$. So, we have $$e_i^T X A Y e_j + e_i^T Y A^T X e_j = e_i^T B e_j$$ which simplifies to $$E_{ij} : A_{ij} X_i Y_j + A_{ji} X_j Y_i = B_{ij}$$ where $$X_i, Y_i$$ are the diagonal elements (for easier typing). Note that $$E_{ij} \equiv E_{ji}$$, so instead of $$n^2$$ equations, we have $$n(n+1)/2$$ equations to solve. But we only have $$2n$$ unknowns, so for $$n > 3$$ the number of equations exceed the number of unknowns.

Regardless, we can construct the following matrix equation: $$\underbrace{\begin{bmatrix} 2 A_{11} & 0 & 0 & 0 & 0 & \dots\\ 0 & A_{12} & A_{21} & 0 & 0 & \dots\\ 0 & 0 & 0 & A_{13} & A_{31} & \dots\\ \vdots & \vdots & \vdots & \vdots & \vdots & \ddots \end{bmatrix}}_{:=M} \underbrace{\begin{bmatrix} X_1 Y_1 \\ X_1 Y_2 \\ X_2 Y_1 \\ X_1 Y_3 \\ X_3 Y_1 \\ \vdots \end{bmatrix}}_{:=x} = \underbrace{\begin{bmatrix} B_{11} \\ B_{12} \\ B_{13} \\ \vdots \end{bmatrix}}_{:=b}$$

So one necessary condition for a solution to exist is $$b \in \operatorname{Im}M$$. However, this is not enough. But you can try to find a solution by gradient descent methods by selecting an initial guess and iterating through a solution after this point.

Edit. Example for $$n=3$$. $$\begin{bmatrix} 2 A_{11} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & A_{12} & A_{21} & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & A_{13} & A_{31} & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 2 A_{22} & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & A_{23} & A_{32} & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 2 A_{33} \end{bmatrix} \begin{bmatrix} X_1 Y_1 \\ X_1 Y_2 \\ X_2 Y_1 \\ X_1 Y_3 \\ X_3 Y_1 \\ X_2 Y_2 \\ X_2 Y_3 \\ X_3 Y_2 \\ X_3 Y_3 \end{bmatrix} = \begin{bmatrix} B_{11} \\ B_{12} \\ B_{13} \\ B_{22} \\ B_{23} \\ B_{33} \end{bmatrix}$$

• Thanks! I think that's a good solution for determining $X$ and $Y$ (if a solution exists). On the right side, isn't the vector $b = [B_{11}, B_{12}, B_{13}, \dots]$?
– Jiro
May 14, 2020 at 18:54
• Actually, with more thought, it's not clear to me how the n(n+1)/2 contributes to determining the $X$ and $Y$. Can you please elaborate more on that?
– Jiro
May 14, 2020 at 20:49
• Yes the right side vector is $b=\begin{bmatrix}B_{11} & B_{12} & \dots \end{bmatrix}^T$. But since $B$ is diagonal most of them are zero. Since the matrices are symmetric, you can only take the elements that are on the upper triangle (including diagonals). Note that $M$ is a matrix of size $n(n+1)/2 \times n^2$. May 15, 2020 at 7:52
• Why do you assume that $B$ is diagonal? I meant it to be symmetric.
– Jiro
May 15, 2020 at 8:28
• I don't know why I thought $B$ is diagonal. Anyway, it is the same idea. See my edit. May 15, 2020 at 8:53

Let $$A=[a_{i,j}],B=[b_{i,j}],X=[v_i],Y=[w_i]$$. When a solution exists, necessarily,

i) for every $$i$$, $$a_{i,i}=0 \implies b_{i,i}=0$$.

ii) for every $$i,j$$, $$a_{i,j}=a_{j,i}=0 \implies b_{i,j}=b_{j,i}=0$$.

Assume that we consider generic matrices $$A,B$$ (for example, randomly choose them). Then, with probability $$1$$, the $$(a_{i,j}),(b_{i,j})$$ are non-zero and the conditions i),ii) are fulfilled. In particular, the $$(v_i),(w_i)$$ are non zero and we can put $$v_1=1$$ (if $$(X,Y)$$ is a solution, then $$(aX,1/aY)$$ too).

Thus one considers a system of $$\dfrac{n(n+1)}{2}+1$$ equations in the $$2n$$ unknowns $$(a_{i,j}),(b_{i,j})$$. As noted by obareey, the system is overdetermined for $$n$$ large enough.

Indeed, we can show that the equations are algebraically independent (when $$A,B$$ are generic); then (if we want that there is at least one solution) we must find $$r_n=\dfrac{n(n+1)}{2}+1-2n=\dfrac{(n-1)(n-2)}{2}$$ relations linking the $$(a_{i,j}),(b_{i,j})$$.

For example, we randomly choose $$B$$ and $$n^2-r_n$$ entries of the matrix $$A$$. Then one has $$r_n$$ supplementary unknowns.

For $$n=2$$, no problem. For $$n=3$$, $$r_3=1$$. For $$n=4$$, $$r_4=3$$.

That follows is an example for $$n=4$$. $$u1,u2,u3$$ are the $$3$$ supplementary unknowns. One of the solutions for $$u1,u2,u3$$ The associated solution $$(X,Y)$$. Let's assume that the diagonals are nonzero, i.e., $$B_{ii} \neq 0$$ and $$A_{ii} \neq 0$$, else simply $$X_{ii} = Y_{ii} = 0$$. More compact, we can write $$x = diag(X)$$ and $$y = diag(Y)$$

Then with $$x_{i} \circ y_{i} = \frac{1}{2} \frac{B_{ii}}{A_{ii}} = d_{i}.$$ where $$\circ$$ is the Hadamard (=element-wise) product.

We can reformulate to (similar to obareey's answer) $$A_{ij} x_i y_j + A_{ji} x_j y_i = B_{ij}$$ and subsitute $$A_{ij} x_i \frac{d_j}{x_j} + A_{ji} x_j \frac{d_i}{x_i} = B_{ij}$$

By $$z_{ij} = x_i / x_j$$, we have $$A_{ij} d_j z_{ij} + A_{ji} d_i z_{ij}^{-1} = B_{ij}$$

And in matrix form: $$A D \circ Z + D A^{T} \circ Z^{\circ -1} = B$$ where $$^{\circ -1}$$ is an element-wise inverse. Thus $$A D \circ Z^{\circ 2} - B \circ Z + D A^{T} = 0$$ which are $$n^2$$ quadratic functions.

The remaining question is which of the two solutions of each of these $$n^2$$ equations satisfies $$Z = x x^{-T}$$?