What does the matrix-valued solution of $X = A X A^T + \operatorname{Id}$ look like?

Question

Let $A \in \mathbb{R}^{n \times n}$ be an invertible contraction, i.e. all singular values are in $(0,1)$. By reformulating the equation \begin{align*} & X = A X A^T + \operatorname{Id} \tag{1} \end{align*} into a Sylvester-equation $A^{-1} X + X(-A^T) = A^{-1}$, the observation that $A^{-1}$ and $A^T$ share no eigenvalues guarantees that there is a unique solution $X \in \mathbb{R}^{n \times n}$.

If $A$ is symmetric, we can use the spectral decomposition $A = V \operatorname{diag}(a_1,...,a_n) V^T$ to check, that $X = V \operatorname{diag}(\frac{1}{1-a_1^2},...,\frac{1}{1-a_n^2}) V^T$ is the solution.

Question: What does $X$ look like, if $A$ is not symmetric? Is there also an explicit way of writing $X$ (preferably without vectorization)?

Any help is much appreciated!

Answer 1

Consider the linear operator $S:X\mapsto AXA^T$, on the space of matrices $\mathbb R^{n\times n}$, with the spectral norm. Its operator norm is $\|A\|^2<1$ so $1-S$ is invertible, with inverse given by the Neumann series $\sum_{k=0}^\infty S^k$. In particular the given equation has a unique solution given by $(\sum_{k=0}^\infty S^k)(I)=\sum_{k=0}^\infty A^k(A^k)^T$ .

Answer 2

A partial generalization of your explicit formula for the symmetric case in terms of the eigendecomposition is possible in some circumstances. Suppose $A$ is not necessarily symmetric, but diagonalizable, $A V = V \Lambda$ with $V$ nonsingular and $\Lambda = \text{diag}(\lambda_1, \dotsc, \lambda_n)$. Suppose additionally that no product of two eigenvalues is 1 (automatically satisfied in your case). Then $$X = V [F \circ (V^T V)^{-1} ] V^T, \quad F_{ij} = \frac{1}{1-\lambda_i \lambda_j}$$ is a solution to $X = A X A^T + I$, where $\circ$ denotes the Schur/Hadamard product. Note this generalizes your explicit solution for symmetric $A$, where we can take $V^TV = I$.

This follows from vectorizing. The equation can be written $$(I \otimes I - A \otimes A) \text{vec}(X) = \text{vec}(I).$$ Multiply on the left by $V^{-1} \otimes V^{-1}$, using that $V^{-1} A = \Lambda V^{-1}$, $$(I \otimes I - \Lambda \otimes \Lambda) (V^{-1} \otimes V^{-1}) \text{vec}(X) = (V^{-1} \otimes V^{-1}) \text{vec}(I).$$ From this it follows that $$ \text{vec}(V^{-1} X V^{-T}) = (I \otimes I - \Lambda \otimes \Lambda)^{-1} \text{vec}(V^{-1} V^{-T}). $$ "Unvectorizing" then gives the form given above.