\documentclass[12pt]{article}

\usepackage{amsmath}
\usepackage{amsthm}
\usepackage{amssymb}
\usepackage{amsfonts}
\newtheorem{theorem}{Theorem}[section]
\newtheorem{corollary}{Corollary}[section]
\def\Z{\mathbb Z}
\def\CC{\mathbb C}
\def\rank{\hbox{\rm \ rank}}
\begin{document}
%Fassbender
%Heike","h.fassbender@tu-bs.de","
\section*{On the numerical solution of large-scale sparse discrete-time Riccati equations} By {\sl Heike Fa\ss bender and Peter Benner}.  

\medskip  \noindent Inspired by a large-scale sparse discrete-time Riccati equation which arises in a spectral factorization problem the efficient numerical solution of such Riccati equations is studied in this work. Spectral factorization is a crucial step in the solution of linear quadratic estimation and control problems. A variety of methods has been developed over the years for the computation of canonical spectral factors for processes with rational spectral densities, see, e.g., the survey \cite{SayK01}. One approach involves the spectral factorization via a discrete-time Riccati equation.  Whenever possible, we consider the generalized discrete--time algebraic Riccati equation  \begin{eqnarray}  0 ~=~ \mathcal{R}(X) &=& C^TQC + A^T X A - E^T X E \label{dare} \\ &&\;\; - (A^T XB + C^T S) (R + B^TXB)^{-1} (B^T XA + S^T C), \nonumber\end{eqnarray}  where $A, E \in \mathbb{R}^{n \times n}, B \in \mathbb{R}^{n \times m}, C \in \mathbb{R}^{p  \times n}, Q \in \mathbb{R}^{p \times p}, R \in \mathbb{R}^{m \times m},$ and $S  \in \mathbb{R}^{p \times m}.$ Furthermore, $Q$ and $R$ are assumed to be  symmetric and $A$ and $E$ are large and spare. For the particular application above, we have  \[ A = \left[ \begin{array}{cccc}      0 & 1 &  \\      & \ddots & \ddots \\      &&0 & 1\\      &&& 0\end{array}\right].      \] The function $\mathcal{R}(X)$ is a rational matrix function, $\mathcal{R}(X) = 0$ defines a system of nonlinear equations. Newton's method for the numerical solution of DAREs can be formulated as follows\\ \phantom{BBBB} {\bf for} {$k = 0,\,1,\,2,\,\ldots$}\\  \phantom{BBBBB} 1. $K_k \gets K(X_k) = (R + B^T X_k B)^{-1} (B^T X_k A + S^T  C)$.\\  \phantom{BBBBB}  2. $A_k \gets A - B K_k$.\\  \phantom{BBBBB} 3. $\mathcal{R}_k \gets \mathcal{R}(X_k)$.\\  \phantom{BBBBB} 4. Solve for $N_k$ in the Stein equation \begin{equation}\label{stein} A_k^T N_k A_k - E^T N_k E  = -\mathcal{R}_k. \end{equation}  \phantom{BBBBB} 5. $X_{k+1} \gets X_k + N_k.$\\  \phantom{BBBB}{\bf end for}  The computational cost for this algorithm mainly depends upon the cost for the numerical solution of the Stein equation (\ref{stein}). This can be done using the Bartels--Stewart algorithm \cite{BarS72} or an extension to the case $E \not= I$ \cite{GarLAM92,GarWLAM92,Pen97}. The Bartels-Stewart algorithm is the standard direct method for the solution of Stein equations of small to moderate size.  This method requires the computation of a Schur decomposition, and thus is not appropriate for large scale problems. The cost for the solution of the Stein equation is $\approx 73n^3$ flops. Iterative schemes have been developed including the Smith method \cite{Smi68}, the sign-function method \cite{Rob80}, and the alternating direction implicit (ADI) iteration method \cite{Wac88}. Unfortunately, all of these methods compute the solution in dense form and hence require ${\cal O}(n^2)$ storage. In case the solution to the Stein equation has low numerical rank (i.e., the eigenvalues decay rapidly) one can take advantage of this low rank structure to obtain approximate solutions in low rank factored form. If the effective rank is $r \ll n$, then the storage is reduced from ${\cal O}(n^2)$ to ${\cal O}(nr)$. This approach will be discussed here in detail.  \begin{thebibliography}{10}  \bibitem{BarS72} {\sc R.H. Bartels and G.W. Stewart}, {\em Solution of the matrix equation   ${AX}+{XB}={C}$: {A}lgorithm 432}, Comm. ACM, 15 (1972), pp.~820--826.  \bibitem{GarLAM92} {\sc J.D. Gardiner, A.J. Laub, J.J. Amato, and C.B. Moler}, {\em Solution of   the {S}ylvester matrix equation ${AXB}+{CXD}={E}$}, {ACM} Trans. Math.   Software, 18 (1992), pp.~223--231.  \bibitem{GarWLAM92} {\sc J.D. Gardiner, M.R. Wette, A.J. Laub, J.J. Amato, and C.B. Moler}, {\em   Algorithm 705: A {F}ortran-77 software package for solving the {S}ylvester   matrix equation ${AXB^T}+{CXD^T}={E}$}, {ACM} Trans. Math. Software, 18   (1992), pp.~232--238.  \bibitem{Pen97} {\sc T.~Penzl}, {\em Numerical solution of generalized {L}yapunov equations},   Adv. Comp. Math., 8 (1997), pp.~33--48.  \bibitem{Rob80} {\sc J.D. Roberts}, {\em Linear model reduction and solution of the algebraic   {R}iccati equation by use of the sign function}, Internat. J. Control, 32   (1980), pp.~677--687. \newblock (Reprint of Technical Report No. TR-13, CUED/B-Control, Cambridge   University, Engineering Department, 1971).  \bibitem{SayK01} {\sc A.H. Sayed and T.~Kailath}, {\em A survey of spectral factorization   methods}, Num. Lin. Alg. Appl., 8 (2001), pp.~467--496.  \bibitem{Smi68} {\sc R.A. Smith}, {\em Matrix equation {$XA + BX = C$}}, {SIAM} J. Appl. Math.,   16 (1968), pp.~198--201.  \bibitem{Wac88} {\sc E.L. Wachspress}, {\em Iterative solution of the {L}yapunov matrix   equation}, Appl. Math. Letters, 107 (1988), pp.~87--90.  \end{thebibliography}
%discrete-time algebraic Riccati equation, Stein equation, large, sparse, Newton method","15A24","","MInsymposium ""MATRIX FUNCTIONS AND MATRIX EQUATIONS"",","08:36:48","Wed Mar 26 2008","134.169.54.97"
\end{document}
 zzzz
