There has been a flurry of activity in recent years in the area of solution of matrix equations. In particular, a good understanding has been reached on how to approach the solution of large scale Lyapunov equations. An effective way to solve Lyapunov equations of the form , where  and  are , is to use Galerkin projection with appropriate extended or rational Krylov subspaces. These methods work in part because the solution is known to be symmetric positive definite with rapidly decreasing singular values, and therefore it can be approximated by a low rank matrix . Thus the computations are performed usually with storage which is lower rank, i.e., much lower than order of .

Generalized Lyapunov equations have additional terms. In this paper, we concentrate on equations of the following form

Such equations arise for example in stochastic control. Seeveral authors have proposed to use approximations to conjugate gradients or to BiCGStab, appropriately preconditioned, where the basis vectors (matrices) and iterates are “truncated” throughout the algorithm to keep all these elements represented by low-rank matrices.

In the present work, we propose a return to classical iterative methods, and consider instead stationary iterations, where the iterate  is the solution of

where .The classical theory of splittings applies here, and we present a new theorem on the convervegence when the linear system at each step is solved inexactly.

One of the advantages of this classical approach is that only the data and the low-rank factors of the old and new iterates  and  need to be kept in storage. Furthermore, the solutions of the Lyapunov equations at each iteration can be performed with the Galerkin projection methods mentioned above, where the growth of rank can usually be well contained.

Numerical experiments show the competitiveness of the proposed approach.

This is joint work with Stephen D. Shank and Valeria Simoncini.