non negative solution of the matrix equation $A^T U A = U - C$ if C is non-negative Given $A$ a matrix with spectral radius smaller than 1 and a matrix $C$ symmetric.
It can be shown that $U=\sum_{k=0}^\infty (A^T)^k C A^k$ converges, is symmetric and is the solution of the equation above.
Is it possible to show that if $C$ is non-negative also $U$ is non-negative?
 A: Here is an alternative approach to solving the original equation.
The desired equation subject to the constraint $U \ge 0$ (elementwise nonnegativity) can be written as
\begin{equation*}
  U = C + A^TUA,\qquad U \ge 0.
\end{equation*}
Define now the iteration:
\begin{equation*}
  U_{k+1} = [C + A^TU_kA]_+,\qquad k=0,1,\ldots
\end{equation*}
where $[\cdot]_+ \equiv \max\lbrace 0, \cdot\rbrace$ denotes projection onto the nonnegative orthant. Assuming the slightly stronger condition that the operator norm $\|A\| < 1$, we can show that the above iteration converges to the desired solution.
$\newcommand{\Gc}{\mathcal{G}}$
Claim. Let $\Gc$ denote the map $U\mapsto [C + A^TUA]_+$. Then, $\Gc$ is a strict contraction in the spectral norm.

Proof. Let $U, V$ be two matrices. Then,
  \begin{eqnarray*}
\|\Gc(U)-\Gc(V)\| &=& \|[C + A^TUA]_+ - [C + A^TVA]_+\|\\
&\le& \|C + A^TUA - C - A^TVA]\|\\
&\le& \|A^T(U-V)A\|\\
&\le& \|A^T\| \|A\| \|U-V\| \\
&\le& \gamma\|U-V\|,\qquad \text{where}\ \gamma=\|A\|^2 < 1.
\end{eqnarray*}
  The first inequality follows as projection is a non-expansive operator, the rest follow from obvious properties of the spectral norm.

Thus, $\Gc$ is a strict contraction. So by the Banach contraction theorem the iteration $U_{k+1}=\Gc(U_k)$ has a unique solution. This solution is nonnegative by construction and satisfies the desired equation.
