Quality of Approximate Balanced Truncation

Lei-Hong Zhang School of Mathematical Sciences, Soochow University, Suzhou 215006, Jiangsu, China. This work was supported in part by the National Natural Science Foundation of China NSFC-12071332, NSFC-12371380, and Jiangsu Shuangchuang Project (JSSCTD202209). Email: [email protected]. Ren-Cang Li Department of Mathematics, University of Texas at Arlington, Arlington, TX 76019-0408, USA. Supported in part by NSF DMS-1719620 and DMS-2009689. Email: [email protected].

Abstract

Model reduction is a powerful tool in dealing with numerical simulation of large scale dynamic systems for studying complex physical systems. Two major types of model reduction methods for linear time-invariant dynamic systems are Krylov subspace-based methods and balanced truncation-based methods. The methods of the second type are much more theoretically sound than the first type in that there is a fairly tight global error bound on the approximation error between the original system and the reduced one. It is noted that the error bound is established based upon the availability of the exact controllability and observability Gramians. However, numerically, the Gramians are not available and have to be numerically calculated, and for a large scale system, a viable option is to compute low-rank approximations of the Gramians from which an approximate balanced truncation is then performed. Hence, rigorously speaking, the existing global error bound is not applicable to any reduced system obtained via approximate Gramians. The goal of this paper is to address this issue by establishing global error bounds for reduced systems via approximate balanced truncation.

Keywords: model reduction, balanced Truncation, transfer function, controllability Gramian, observability Gramian, Hankel singular value, low-rank approximation, error bound

Mathematics Subject Classification 78M34, 93A15, 93B40

1 Introduction

Model reduction is a powerful tool in dealing with numerical simulation of large scale dynamic systems for studying complex physical systems [2, 11, 17]. In this paper, we are interested in the following continuous linear time-invariant dynamic system


$\displaystyle\boldsymbol{x}^{\prime}(t)$	$\displaystyle=A\,\boldsymbol{x}(t)+B\,\boldsymbol{u}(t),\quad\mbox{given $\boldsymbol{x}(0)=\boldsymbol{x}_{0}$},$	(1.1a)
$\displaystyle\boldsymbol{y}(t)$	$\displaystyle=C^{\operatorname{T}}\,\boldsymbol{x}(t),$	(1.1b)

where $\boldsymbol{x}\,:\,t\in[0,\infty)\to\mathbb{R}^{n}$ is the state vector, and $\boldsymbol{u}\,:\,t\in[0,\infty)\to\mathbb{R}^{m}$ is the input, $\boldsymbol{y}\,:\,t\in[0,\infty)\to\mathbb{R}^{p}$ is the output, and $A\in\mathbb{R}^{n\times n},\,B\in\mathbb{R}^{n\times m},\,C\in\mathbb{R}^{n\times p}$ are constant matrices that define the dynamic system. In today’s applications of interests, such as very large scale integration (VLSI) circuit designs and structural dynamics, $n$ can be up to millions [3, 5, 11], but usually the dimensions of input and output vectors are much smaller, i.e., $p,\,m\ll n$ . Large $n$ can be an obstacle in practice both computationally and in memory usage. Model reduction is then called for to overcome the obstacle.

In a nutshell, model reduction for dynamic system (1.1) seeks two matrices $X,\,Y\in\mathbb{R}^{n\times r}$ such that $Y^{\operatorname{T}}X=I_{r}$ to reduce the system (1.1) to


$\displaystyle\hat{\boldsymbol{x}}_{r}^{\prime}(t)$	$\displaystyle=A_{r}\,\hat{\boldsymbol{x}}_{r}(t)+B_{r}\,\boldsymbol{u}(t),\quad\mbox{given $\hat{\boldsymbol{x}}_{r}(0)=Y^{\operatorname{T}}\boldsymbol{x}_{0}$},$	(1.2a)
$\displaystyle\boldsymbol{y}(t)$	$\displaystyle=C_{r}^{\operatorname{T}}\,\hat{\boldsymbol{x}}_{r}(t),$	(1.2b)

where $A_{r},\,B_{r},\,C_{r}$ are given by

A_{r}:=Y^{\operatorname{T}}AX\in\mathbb{R}^{r\times r},\quad B_{r}:=Y^{\operatorname{T}}B\in\mathbb{R}^{r\times m},\quad C_{r}:=X^{\operatorname{T}}C\in\mathbb{R}^{r\times p}.

(1.3)

Intuitively, this reduced system (1.2) may be thought of obtaining from (1.1) by letting $\boldsymbol{x}=X\hat{\boldsymbol{x}}_{r}$ and performing Galerkin projection with $Y$ . The new state vector $\hat{\boldsymbol{x}}_{r}$ is now in $\mathbb{R}^{r}$ , a much smaller space in dimension than $\mathbb{R}^{n}$ . In practice, for the reduced system to be of any use, the two systems must be “close” in some sense.

Different model reduction methods differ in their choosing $X$ and $Y$ , the projection matrices. There are two major types: Krylov subspace-based methods [5, 11, 17] and balanced truncation-based methods [2, 12]. The methods of the first type are computationally more efficient for large scale systems and reduced models are accurate around points where Krylov subspaces are built, while those of the second type are theoretically sound in that fairly tight global error bounds are known but numerically much more expensive in that controllability and observability Gramians which are provably positive definite have to be computed at cost of $O(n^{3})$ complexity.

Modern balanced truncation-based methods have improved, thanks to the discovery that the Gramians are usually numerically low-rank [4, 20, 21] and methods that compute their low-rank approximations in the factor form [15, 9]. The low-rank factors are then naturally used to compute an approximate balanced truncation. Moments ago, we pointed out the advantage of balanced truncation-based methods in their sound global approximations guaranteed by tight global error bounds, but these bounds are established based on exact Gramians and hence the exiting global error bounds, though suggestive, are no longer valid. To the best of our knowledge, there is no study as to the quality of reduced models by modern balanced truncation-based methods that use the low-rank approximate Gramians. Our aim in this paper is to address the void.

The rest of this paper is organized as follows. Section 2 reviews the basics of balanced truncation methods. Section 3 explains approximate balanced truncation, when some low-rank approximations of controllability and observability Gramians, not the exact Gramians themselves, are available. In Section 4, we establish our main results to quantify the accuracy of the reduced model by approximate balanced reduction. We draw our conclusions and makes some remarks. Some preliminary material on subspaces of $\mathbb{R}^{n}$ and perturbation for Lyapunov equation are discussed in appendixes.

Notation. $\mathbb{R}^{m\times n}$ is the set of $m\times n$ real matrices, $\mathbb{R}^{n}=\mathbb{R}^{n\times 1}$ , and $\mathbb{R}=\mathbb{R}^{1}$ . $I_{n}\in\mathbb{R}^{n\times n}$ is the identity matrix or simply $I$ if its size is clear from the context, and $\boldsymbol{e}_{j}$ is the $j$ th column of $I$ of apt size. $B^{\operatorname{T}}$ stands for the transpose of a matrix/vector. ${\cal R}(B)$ is the column subspace of $B$ , spanned by its columns. For $B\in\mathbb{R}^{m\times n}$ , its singular values are

\sigma_{1}(B)\geq\sigma_{2}(B)\geq\cdots\geq\sigma_{k}(B)\geq 0,

where $k=\min\{m,n\}$ , and $\sigma_{\max}(B)=\sigma_{1}(B)$ and $\sigma_{\min}(B)=\sigma_{k}(B)$ . $\|B\|_{2}$ , $\|B\|_{\operatorname{F}}$ , and $\|B\|_{\operatorname{ui}}$ are its spectral and Frobenius norms:

\|B\|_{2}=\sigma_{1}(B),\,\,\|B\|_{\operatorname{F}}=\Big{(}\sum_{i=1}^{k}[\sigma_{i}(B)]^{2}\Big{)}^{1/2},\,\,

respectively. $\|B\|_{\operatorname{ui}}$ is some unitarily invariant norm of $B$ [18, 23]. For a matrix $A\in\mathbb{R}^{n\times n}$ that is known to have real eigenvalues only, $\operatorname{eig}(A)=\{\lambda_{i}(A)\}_{i=1}^{n}$ denotes the set of its eigenvalues (counted by multiplicities) arranged in the decreasing order, and $\lambda_{\max}(A)=\lambda_{1}(A)$ and $\lambda_{\min}(A)=\lambda_{n}(A)$ . $A\succ 0\,(\succeq 0)$ means that it is symmetric and positive definite (semi-definite), and accordingly $A\prec 0\,(\preceq 0)$ if $-A\succ 0\,(\succeq 0)$ . MATLAB-like notation is used to access the entries of a matrix: $X_{(i:j,k:\ell)}$ to denote the submatrix of a matrix $X$ , consisting of the intersections of rows $i$ to $j$ and columns $k$ to $\ell$ , and when $i:j$ is replaced by $:$ , it means all rows, similarly for columns.

2 Review of balanced truncation

In this section, we will review the balanced truncation, minimally to the point to serve our purpose in this paper. The reader is referred to [2] for a more detailed exposition.

Consider continuous linear time-invariant dynamic system (1.1) and suppose that it is stable, observable and controllable [1, 25].

2.1 Quality of a reduced order model

Suppose initially $\boldsymbol{x}_{0}=0$ . Applying the Laplacian transformation to (1.1) yields

\boldsymbol{Y}(s)=\underbrace{C^{\operatorname{T}}(sI_{n}-A)^{-1}B}_{:=H(s)}\boldsymbol{U}(s),\quad s\in\mathbb{C},

where $\boldsymbol{U}(s)$ and $\boldsymbol{Y}(s)$ are the Laplacian transformations of $\boldsymbol{u}$ and $\boldsymbol{y}$ , respectively, and $H(s)\in\mathbb{C}^{p\times m}$ is the so-called transfer function of system (1.1). Conveniently, we will adopt the notation to denote the system (1.1) symbolically by

\mathscr{S}=\left(\begin{array}[]{c|c}A&B\\ \hline\cr\vphantom{\vrule height=12.80365pt,width=0.9pt,depth=2.84544pt}C^{\operatorname{T}}&\end{array}\right)

with the round bracket to distinguish it from the square bracket for matrices. The infinity Hankel norm of the system $\mathscr{S}$ , also known as the infinity Hankel norm of $H(\cdot)$ , is defined as

\|\mathscr{S}\|_{{\cal H}_{\infty}}=\|H(\cdot)\|_{{\cal H}_{\infty}}:=\sup_{\omega\in\mathbb{R}}\|H(\iota\omega)\|_{2}=\sup_{\omega\in\mathbb{R}}\sigma_{\max}(H(\iota\omega)),

(2.1)

where $\|\cdot\|_{2}$ is the spectral norm of a matrix, and $\iota=\sqrt{-1}$ is the imaginary unit.

In Section 1, we introduced the framework of model reduction with two matrices $X,\,Y\in\mathbb{R}^{n\times r}$ such that $Y^{\operatorname{T}}X=I_{r}$ . For the ease of our presentation going forward, we shall rename them as $X_{1},\,Y_{1}\in\mathbb{R}^{n\times r}$ and $Y_{1}^{\operatorname{T}}X_{1}=I_{r}$ . Next we look for $X_{2},\,Y_{2}\in\mathbb{R}^{n\times(n-r)}$ such that

I_{n}=[Y_{1},Y_{2}]^{\operatorname{T}}[X_{1},X_{2}]=\begin{bmatrix}Y_{1}^{\operatorname{T}}\\ Y_{2}^{\operatorname{T}}\end{bmatrix}[X_{1},X_{2}]=\begin{bmatrix}Y_{1}^{\operatorname{T}}X_{1}&Y_{1}^{\operatorname{T}}X_{2}\\ Y_{2}^{\operatorname{T}}X_{1}&Y_{2}^{\operatorname{T}}X_{2}\end{bmatrix}.

(2.2)

Such $X_{2},\,Y_{2}\in\mathbb{R}^{n\times(n-r)}$ always exist by Lemmas A.1 and A.2 if $\|\sin\Theta({\cal R}(X_{1}),{\cal R}(Y_{1}))\|_{2}<1$ . In any practical model reduction method, only $X_{1},\,Y_{1}$ need to be produced. That $X_{2},\,Y_{2}$ are introduced here is only for our analysis later. Denote by

T=[Y_{1},Y_{2}]^{\operatorname{T}},\quad T^{-1}=[X_{1},X_{2}],

(2.3)

which are consistent because of (2.2). To the original system (1.1), perform transformation: $\boldsymbol{x}(t)=T\hat{\boldsymbol{x}}(t)$ , to get


	$\displaystyle\hat{\boldsymbol{x}}^{\prime}(t)$	$\displaystyle=\widehat{A}\,\hat{\boldsymbol{x}}(t)+\widehat{B}\,\boldsymbol{u}(t),\quad\hat{\boldsymbol{x}}(0)=T^{-1}\boldsymbol{x}_{0},$	(2.4a)
	$\displaystyle\boldsymbol{y}(t)$	$\displaystyle=\widehat{C}^{\operatorname{T}}\,\hat{\boldsymbol{x}}(t),$	(2.4b)
where
	$\widehat{A}=TAT^{-1},\quad\widehat{B}=TB,\quad\widehat{C}=T^{-\operatorname{T}}C,$		(2.4c)

naturally partitioned as

\widehat{A}=\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle n-r\\\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\widehat{A}_{11}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\widehat{A}_{12}\\\scriptstyle n-r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\widehat{A}_{21}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\widehat{A}_{22}$\hfil\kern 5.0pt\crcr}}}}\right]$}},\quad\widehat{B}=\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle\scriptstyle\\\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\widehat{B}_{1}\\\scriptstyle n-r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\widehat{B}_{2}$\hfil\kern 5.0pt\crcr}}}}\right]$}},\quad\widehat{C}=\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle\scriptstyle\\\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\widehat{C}_{1}\\\scriptstyle n-r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\widehat{C}_{2}$\hfil\kern 5.0pt\crcr}}}}\right]$}}.

In particular,

\widehat{A}_{11}:=Y_{1}^{\operatorname{T}}AX_{1}\in\mathbb{R}^{r\times r},\quad\widehat{B}_{1}:=Y_{1}^{\operatorname{T}}B\in\mathbb{R}^{r\times m},\quad\widehat{C}_{1}:=X_{1}^{\operatorname{T}}C\in\mathbb{R}^{r\times p}.

(2.5)

One can verify that the transfer functions of (1.1) and (2.4) are exactly the same.

In current notation, the reduced system (1.2) in Section 1 takes the form


$\displaystyle\hat{\boldsymbol{x}}_{r}^{\prime}(t)$	$\displaystyle=\widehat{A}_{11}\,\hat{\boldsymbol{x}}_{r}(t)+\widehat{B}_{1}\,\boldsymbol{u}(t),\quad\hat{\boldsymbol{x}}_{r}(0)=Y_{1}^{\operatorname{T}}\boldsymbol{x}_{0},$	(2.6a)
$\displaystyle\boldsymbol{y}(t)$	$\displaystyle=\widehat{C}_{1}^{\operatorname{T}}\,\hat{\boldsymbol{x}}_{r}(t),$	(2.6b)

which will be denoted in short by $\mathscr{S}_{\operatorname{rd}}=\left(\begin{array}[]{c|c}\widehat{A}_{11}&\widehat{B}_{1}\\ \hline\cr\vphantom{\vrule height=12.80365pt,width=0.9pt,depth=2.84544pt}\widehat{C}_{1}^{\operatorname{T}}&\end{array}\right)$ . Its transfer function is given by

H_{\operatorname{rd}}(s)=\widehat{C}_{1}^{\operatorname{T}}(sI-\widehat{A}_{11})^{-1}\widehat{B}_{1},\quad s\in\mathbb{C}.

(2.7)

Naturally, we would like that the full system (2.4) and its reduced one (2.6) are “close”. One way to measure the closeness is the ${\cal H}_{\infty}$ -norm of the difference between the two transfer functions [2, 25]:

\|H(\cdot)-H_{\operatorname{rd}}(\cdot)\|_{{\cal H}_{\infty}}:=\sup_{\omega\in\mathbb{R}}\|H(\iota\omega)-H_{\operatorname{rd}}(\iota\omega)\|_{2},

assuming both systems are stable, observable and controllable [25]. Another way is by ${\cal H}_{2}$ -norm which we will get to later. It turns out that $H_{\operatorname{err}}(s)=H(s)-H_{\operatorname{rd}}(s)$ is the transfer function of an expanded dynamic system:


$\displaystyle\begin{bmatrix}\hat{\boldsymbol{x}}^{\prime}(t)\\ \hat{\boldsymbol{x}}_{r}^{\prime}(t)\end{bmatrix}$	$\displaystyle=\begin{bmatrix}\widehat{A}&\\ &\widehat{A}_{11}\end{bmatrix}\begin{bmatrix}\hat{\boldsymbol{x}}(t)\\ \hat{\boldsymbol{x}}_{r}(t)\end{bmatrix}+\begin{bmatrix}\widehat{B}\\ \widehat{B}_{1}\end{bmatrix}\boldsymbol{u}(t),\quad\begin{bmatrix}\hat{\boldsymbol{x}}(0)\\ \hat{\boldsymbol{x}}_{r}(0)\end{bmatrix}=\begin{bmatrix}T^{-1}\boldsymbol{x}_{0}\\ Y_{1}^{\operatorname{T}}\boldsymbol{x}_{0}\end{bmatrix},$	(2.8a)
$\displaystyle\hat{\boldsymbol{y}}(t)$	$\displaystyle=\begin{bmatrix}\widehat{C}\\ -\widehat{C}_{1}\end{bmatrix}^{\operatorname{T}}\begin{bmatrix}\hat{\boldsymbol{x}}(t)\\ \hat{\boldsymbol{x}}_{r}(t)\end{bmatrix},$	(2.8b)

or in the short notation

\mathscr{S}_{\operatorname{err}}=\left(\begin{array}[]{cc|c}\widehat{A}_{11}&&\widehat{B}\\ &\widehat{A}_{11}&\widehat{B}_{1}\\ \hline\cr\vphantom{\vrule height=12.80365pt,width=0.9pt,depth=2.84544pt}\widehat{C}^{\operatorname{T}}&-\widehat{C}_{1}^{\operatorname{T}}&\end{array}\right).

The key that really determines the quality of a reduced system is the subspaces ${\cal X}_{1}:={\cal R}(X_{1})$ and ${\cal Y}_{1}:={\cal R}(Y_{1})$ as far as the transfer function (2.7) is concerned, as guaranteed by the next theorem.

Theorem 2.1.

Given the subspaces ${\cal X}_{1}$ and ${\cal Y}_{1}$ of dimension $r$ such that $\|\sin\Theta({\cal X}_{1},{\cal Y}_{1})\|_{2}<1$ , any realizations of their basis matrices $X_{1},\,Y_{1}\in\mathbb{R}^{n\times r}$ satisfying $Y_{1}^{\operatorname{T}}X_{1}=I_{r}$ , respectively, do not affect the transfer function (2.7) of reduced system (2.6).

Proof.

Fix a pair of basis matrices $X_{1},\,Y_{1}\in\mathbb{R}^{n\times r}$ of ${\cal X}_{1}$ and ${\cal Y}_{1}$ , respectively, such that $Y_{1}^{\operatorname{T}}X_{1}=I_{r}$ . Consider any other two basis matrices $\check{X}_{1},\,\check{Y}_{1}\in\mathbb{R}^{n\times r}$ of ${\cal X}_{1}$ and ${\cal Y}_{1}$ , respectively, such that $\check{Y}_{1}^{\operatorname{T}}\check{X}_{1}=I_{r}$ . Then $\check{X}_{1}=X_{1}Z$ and $\check{Y}_{1}=Y_{1}W$ for some nonsingular $Z,\,W\in\mathbb{R}^{r\times r}$ . We have

I_{r}=\check{Y}_{1}^{\operatorname{T}}\check{X}_{1}=(Y_{1}W)^{\operatorname{T}}(X_{1}Z)=W^{\operatorname{T}}(Y_{1}^{\operatorname{T}}X_{1})Z=W^{\operatorname{T}}Z,

implying $W^{\operatorname{T}}=Z^{-1}$ , and

\check{Y}_{1}^{\operatorname{T}}A\check{X}_{1}=Z^{-1}\big{(}Y_{1}^{\operatorname{T}}AX_{1}\big{)}Z,\,\,\check{Y}_{1}^{\operatorname{T}}B=Z^{-1}\big{(}Y_{1}^{\operatorname{T}}B\big{)},\,\,\check{X}_{1}^{\operatorname{T}}C=Z^{\operatorname{T}}\big{(}X_{1}^{\operatorname{T}}C\big{)}.

The transfer function associated with $\check{X}_{1},\,\check{Y}_{1}$ is

	$\displaystyle\big{(}\check{X}_{1}^{\operatorname{T}}C\big{)}^{\operatorname{T}}\big{(}sI_{r}-\check{Y}_{1}^{\operatorname{T}}A\check{X}_{1}\big{)}^{-1}\big{(}\check{Y}_{1}^{\operatorname{T}}B\big{)}$	$\displaystyle=\big{(}X_{1}^{\operatorname{T}}C\big{)}^{\operatorname{T}}Z\big{[}Z^{-1}\big{(}sI_{r}-Y_{1}^{\operatorname{T}}AX_{1}\big{)}Z\big{]}^{-1}Z^{-1}\big{(}Y_{1}^{\operatorname{T}}B\big{)}$
		$\displaystyle=\big{(}X_{1}^{\operatorname{T}}C\big{)}^{\operatorname{T}}\big{(}sI_{r}-Y_{1}^{\operatorname{T}}AX_{1}\big{)}^{-1}\big{(}Y_{1}^{\operatorname{T}}B\big{)},$

having nothing to do with $Z$ and $W$ , as was to be shown. ∎

2.2 Balanced truncation

Balanced truncation fits into the general framework of model reduction, and thus it suffices for us to define $X_{1},\,Y_{1}\in\mathbb{R}^{n\times r}$ and $Y_{1}^{\operatorname{T}}X_{1}=I_{r}$ for balanced truncation accordingly.

The controllability and observability Gramians $P$ and $Q$ are defined as the solutions to the Lyapunov equations:


$\displaystyle AP+PA^{\operatorname{T}}$	$\displaystyle+BB^{\operatorname{T}}=0,$	(2.9a)
$\displaystyle A^{\operatorname{T}}Q+QA$	$\displaystyle+CC^{\operatorname{T}}=0,$	(2.9b)

respectively. Under the assumption that dynamic system (1.1) is stable, observable and controllable, the Lyapunov equations have unique solutions that are positive definite, i.e., $P\succ 0$ and $Q\succ 0$ . The model order reduction based on balanced truncation [2, 10] starts with a balanced transformation to dynamic system (1.1) such that both Gramians are the same and diagonal with diagonal entries being the system’s invariants, known as the Hankel singular values of the system.

Balanced truncation is classically introduced in the literature through some full-rank decompositions of $P$ and $Q$ :

P=SS^{\operatorname{T}}\quad\mbox{and}\quad Q=RR^{\operatorname{T}},

(2.10)

where $S,\,R\in\mathbb{R}^{n\times n}$ and are nonsingular because $P\succ 0$ and $Q\succ 0$ . But that is not necessary in theory, namely $S,\,R$ do not have to be square, in which case both will have no fewer than $n$ columns because the equalities in (2.10) ensure $\operatorname{rank}(S)=\operatorname{rank}(P)$ and $\operatorname{rank}(R)=\operatorname{rank}(Q)$ . Later in Theorem 2.3, we will show that balanced truncation is invariant with respect to how the decompositions in (2.10) are done, including non-square $S$ and $R$ . Such an invariance property is critical to our analysis.

Suppose that we have (2.10) with

S\in\mathbb{R}^{n\times m_{1}}\quad\mbox{and}\quad R\in\mathbb{R}^{n\times m_{2}}.

(2.11)

Without loss of generality, we may assume $m_{1}\geq m_{2}$ . Let the SVD of $S^{\operatorname{T}}R\in\mathbb{R}^{m_{1}\times m_{2}}$ be

	$S^{\operatorname{T}}R=U\Sigma V^{\operatorname{T}}\equiv\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle m_{1}-r\\$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle U_{1}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle U_{2}$\hfil\kern 5.0pt\crcr}}}}\right]$}}\times\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle m_{2}-r\\\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\Sigma_{1}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\\\scriptstyle m_{1}-r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\Sigma_{2}$\hfil\kern 5.0pt\crcr}}}}\right]$}}\times\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\\\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle V_{1}^{\operatorname{T}}\\\scriptstyle m_{2}-r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle V_{2}^{\operatorname{T}}$\hfil\kern 5.0pt\crcr}}}}\right]$}},$	(2.12a)
where

	$\displaystyle\Sigma_{1}=\operatorname{diag}(\sigma_{1},\ldots,\sigma_{r}),\quad\Sigma_{2}=\begin{bmatrix}\operatorname{diag}(\sigma_{r+1},\ldots,\sigma_{m_{2}})\\ 0_{(m_{1}-m_{2})\times(m_{2}-r)}\end{bmatrix},$	(2.12b)
	$\displaystyle\sigma_{1}\geq\sigma_{2}\geq\cdots\geq\sigma_{m_{2}}\geq 0.$	(2.12c)

Only $\sigma_{i}$ for $1\leq i\leq n$ are positive and the rest are $0$ . Those $\sigma_{i}$ for $1\leq i\leq n$ are the so-called Hankel singular values of the system, and they are invariant with respect to different ways of decomposing $P$ and $Q$ in (2.10) with (2.11), and, in fact, they are the square roots of the eigenvalues of $PQ$ , which are real and positive. To see this, we note $\{\sigma_{i}^{2}\}$ are the eigenvalues of $(S^{\operatorname{T}}R)^{\operatorname{T}}(S^{\operatorname{T}}R)=R^{\operatorname{T}}SS^{\operatorname{T}}R=R^{\operatorname{T}}PR$ whose nonzero eigenvalues are the same as those of $PRR^{\operatorname{T}}=PQ$ .

Define

T=(\Sigma_{(1:n,1:n)})^{-1/2}V_{(:,1:n)}^{\operatorname{T}}R^{\operatorname{T}}.

(2.13)

It can be verified that $T^{-1}=SU_{(:,1:n)}(\Sigma_{(1:n,1:n)})^{-1/2}$ because

	$\displaystyle\big{[}(\Sigma_{(1:n,1:n)})^{-1/2}V_{(:,1:n)}^{\operatorname{T}}R^{\operatorname{T}}\big{]}$	$\displaystyle\big{[}SU_{(:,1:n)}(\Sigma_{(1:n,1:n)})^{-1/2}\big{]}$
		$\displaystyle=(\Sigma_{(1:n,1:n)})^{-1/2}V_{(:,1:n)}^{\operatorname{T}}(R^{\operatorname{T}}S)U_{(:,1:n)}(\Sigma_{(1:n,1:n)})^{-1/2}$
		$\displaystyle=(\Sigma_{(1:n,1:n)})^{-1/2}V_{(:,1:n)}^{\operatorname{T}}(V\Sigma U^{\operatorname{T}})U_{(:,1:n)}(\Sigma_{(1:n,1:n)})^{-1/2}$
		$\displaystyle=(\Sigma_{(1:n,1:n)})^{-1/2}\Sigma_{(1:n,1:n)}(\Sigma_{(1:n,1:n)})^{-1/2}$
		$\displaystyle=I_{n}.$

With $T$ and $T^{-1}$ , we define $\widehat{A}$ , $\widehat{B}$ , and $\widehat{C}$ according to (2.4c), and, as a result, the transformed system (2.4). In turn, we have $A=T^{-1}\widehat{A}T$ , $B=T^{-1}\widehat{B}$ , and $C=T^{\operatorname{T}}\widehat{C}$ . Plug these relations into (2.9) to get, after simple re-arrangements,


$\displaystyle\widehat{A}(TPT^{\operatorname{T}})+(TPT^{-1})\widehat{A}^{\operatorname{T}}$	$\displaystyle+\widehat{B}\widehat{B}^{\operatorname{T}}=0,$	(2.14a)
$\displaystyle\widehat{A}^{\operatorname{T}}(T^{-\operatorname{T}}QT^{-1})+(T^{-\operatorname{T}}QT^{-1})\widehat{A}$	$\displaystyle+\widehat{C}\widehat{C}^{\operatorname{T}}=0,$	(2.14b)

which are precisely the Lyapunov equations for the Gramians

\widehat{P}=TPT^{\operatorname{T}},\quad\widehat{Q}=T^{-\operatorname{T}}QT^{-1},

(2.15)

of the transformed system (2.4). With the help of (2.10), (2.12) and (2.13), it is not hard to verify that

\widehat{P}=\widehat{Q}=\Sigma_{(1:n,1:n)},

balancing out the Gramians.

Given integer $1\leq r\leq n$ (usually $r\ll n$ ), according to the partitions of $U$ , $\Sigma$ , and $V$ in (2.12), we write


$\displaystyle T^{-1}$	$\displaystyle=\Big{[}SU_{1}\Sigma_{1}^{-1/2},SU_{2}(\Sigma_{2})_{(1:n-r,1:n-r)}^{-1/2}\Big{]}=:[X_{1},X_{2}],$	(2.16a)
$\displaystyle T$	$\displaystyle=\begin{bmatrix}\Sigma_{1}^{-1/2}V_{1}^{\operatorname{T}}R^{\operatorname{T}}\\ (\Sigma_{2})_{(1:n-r,1:n-r)}^{-1/2}V_{2}^{\operatorname{T}}R^{\operatorname{T}}\end{bmatrix}=:\begin{bmatrix}Y_{1}^{\operatorname{T}}\\ Y_{2}^{\operatorname{T}}\end{bmatrix},$	(2.16b)

leading to the reduced system (2.6) in form but with newly defined $X_{1},\,Y_{1}\in\mathbb{R}^{n\times r}$ by (2.16). In the rest of this section, we will adopt the notations in Section 2.1 but with $X_{1},\,Y_{1}$ given by (2.16).

Balanced truncation as stated is a very expensive procedure that generates (2.6) computationally. The computations of $P$ and $Q$ fully costs $O(n^{3})$ each, by, e.g., the Bartels-Stewart algorithm [7], decompositions $P=SS^{\operatorname{T}}$ and $Q=RR^{\operatorname{T}}$ costs $O(n^{3})$ each, and so does computing SVD of $S^{\operatorname{T}}R$ , not to mention $O(n^{2})$ storage requirements. However, it is a well-understood method in that the associated reduced system (1.2) inherits most important system properties of the original system: being stable, observable and controllable, and also there is a global error bound that guarantees the overall quality of the reduced system.

In terms of Gramians, the ${\cal H}_{\infty}$ - and ${\cal H}_{2}$ -norms of $H(\cdot)$ previously defined in (2.1) are given by (e.g., [2, Section 5.4.2])

	$\displaystyle\\|H(\cdot)\\|_{{\cal H}_{\infty}}$	$\displaystyle=\sqrt{\lambda_{\max}(PQ)}=\sigma_{1},$
	$\displaystyle\\|H(\cdot)\\|_{{\cal H}_{2}}$	$\displaystyle=\sqrt{\operatorname{tr}(B^{\operatorname{T}}QB)}=\sqrt{\operatorname{tr}(C^{\operatorname{T}}PC)},$

where $\sigma_{1}$ is the largest Hankel singular value in (2.12c). We remark that the transformations on $P$ and $Q$ as in (2.15) for any nonsingular $T$ , not necessarily the one in (2.13), preserve eigenvalues of $PQ$ because

\widehat{P}\widehat{Q}=(TPT^{\operatorname{T}})(T^{-\operatorname{T}}QT^{-1})=T(PQ)T^{-1}.

For the ease of future reference, we will denote by $H_{\operatorname{bt}}(s)$ :

H_{\operatorname{bt}}(s):=\widehat{C}_{1}^{\operatorname{T}}(sI_{r}-\widehat{A}_{11})^{-1}\widehat{B}_{1},

(2.17)

the transfer function of the reduced system (2.6) with $X_{1},\,Y_{1}\in\mathbb{R}^{n\times r}$ as in (2.16) by the balanced truncation.

The next theorem is well-known; see, e.g., [2, Theorem 7.9], [25, Theorem 8.16].

Theorem 2.2 ([2, 25]).

For $X_{1}$ and $Y_{1}$ from the balanced truncation as in (2.16), we have

\sigma_{r+1}\leq\|H(\cdot)-H_{\operatorname{bt}}(\cdot)\|_{{\cal H}_{\infty}}\leq 2\sum_{j=r+1}^{n}\sigma_{j},

(2.18)

where $\sigma_{1}\geq\sigma_{2}\geq\cdots\geq\sigma_{n}$ are the Hankel singular values of the system, i.e., the first $n$ singular values of $S^{\operatorname{T}}R$ .

Remark 2.1.

The left inequality in (2.18) actually holds for any reduced system of order $r$ , not necessarily from balanced truncation. In fact, it is known that (see e.g., [2, Proposition 8.3] and [25, Lemma 8.5])

\sigma_{r+1}\leq\|H(\cdot)-H_{\operatorname{rd}}(\cdot)\|_{{\cal H}_{\infty}},

where $H_{\operatorname{rd}}(s)$ is the transfer function (2.6) of reduced system (2.6) by any $X_{1},\,Y_{1}\in\mathbb{R}^{n\times r}$ such that $Y_{1}^{\operatorname{T}}X_{1}=I_{r}$ .

One thing that is not clear yet and hasn’t been drawn much attention in the literature is whether the reduced system by the balanced truncation of order $r$ varies with the decompositions $P=SS^{\operatorname{T}}$ and $Q=RR^{\operatorname{T}}$ which are not unique, including $S$ and $R$ that may not necessarily be square. This turns out to be an easy question to answer.

Theorem 2.3.

If $\sigma_{r}>\sigma_{r+1}$ , then the transfer function of the reduced system (2.6) by the balanced truncation of order $r$ is unique, regardless of any variations in the decompositions in (2.10).

Proof.

We will show that the projection matrices $X_{1}$ and $Y_{1}$ defined in (2.16) are invariant with respect to any choices of decompositions for $P$ and $Q$ of the said kind. Suppose we have two different decompositions for each one of $P$ and $Q$

	$\displaystyle P$	$\displaystyle=SS^{\operatorname{T}}=\check{S}\check{S}^{\operatorname{T}}\quad\mbox{with}\,\,S\in\mathbb{R}^{n\times n},\,\check{S}\in\mathbb{R}^{n\times\check{n}_{1}},$
	$\displaystyle Q$	$\displaystyle=RR^{\operatorname{T}}=\check{R}\check{R}^{\operatorname{T}}\quad\mbox{with}\,\,R\in\mathbb{R}^{n\times n},\,\check{R}\in\mathbb{R}^{n\times\check{n}_{2}}.$

The idea is to show that after fixing one pair of decompositions $P=SS^{\operatorname{T}}$ and $Q=RR^{\operatorname{T}}$ , $X_{1}$ and $Y_{1}$ constructed from any other decompositions $P=\check{S}\check{S}^{\operatorname{T}}$ and $Q=\check{R}\check{R}^{\operatorname{T}}$ , including nonsquare $\check{S}$ and $\check{R}$ , remain the same. Evidentally $\check{n}_{1},\,\check{n}_{2}\geq n$ .

Without loss of generality, we may assume $\check{n}_{1}\geq\check{n}_{2}$ ; otherwise we can append some columns of $0$ to $\check{S}$ from the right.

Since ${\cal R}(P)={\cal R}(S)={\cal R}(\check{S})$ and ${\cal R}(Q)={\cal R}(R)={\cal R}(\check{R})$ , there exist $W\in\mathbb{R}^{\check{n}_{1}\times n}$ and $Z\in\mathbb{R}^{\check{n}_{2}\times n}$ such that

\check{S}=SW,\quad\check{R}=RZ.

It can be verified that $WW^{\operatorname{T}}=I_{n}$ and $ZZ^{\operatorname{T}}=I_{n}$ , i.e., both $W\in\mathbb{R}^{n\times\check{n}_{1}}$ and $Z\in\mathbb{R}^{n\times\check{n}_{2}}$ have orthonormal rows. Suppose we already have the SVD of $S^{\operatorname{T}}R$ as in (2.12) with $m_{1}=m_{2}=n$ . Both $W^{\operatorname{T}}U\in\mathbb{R}^{\check{n}_{1}\times n}$ and $Z^{\operatorname{T}}V\in\mathbb{R}^{\check{n}_{1}\times n}$ have orthonormal columns. There exist $\check{U}_{3}\in\mathbb{R}^{\check{n}_{1}\times(\check{n}_{1}-n)}$ and $\check{V}_{3}\in\mathbb{R}^{\check{n}_{2}\times(\check{n}_{2}-n)}$ such that

[W^{\operatorname{T}}U,\check{U}_{3}]\in\mathbb{R}^{\check{n}_{1}\times\check{n}_{1}}\quad\mbox{and}\quad[Z^{\operatorname{T}}V,\check{V}_{3}]\in\mathbb{R}^{\check{n}_{1}\times\check{n}_{1}}

are orthogonal matrices. We have

	$\displaystyle\check{S}^{\operatorname{T}}\check{R}=W^{\operatorname{T}}S^{\operatorname{T}}RZ$	$\displaystyle=W^{\operatorname{T}}[U_{1},U_{2}]\begin{bmatrix}\Sigma_{1}&\\ &\Sigma_{2}\end{bmatrix}\begin{bmatrix}V_{1}^{\operatorname{T}}\\ V_{2}^{\operatorname{T}}\end{bmatrix}Z$
		$\displaystyle=[W^{\operatorname{T}}U_{1},W^{\operatorname{T}}U_{2},\check{U}_{3}]\begin{bmatrix}\Sigma_{1}&&\\ &\Sigma_{2}&\\ &&0_{(\check{n}_{1}-n)\times(\check{n}_{2}-n)}\end{bmatrix}\begin{bmatrix}(Z^{\operatorname{T}}V_{1})^{\operatorname{T}}\\ (Z^{\operatorname{T}}V_{2})^{\operatorname{T}}\\ \check{V}_{3}^{\operatorname{T}}\end{bmatrix},$

yielding an SVD of $\check{S}^{\operatorname{T}}\check{R}$ , for which the corresponding projection matrices from $P=\check{S}\check{S}^{\operatorname{T}}$ and $Q=\check{R}\check{R}^{\operatorname{T}}$ are given by

\check{S}(W^{\operatorname{T}}U_{1})\Sigma_{1}^{-1/2}=SW(W^{\operatorname{T}}U_{1})\Sigma_{1}^{-1/2}=SU_{1}\Sigma_{1}^{-1/2}

and, similarly, $\check{R}(Z^{\operatorname{T}}V_{1})\Sigma_{1}^{-1/2}=RV_{1}\Sigma_{1}^{-1/2}$ , yielding the same projection matrices as $X_{1}$ and $Y_{1}$ in (2.16) from $P=SS^{\operatorname{T}}$ and $Q=RR^{\operatorname{T}}$ , which in turn leads to the same reduced system (2.6) and hence the same transfer function. Now let

\check{S}^{\operatorname{T}}\check{R}=[\check{U}_{1},\check{U}_{2},\check{U}_{3}]\begin{bmatrix}\Sigma_{1}&&\\ &\Sigma_{2}&\\ &&0_{(\check{n}_{1}-n)\times(\check{n}_{2}-n)}\end{bmatrix}\begin{bmatrix}\check{V}_{1}^{\operatorname{T}}\\ \check{V}_{2}^{\operatorname{T}}\\ \check{V}_{3}^{\operatorname{T}}\end{bmatrix}

be another SVD of $\check{S}^{\operatorname{T}}\check{R}$ subject to the inherent freedom in SVD, where $\check{U}_{1}\in\mathbb{R}^{\check{n}_{1}\times r}$ and $\check{V}_{1}\in\mathbb{R}^{\check{n}_{2}\times r}$ . Since $\sigma_{r}>\sigma_{r+1}$ , by the uniqueness of singular subspaces, we know ${\cal R}(\check{U}_{1})={\cal R}(W^{\operatorname{T}}U_{1})$ and ${\cal R}(\check{V}_{1})={\cal R}(Z^{\operatorname{T}}V_{1})$ . Therefore

{\cal R}(\check{S}\check{U}_{1}\Sigma_{1}^{-1/2})={\cal R}(\check{S}\check{U}_{1})={\cal R}(\check{S}W^{\operatorname{T}}U_{1})={\cal R}(\check{S}W^{\operatorname{T}}U_{1}\Sigma_{1}^{-1/2})={\cal R}(SU_{1}\Sigma_{1}^{-1/2}),

and similarly, ${\cal R}(\check{R}\check{V}_{1}\Sigma_{1}^{-1/2})={\cal R}(RV_{1}\Sigma_{1}^{-1/2})$ , implying the same transfer function regardless of whether the reduced system is obtained by the projection matrix pair $(X_{1},Y_{1})$ or by the pair $(\check{S}\check{U}_{1}\Sigma_{1}^{-1/2},\check{R}\check{V}_{1}\Sigma_{1}^{-1/2})$ by Theorem 2.1. ∎

2.3 A variant of balanced truncation

A distinguished feature of the transformation $T$ in (2.16) is that it makes the transformed system (2.4) balanced, i.e., both controllability and observability Gramians are the same and diagonal, and so the reduced system (2.6) is balanced, too. But as far as just the transfer function of the reduced system is concerned, there is no need to have $X_{1}$ and $Y_{1}$ precisely the same as the ones in (2.16) because of Theorem 2.1. In fact, all that we need is to make sure ${\cal R}(X_{1})={\cal R}(SU_{1})$ and ${\cal R}(Y_{1})={\cal R}(RV_{1})$ , besides $Y_{1}^{\operatorname{T}}X_{1}=I_{r}$ . Specifically, we have by Theorem 2.1

Corollary 2.1.

Let $\check{X}_{1},\,\check{Y}_{1}\in\mathbb{R}^{n\times r}$ such that

{\cal R}(\check{X}_{1})={\cal R}(SU_{1}),\,\,{\cal R}(\check{Y}_{1})={\cal R}(RV_{1}),\,\,\check{Y}_{1}^{\operatorname{T}}\check{X}_{1}=I_{r}.

(2.19)

Then $H_{\operatorname{bt}}(s)\equiv\big{(}\check{X}_{1}^{\operatorname{T}}C\big{)}^{\operatorname{T}}\big{(}sI_{r}-\check{Y}_{1}^{\operatorname{T}}A\check{X}_{1}\big{)}^{-1}\big{(}\check{Y}_{1}^{\operatorname{T}}B\big{)}$ , i.e., the reduced system (2.6) with (2.5) obtained by replacing $X_{1},\,Y_{1}$ from (2.16) with $\check{X}_{1},\,\check{Y}_{1}$ satisfying $\check{Y}_{1}^{\operatorname{T}}\check{X}_{1}=I_{r}$ has the same transfer function as the one from the true balanced truncation.

$X_{1},\,Y_{1}\in\mathbb{R}^{n\times r}$ defined by (2.16) for balanced truncation are difficult to work with in analyzing the quality of balanced truncation. Luckily, the use of transfer function for analysis allows us to focus on the subspaces ${\cal R}(X_{1})$ and ${\cal R}(Y_{1})$ . Later, instead of the concrete forms of $X_{1}$ and $Y_{1}$ in (2.16), we will work with the reduced system (2.6) with

X_{1}=SU_{1},\,\,Y_{1}=RV_{1}\Sigma_{1}^{-1}.

(2.20)

It is not hard to verify that ${\cal R}(X_{1})={\cal R}(SU_{1})$ , ${\cal R}(Y_{1})={\cal R}(RV_{1})$ , and $Y_{1}^{\operatorname{T}}X_{1}=I_{r}$ . Effectively, in the notations of Section 2.2 up to SVD (2.12), this relates to transform the original system (1.1) to (2.4) with

T^{-1}=\Big{[}SU_{1},SU_{2}\Big{]},\quad T=\begin{bmatrix}\Sigma_{1}^{-1}V_{1}^{\operatorname{T}}R^{\operatorname{T}}\\ (\Sigma_{2})_{(1:n-r,1:n-r)}^{-1}V_{2}^{\operatorname{T}}R^{\operatorname{T}}\end{bmatrix}.

(2.21)

Accordingly, the Gramians for the reduced system, by (2.15), are

\widehat{P}=TPT^{\operatorname{T}}=I_{n},\quad\widehat{Q}=T^{-\operatorname{T}}QT^{-1}=\Sigma_{(1:n,1:n)}^{2},

(2.22)

which are not balanced, but the reduced system has the same transfer function as by the balanced truncation with (2.16) nonetheless.

3 Approximate balanced truncation

When $n$ is large, balanced truncation as stated is a very expensive procedure both computationally and in storage usage. Fortunately, $P$ and $Q$ are usually numerically low-rank [8, 20, 6, 21, 4], which means, $P$ and $Q$ can be very well approximated by $\widetilde{P}=\widetilde{S}\widetilde{S}^{\operatorname{T}}$ and $\widetilde{Q}=\widetilde{R}\widetilde{R}^{\operatorname{T}}$ , respectively, where $\widetilde{S}\in\mathbb{R}^{n\times\widetilde{r}_{1}}$ and $\widetilde{R}\in\mathbb{R}^{n\times\widetilde{r}_{2}}$ with $\widetilde{r}_{1},\,\widetilde{r}_{2}\ll n$ . Naturally, we will use $\widetilde{S}$ and $\widetilde{R}$ to play the roles of $S$ and $R$ in Section 2.2. Specifically, a model order reduction by approximate balanced truncation goes as follows.

1.

compute some low-rank approximations to $P$ and $Q$ in the product form

$P\approx\widetilde{P}=\widetilde{S}\widetilde{S}^{\operatorname{T}},\quad Q\approx\widetilde{Q}=\widetilde{R}\widetilde{R}^{\operatorname{T}},$ (3.1)

where $\widetilde{S}\in\mathbb{R}^{n\times\widetilde{r}_{1}}$ and $\widetilde{R}\in\mathbb{R}^{n\times\widetilde{r}_{2}}$ . Without loss of generality, assume $\widetilde{r}_{1}\geq\widetilde{r}_{2}$ , for our presentation.

compute SVD

\widetilde{S}^{\operatorname{T}}\widetilde{R}=\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle\widetilde{r}_{1}-r\\$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\widetilde{U}_{1}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\widetilde{U}_{2}$\hfil\kern 5.0pt\crcr}}}}\right]$}}\times\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle\widetilde{r}_{2}-r\\\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\widetilde{\Sigma}_{1}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\\\scriptstyle\widetilde{r}_{1}-r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\widetilde{\Sigma}_{2}$\hfil\kern 5.0pt\crcr}}}}\right]$}}\times\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\\\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\widetilde{V}_{1}^{\operatorname{T}}\\\scriptstyle\widetilde{r}_{2}-r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\widetilde{V}_{2}^{\operatorname{T}}$\hfil\kern 5.0pt\crcr}}}}\right]$}},

where $\widetilde{\Sigma}_{1}=\operatorname{diag}(\widetilde{\sigma}_{1},\widetilde{\sigma}_{2},\ldots,\widetilde{\sigma}_{r})$ , and $\widetilde{\Sigma}_{2}=\begin{bmatrix}\operatorname{diag}(\widetilde{\sigma}_{r+1},\widetilde{\sigma}_{2},\ldots,\widetilde{\sigma}_{\widetilde{r}_{2}})\\ 0_{(\widetilde{r}_{1}-\widetilde{r}_{2})\times(\widetilde{r}_{2}-r)}\\ \end{bmatrix}$ with these $\widetilde{\sigma}_{i}$ arranged in the decreasing order, as in (2.12c) for $\sigma_{i}$ .

finally, $A$ , $B$ , and $C$ are reduced to

\widetilde{A}_{11}:=\widetilde{Y}_{1}^{\operatorname{T}}A\widetilde{X}_{1},\quad\widetilde{B}_{1}:=\widetilde{Y}_{1}^{\operatorname{T}}B,\quad\widetilde{C}_{1}:=\widetilde{X}_{1}^{\operatorname{T}}C,

(3.2)

where

\widetilde{X}_{1}=\widetilde{S}\widetilde{U}_{1}\widetilde{\Sigma}_{1}^{-1/2},\quad\widetilde{Y}_{1}=\widetilde{R}\widetilde{V}_{1}\widetilde{\Sigma}_{1}^{-1/2}.

(3.3)

It can be verified that $\widetilde{Y}^{\operatorname{T}}\widetilde{X}=I_{r}$ . Accordingly, we will have a reduced system


$\displaystyle\widetilde{\boldsymbol{x}}_{r}^{\prime}(t)$	$\displaystyle=\widetilde{A}_{11}\,\widetilde{\boldsymbol{x}}_{r}(t)+\widetilde{B}_{1}\,\boldsymbol{u}(t),\quad\mbox{given $\widetilde{\boldsymbol{x}}_{r}(0)=\widetilde{Y}^{\operatorname{T}}\boldsymbol{x}_{0}$},$	(3.4a)
$\displaystyle\widetilde{\boldsymbol{y}}(t)$	$\displaystyle=\widetilde{C}_{1}^{\operatorname{T}}\,\widetilde{\boldsymbol{x}}_{r}(t),$	(3.4b)

which will not be quite the same as (1.2) with $X_{1}$ and $Y_{1}$ in (2.16) from the (exact) balanced truncation. The transfer function of (3.4) is

\widetilde{H}_{\operatorname{bt}}(s)=\widetilde{C}_{1}^{\operatorname{T}}(sI-\widetilde{A}_{11})^{-1}\widetilde{B}_{1},\quad s\in\mathbb{C}.

(3.5)

One lingering question that has not been addressed in the literature is how good reduced system (3.4) is, compared to the true reduced system of balanced truncation. The seemingly convincing argument that if $\widetilde{P}=\widetilde{S}\widetilde{S}^{\operatorname{T}}$ and $\widetilde{Q}=\widetilde{R}\widetilde{R}^{\operatorname{T}}$ are sufficiently accurate then $\widetilde{S}^{\operatorname{T}}\widetilde{R}$ should approximate $S^{\operatorname{T}}R$ well could be doubtful because usually $\widetilde{r}_{1},\,\widetilde{r}_{2}\ll n$ . A different argument may say otherwise. In order for $\widetilde{P}=\widetilde{S}\widetilde{S}^{\operatorname{T}}$ and $\widetilde{Q}=\widetilde{R}\widetilde{R}^{\operatorname{T}}$ to approximate $P$ and $Q$ well, respectively, both $\widetilde{S}$ and $\widetilde{R}$ must approximate the dominant components of the factors $S$ and $R$ of $P$ and $Q$ well. The problem is $\widetilde{r}_{1},\,\widetilde{r}_{2}\ll n$ here while it is possible that the dominant components of $S$ and $R$ could mismatch in forming $S^{\operatorname{T}}R$ , i.e., in the unlucky scenario, the dominant components of $S$ match the least dominant components of $R$ in forming $S^{\operatorname{T}}R$ and simply extracting out the dominant components of $S$ and $R$ is not enough. Hence it becomes critically important to provide theoretical analysis that shows the quality of approximate balanced truncation derived from $\widetilde{P}=\widetilde{S}\widetilde{S}^{\operatorname{T}}$ and $\widetilde{Q}=\widetilde{R}\widetilde{R}^{\operatorname{T}}$ , assuming $\|P-\widetilde{P}\|$ and $\|Q-\widetilde{Q}\|$ are tiny.

By the same reasoning as we argue in Subsection 2.2, the transfer function $\widetilde{H}_{\operatorname{bt}}(\cdot)$ stays the same for any $\widetilde{X}_{1},\,\widetilde{Y}_{1}\in\mathbb{R}^{n\times r}$ that satisfy

{\cal R}(\widetilde{X}_{1})={\cal R}(\widetilde{S}\widetilde{U}_{1}),\quad{\cal R}(\widetilde{Y}_{1})={\cal R}(\widetilde{R}\widetilde{V}_{1})\quad\mbox{such that}\quad\widetilde{Y}^{\operatorname{T}}\widetilde{X}=I_{r},

(3.6)

and the pair $(\widetilde{X}_{1},\widetilde{Y}_{1})$ in (3.3) is just one of many concrete pairs that satisfy (3.6). Again $\widetilde{X}_{1}$ and $\widetilde{Y}_{1}$ in (3.3) for approximate balanced truncation are difficult to work with in our later analysis. Luckily, we can again focus on the subspaces ${\cal R}(\widetilde{X}_{1})$ and ${\cal R}(\widetilde{Y}_{1})$ because of Theorem 2.1. Precisely what $\widetilde{X}_{1},\,\widetilde{Y}_{1}\in\mathbb{R}^{n\times r}$ to use will be specified later in Section 4 so that they will be close to $X_{1}$ and $Y_{1}$ in (2.20), respectively.

We reiterate our notations for the reduced models going forward.

•

$(\widehat{A}_{11},\widehat{B}_{1},\widehat{C}_{1})$ stands for the matrices for the reduced model (2.6) by balanced truncation with $X_{1}$ and $Y_{1}$ in (2.20). It is different from the one in the literature we introduced earlier with $X_{1}$ and $Y_{1}$ in (2.16), but both share the same transfer function denoted by $H_{\operatorname{bt}}(\cdot)$ .
•

$(\widetilde{A}_{11},\widetilde{B}_{1},\widetilde{C}_{1})$ stands for the matrices for the reduced model (3.4) by approximate balanced truncation with $\widetilde{X}_{1}$ and $\widetilde{Y}_{1}$ specified later in (4.21). It is different from the one in the literature we introduced earlier with $\widetilde{X}_{1}$ and $\widetilde{Y}_{1}$ in (3.3), but both share the same transfer function denoted by $\widetilde{H}_{\operatorname{bt}}(\cdot)$ .

In the rest of this paper, assuming $\big{\|}P-\widetilde{P}\big{\|}_{2},\,\big{\|}Q-\widetilde{Q}\big{\|}_{2}\leq\epsilon$ , we will

(i)

bound $\widehat{A}_{11}-\widetilde{A}_{11}$ , $\widehat{B}_{1}-\widetilde{B}_{1}$ , $\widehat{C}_{1}-\widetilde{C}_{1}$ in terms of $\epsilon$ , where $\widehat{A}_{11}$ , $\widehat{B}_{1}$ , and $\widehat{C}_{1}$ are from exact balanced truncation as in (2.5) with $X_{1},\,Y_{1}$ given by (2.20), while $\widetilde{A}_{11}$ , $\widetilde{B}_{1}$ , and $\widetilde{C}_{1}$ are from the approximate balanced truncation as in (3.2) with $\widetilde{X}_{1},\,\widetilde{Y}_{1}\in\mathbb{R}^{n\times r}$ to be specified;
(ii)

bound $\big{\|}H_{\operatorname{bt}}(\cdot)-\widetilde{H}_{\operatorname{bt}}(\cdot)\big{\|}$ and $\big{\|}H(\cdot)-\widetilde{H}_{\operatorname{bt}}(\cdot)\big{\|}$ in terms of $\epsilon$ for both $\|\cdot\|_{{\cal H}_{\infty}}$ and $\|\cdot\|_{{\cal H}_{2}}$ .

4 Quality of the approximate balanced reduction

The true balanced truncation requires computing the controllability and observability Gramians $P$ and $Q$ to the working precision, performing their full-rank decompositions (such as the Cholesky decomposition) and an SVD, each of which costs $O(n^{3})$ flops. It is infeasible for large scale dynamic systems. Luckily, the numbers of columns in $B$ and $C$ are usually of $O(1)$ and $P$ and $Q$ numerically have extremely low ranks. In practice, due to the fast decay of the Hankel singular values $\sigma_{i}$ [4, 6, 10, 8, 20], and the fact that solving the Lyapunov equations in (2.9) for the full Gramians is too expensive and storing the full Gramians takes too much space, we can only afford to compute low-rank approximations to $P$ and $Q$ in the product form as in (3.1) [14, 21, 22]. More than that, $\widetilde{P}$ and $\widetilde{Q}$ approach $P$ and $Q$ from below, i.e.,

	$0\preceq\widetilde{P}=\widetilde{S}\widetilde{S}^{\operatorname{T}}\preceq P,~{}~{}0\preceq\widetilde{Q}=\widetilde{R}\widetilde{R}^{\operatorname{T}}\preceq Q,$	(4.1a)
where $\widetilde{S}\in\mathbb{R}^{n\times\widetilde{r}_{1}}$ and $\widetilde{R}\in\mathbb{R}^{n\times\widetilde{r}_{2}}$ . This is what we will assume about $\widetilde{P}$ amd $\widetilde{Q}$ in the rest of this paper, besides
	$\\|P-\widetilde{P}\\|_{2}\leq\epsilon_{1},\quad\\|Q-\widetilde{Q}\\|_{2}\leq\epsilon_{2}$	(4.1b)

for some sufficiently tiny $\epsilon_{1}$ and $\epsilon_{2}$ . Except their existences, exactly what $P$ , $Q$ and their full-rank factors $S$ and $R$ are not needed in our analysis. Because of (4.1a), we may write


	$\displaystyle P$	$\displaystyle=\widetilde{P}+EE^{\operatorname{T}}=[\widetilde{S},E][\widetilde{S},E]^{\operatorname{T}}=SS^{\operatorname{T}},$	(4.2a)
	$\displaystyle Q$	$\displaystyle=\widetilde{Q}+FF^{\operatorname{T}}=[\widetilde{R},F][\widetilde{R},F]^{\operatorname{T}}=RR^{\operatorname{T}},$	(4.2b)
where $E\in\mathbb{R}^{n\times p_{1}}$ and $F\in\mathbb{R}^{n\times p_{2}}$ are unknown, and neither are
	$S=[\widetilde{S},E]\in\mathbb{R}^{n\times m_{1}},\quad R=[\widetilde{R},F]\in\mathbb{R}^{n\times m_{2}},$		(4.2c)

$m_{1}=\widetilde{r}_{1}+p_{1}$ and $m_{2}=\widetilde{r}_{2}+p_{2}$ . Without loss of generality, we may assume

m_{1}\geq m_{2};

otherwise, we simply append a few columns of $0$ to $E$ . Let


	$\displaystyle G$	$\displaystyle:=S^{\operatorname{T}}R=[\widetilde{S},E]^{\operatorname{T}}[\widetilde{R},F]=\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle\widetilde{r}_{2}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle p_{2}\\\scriptstyle\widetilde{r}_{1}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\widetilde{S}^{\operatorname{T}}\widetilde{R}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\widetilde{S}^{\operatorname{T}}F\\\scriptstyle p_{1}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle E^{\operatorname{T}}\widetilde{R}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle E^{\operatorname{T}}F$\hfil\kern 5.0pt\crcr}}}}\right]$}},\\\widetilde{G}$	$\displaystyle:=\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle\widetilde{r}_{2}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle p_{2}\\\scriptstyle\widetilde{r}_{1}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\widetilde{S}^{\operatorname{T}}\widetilde{R}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle 0\\\scriptstyle p_{1}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle 0$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle 0$\hfil\kern 5.0pt\crcr}}}}\right]$}}=G-\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle\widetilde{r}_{2}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle p_{2}\\\scriptstyle\widetilde{r}_{1}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle 0$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\widetilde{S}^{\operatorname{T}}F\\\scriptstyle p_{1}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle E^{\operatorname{T}}\widetilde{R}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle E^{\operatorname{T}}F$\hfil\kern 5.0pt\crcr}}}}\right]$}}.$		(4.3d)

It is reasonable to require

\widetilde{r}_{i}\geq r\quad\mbox{for $i=1,2$},

because we are looking for balanced truncation of order $r$ . Lemma 4.1 provides some basic inequalities we need in the rest of this paper.

Lemma 4.1.

Suppose that (4.1) holds. Then

	$\displaystyle\\|\widetilde{S}\\|_{2}=\sqrt{\\|\widetilde{P}\\|_{2}}$	$\displaystyle\leq\\|S\\|_{2}=\sqrt{\\|P\\|_{2}},$	$\displaystyle\quad\big{\\|}\widetilde{R}\big{\\|}_{2}=\sqrt{\\|\widetilde{Q}\\|_{2}}$	$\displaystyle\leq\\|R\\|_{2}=\sqrt{\\|Q\\|_{2}},$		(4.4)
	$\displaystyle\\|E\\|_{2}$	$\displaystyle\leq\sqrt{\epsilon_{1}},$	$\displaystyle\quad\\|F\\|_{2}$	$\displaystyle\leq\sqrt{\epsilon_{2}},$		(4.5)

and

	$\displaystyle\\|\widetilde{G}-G\\|_{2}$	$\displaystyle=\left\\|\begin{bmatrix}0&\widetilde{S}^{\operatorname{T}}F\\ E^{\operatorname{T}}\widetilde{R}&E^{\operatorname{T}}F\end{bmatrix}\right\\|_{2}$
		$\displaystyle\leq\max\left\{\sqrt{\\|P\\|_{2}\epsilon_{2}},\sqrt{\\|Q\\|_{2}\epsilon_{1}}\right\}+\sqrt{\epsilon_{1}\epsilon_{2}}=:\varepsilon.$		(4.6)

Proof.

We have $\|\widetilde{S}\|_{2}^{2}=\|\widetilde{S}\widetilde{S}^{\operatorname{T}}\|_{2}=\|\widetilde{P}\|_{2}\leq\|P\|_{2}=\|SS^{\operatorname{T}}\|_{2}=\|S\|_{2}^{2}$ , proving the first relation in (4.4). It follows from $P-\widetilde{P}=EE^{\operatorname{T}}$ that $\|P-\widetilde{P}\|_{2}=\|E\|_{2}^{2}$ , yielding the first inequality in (4.5) upon using (4.1). Similarly, we will have the second relation in (4.4) and the second inequality in (4.5).

For (4.6), we have

	$\displaystyle\left\\|\begin{bmatrix}0&\widetilde{S}^{\operatorname{T}}F\\ E^{\operatorname{T}}\widetilde{R}&E^{\operatorname{T}}F\end{bmatrix}\right\\|_{2}$	$\displaystyle\leq\left\\|\begin{bmatrix}0&\widetilde{S}^{\operatorname{T}}F\\ E^{\operatorname{T}}\widetilde{R}&0\end{bmatrix}\right\\|_{2}+\left\\|\begin{bmatrix}0&0\\ 0&E^{\operatorname{T}}F\end{bmatrix}\right\\|_{2}$
		$\displaystyle=\max\left\{\\|\widetilde{S}^{\operatorname{T}}F\\|_{2},\\|E^{\operatorname{T}}\widetilde{R}\\|_{2}\right\}+\\|E^{\operatorname{T}}F\\|_{2}$
		$\displaystyle\leq\max\left\{\sqrt{\\|P\\|_{2}\epsilon_{2}},\sqrt{\\|Q\\|_{2}\epsilon_{1}}\right\}+\sqrt{\epsilon_{1}\epsilon_{2}},$

as was to be shown. ∎

Remark 4.1.

Besides the spectral norm $\|\cdot\|_{2}$ , the Frobenius norm is another commonly used matrix norm, too. Naturally, we are wondering if we could have Frobenius-norm versions of (4.1b) and Lemma 4.1. Theoretically, it can be done, but there is one potential problem which is that matrix dimension $n$ will show up. Here is why:

\|E\|_{\operatorname{F}}^{2}\leq{\sqrt{\operatorname{rank}(E)}}\,\|EE^{\operatorname{T}}\|_{\operatorname{F}}=\sqrt{\operatorname{rank}(E)}\,\|P-\widetilde{P}\|_{\operatorname{F}},

and this inequality becomes an equality if all singular values of $E$ are the same. Although $\operatorname{rank}(E)\leq n$ always, potentially $\operatorname{rank}(E)=n$ , bringing $n$ into the estimates here and forward. That can be an unfavorable thing to have for huge $n$ .

4.1 Associated SVDs

Let the SVD of $G$ in (4.3d) be

	$G=U\Sigma V^{\operatorname{T}}\equiv\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle m_{1}-r\\$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle U_{1}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle U_{2}$\hfil\kern 5.0pt\crcr}}}}\right]$}}\times\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle m_{2}-r\\\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\Sigma_{1}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\\\scriptstyle m_{1}-r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\Sigma_{2}$\hfil\kern 5.0pt\crcr}}}}\right]$}}\times\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\\\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle V_{1}^{\operatorname{T}}\\\scriptstyle m_{2}-r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle V_{2}^{\operatorname{T}}$\hfil\kern 5.0pt\crcr}}}}\right]$}},$	(4.7a)
where

	$\displaystyle\Sigma_{1}=\operatorname{diag}(\sigma_{1},\ldots,\sigma_{r}),\quad\Sigma_{2}=\begin{bmatrix}\operatorname{diag}(\sigma_{r+1},\ldots,\sigma_{m_{2}})\\ 0_{(m_{1}-m_{2})\times(m_{2}-r)}\end{bmatrix},$	(4.7b)
	$\displaystyle\sigma_{1}\geq\sigma_{2}\geq\cdots\geq\sigma_{m_{2}}.$	(4.7c)

Despite of its large size, $G$ still has only $n$ nonzero singular values, namely $\{\sigma_{i}\}_{i=1}^{n}$ , which are the Hankel singular values of the system, and the rest of its singular values $\sigma_{i}=0$ for $i=n+1,\ldots,m_{2}$ .

Lemma 4.2.

Suppose that (4.1a) holds, and let the singular values of $\widetilde{G}$ be

\widetilde{\sigma}_{1}\geq\widetilde{\sigma}_{2}\geq\cdots\geq\widetilde{\sigma}_{m_{2}}.

Then $\widetilde{\sigma}_{i}\leq\sigma_{i}$ for $i=1,2,\ldots,m_{2}$ . As a corollary, $\|\widetilde{G}\|_{2}=\widetilde{\sigma}_{1}\leq\sigma_{1}$ .

Proof.

The nonzero singular values of $\widetilde{G}$ are given by those of $\widetilde{S}^{\operatorname{T}}\widetilde{R}$ . It suffices to show $\widetilde{\sigma}_{i}^{2}\leq\sigma_{i}^{2}$ for $i=1,2,\ldots,\min\{\widetilde{r}_{1},\widetilde{r}_{2}\}$ . Note $\widetilde{\sigma}_{i}^{2}$ for $i=1,2,\ldots,m_{1}$ are the eigenvalues of

\widetilde{G}\widetilde{G}^{\operatorname{T}}=\widetilde{S}^{\operatorname{T}}\widetilde{R}\widetilde{R}^{\operatorname{T}}\widetilde{S}=\widetilde{S}^{\operatorname{T}}\widetilde{Q}\widetilde{S}\preceq\widetilde{S}^{\operatorname{T}}Q\widetilde{S},

whose nonzero eigenvalues are the same as those of

Q\widetilde{S}\widetilde{S}^{\operatorname{T}}=Q\widetilde{P}=RR^{\operatorname{T}}\widetilde{P},

whose nonzero eigenvalues are the same as those of

R^{\operatorname{T}}\widetilde{P}R\preceq R^{\operatorname{T}}PR,

whose nonzero eigenvalues are $\sigma_{i}^{2}$ for $i=1,2,\ldots,n$ . ∎

Partition

U^{\operatorname{T}}(\widetilde{G}-G)V=-U^{\operatorname{T}}\begin{bmatrix}0&\widetilde{S}^{\operatorname{T}}F\\ E^{\operatorname{T}}\widetilde{R}&E^{\operatorname{T}}F\end{bmatrix}V=\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle m_{2}-r\\\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle E_{11}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle E_{12}\\\scriptstyle m_{1}-r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle E_{21}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle E_{22}$\hfil\kern 5.0pt\crcr}}}}\right]$}}.

By Lemma 4.1, we find

\|E_{ij}\|_{2}\leq\|\widetilde{G}-G\|_{2}\leq\varepsilon\,\,\mbox{for $i,j\in\{1,2\}$},

(4.8)

where $\varepsilon$ is defined in (4.6). Now we will apply [19, Theorem 3.1] to $G,\,\widetilde{G}$ to yield an almost SVD decomposition of $\widetilde{G}$ :

\widetilde{G}=\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle m_{1}-r\\$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\check{U}_{1}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\check{U}_{2}$\hfil\kern 5.0pt\crcr}}}}\right]$}}\times\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle m_{2}-r\\\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\check{\Sigma}_{1}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle 0\\\scriptstyle m_{1}-r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle 0$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\check{\Sigma}_{2}$\hfil\kern 5.0pt\crcr}}}}\right]$}}\times\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\\\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\check{V}_{1}^{\operatorname{T}}\\\scriptstyle m_{2}-r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\check{V}_{2}^{\operatorname{T}}$\hfil\kern 5.0pt\crcr}}}}\right]$}},

(4.9)

where


	$\displaystyle\check{U}$	$\displaystyle\equiv\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle m_{1}-r\\$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\check{U}_{1}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\check{U}_{2}$\hfil\kern 5.0pt\crcr}}}}\right]$}}=[U_{1},U_{2}]\begin{bmatrix}I_{r}&\Gamma^{\operatorname{T}}\\ -\Gamma&I_{m_{1}-r}\end{bmatrix}\begin{bmatrix}(I+\Gamma^{\operatorname{T}}\Gamma)^{-1/2}&0\\ 0&(I+\Gamma\Gamma^{\operatorname{T}})^{-1/2}\end{bmatrix},\\\check{V}$	$\displaystyle\equiv\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle m_{2}-r\\$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\check{V}_{1}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\check{V}_{2}$\hfil\kern 5.0pt\crcr}}}}\right]$}}=[V_{1},V_{2}]\begin{bmatrix}I_{r}&-\Omega^{\operatorname{T}}\\ \Omega&I_{m_{2}-r}\end{bmatrix}\begin{bmatrix}(I+\Omega^{\operatorname{T}}\Omega)^{-1/2}&0\\ 0&(I+\Omega\Omega^{\operatorname{T}})^{-1/2}\end{bmatrix}$		(4.10c)

are two orthogonal matrices, $\Omega\in\mathbb{R}^{(m_{2}-r)\times r}$ and $\Gamma\in\mathbb{R}^{(m_{1}-r)\times r}$ .

Theorem 4.1.

Let $\varepsilon$ be as in (4.6), and let

\delta=\sigma_{r}-\sigma_{r+1},\,\,\underaccent{\bar}{\delta}=\delta-2\varepsilon,\,\,\underaccent{\bar{}}{\sigma}_{r}=\sigma_{r}-\varepsilon.

\underaccent{\bar}{\delta}=\delta-2\varepsilon>0\quad\mbox{and}\quad\frac{\varepsilon^{2}}{\underaccent{\bar}{\delta}^{2}}<\frac{1}{4},

then the following statements hold:

(a)

there exist $\Omega\in\mathbb{R}^{(m_{2}-r)\times r}$ and $\Gamma\in\mathbb{R}^{(m_{1}-r)\times r}$ satisfying

$\max\{\|\Omega\|_{2},\|\Gamma\|_{2}\}\leq\frac{2\varepsilon}{\underaccent{\bar}{\delta}}$ (4.11)

such that $\widetilde{G}$ admits decomposition (4.9) with (4.10);

(b)

the singular values of $\widetilde{G}$ is the multiset union of


$\displaystyle\check{\Sigma}_{1}$	$\displaystyle=\check{U}_{1}^{\operatorname{T}}\widetilde{G}\check{V}_{1}$
	$\displaystyle=(I+\Gamma^{\operatorname{T}}\Gamma)^{1/2}(\Sigma_{1}+E_{11}+E_{12}\Omega)(I+\Omega^{\operatorname{T}}\Omega)^{-1/2}$	(4.12a)
	$\displaystyle=(I+\Gamma^{\operatorname{T}}\Gamma)^{-1/2}(\Sigma_{1}+E_{11}+\Gamma^{\operatorname{T}}E_{21})(I+\Omega^{\operatorname{T}}\Omega)^{1/2},$	(4.12b)

and


$\displaystyle\check{\Sigma}_{2}$	$\displaystyle=\check{U}_{2}^{\operatorname{T}}\widetilde{G}\check{V}_{2}$
	$\displaystyle=(I+\Gamma\Gamma^{\operatorname{T}})^{1/2}(\Sigma_{2}+E_{22}-E_{21}\Omega^{\operatorname{T}})(I+\Omega\Omega^{\operatorname{T}})^{-1/2}$	(4.13a)
	$\displaystyle=(I+\Gamma\Gamma^{\operatorname{T}})^{-1/2}(\Sigma_{2}+E_{22}-\Gamma E_{12})(I+\Omega\Omega^{\operatorname{T}})^{1/2};$	(4.13b)

(c)

we have

\sigma_{\min}(\check{\Sigma}_{1})\geq\sigma_{r}-\varepsilon-\frac{2\varepsilon^{2}}{\underaccent{\bar}{\delta}},\quad\sigma_{\max}(\check{\Sigma}_{2})\leq\sigma_{r+1}+\varepsilon+\frac{2\varepsilon^{2}}{\underaccent{\bar}{\delta}},

(4.14)

where $\sigma_{\min}(\check{\Sigma}_{1})$ and $\sigma_{\max}(\check{\Sigma}_{2})$ are the smallest singular value of $\check{\Sigma}_{1}$ and the largest singular value of $\check{\Sigma}_{2}$ , respectively;

(d)

if also $\varepsilon/\underaccent{\bar}{\delta}<1/2$ , then the top $r$ singular values of $\widetilde{G}$ are exactly the $r$ singular values of $\check{\Sigma}_{1}$ ,

\sigma_{\min}(\check{\Sigma}_{1})\geq\sigma_{r}-\varepsilon=\underaccent{\bar}{\sigma}_{r},

(4.15)

and the dominant left and right singular subspaces are spanned by the columns of


$\displaystyle\check{U}_{1}$	$\displaystyle=(U_{1}-U_{2}\Gamma)(I+\Gamma^{\operatorname{T}}\Gamma)^{-1/2},$	(4.16a)
$\displaystyle\check{V}_{1}$	$\displaystyle=(V_{1}+V_{2}\Omega)(I+\Omega^{\operatorname{T}}\Omega)^{-1/2},$	(4.16b)

respectively. In particular,


$\displaystyle\\|\check{U}_{1}-U_{1}\\|_{2}$	$\displaystyle=\frac{\sqrt{2}\,\\|\Gamma\\|_{2}}{\big{[}\sqrt{1+\\|\Gamma\\|_{2}^{2}}(\sqrt{1+\\|\Gamma\\|_{2}^{2}}+1)\big{]}^{1/2}}\leq\\|\Gamma\\|_{2}\leq\frac{2\varepsilon}{\underaccent{\bar}{\delta}},$	(4.17a)
$\displaystyle\\|\check{V}_{1}-V_{1}\\|_{2}$	$\displaystyle=\frac{\sqrt{2}\,\\|\Omega\\|_{2}}{\big{[}\sqrt{1+\\|\Omega\\|_{2}^{2}}(\sqrt{1+\\|\Omega\\|_{2}^{2}}+1)\big{]}^{1/2}}\leq\\|\Omega\\|_{2}\leq\frac{2\varepsilon}{\underaccent{\bar}{\delta}},$	(4.17b)

and¹¹1 It is tempting to wonder if $\|\check{\Sigma}_{1}-\Sigma_{1}\|_{2}\leq\varepsilon$ , considering the standard perturbation result of singular values [23, p.204], [16, p.21-7]. Unfortunately, $\check{\Sigma}_{1}$ is unlikely diagonal. Another set of two inequalities for the same purpose as (4.18) can be obtained as outlined in Remark 4.2.


$\displaystyle\\|\check{\Sigma}_{1}-\Sigma_{1}\\|_{2}$	$\displaystyle\leq\left(1+\frac{4\sigma_{1}}{\underaccent{\bar}{\delta}}\right)\varepsilon,$	(4.18a)
$\displaystyle\\|\check{\Sigma}_{1}^{-1}-\Sigma_{1}^{-1}\\|_{2}$	$\displaystyle\leq\frac{1}{\sigma_{r}\underaccent{\bar}{\sigma}_{r}}\left(1+\frac{4\sigma_{1}}{\underaccent{\bar}{\delta}}\right){\varepsilon}.$	(4.18b)

Proof.

Recall (4.8). Apply [19, Theorem 3.1] to $G,\,\widetilde{G}$ with $\varepsilon_{ij}=\varepsilon$ , $\delta$ and $\underaccent{\bar}{\delta}$ here to yield all conclusions of the theorem, except (4.17) and (4.18), which we will now prove.

To prove (4.17a), we have

	$\displaystyle\check{U}_{1}-U_{1}$	$\displaystyle=U_{1}\big{[}(I+\Gamma^{\operatorname{T}}\Gamma)^{-1/2}-I\big{]}-U_{2}\Gamma(I+\Gamma^{\operatorname{T}}\Gamma)^{-1/2}$
		$\displaystyle=-[U_{1},U_{2}]\begin{bmatrix}I-(I+\Gamma^{\operatorname{T}}\Gamma)^{-1/2}\\ \Gamma(I+\Gamma^{\operatorname{T}}\Gamma)^{-1/2}\end{bmatrix}.$

Let $\Gamma=Z\Xi W^{\operatorname{T}}$ be the SVD of $\Gamma$ . We find

\begin{bmatrix}I-(I+\Gamma^{\operatorname{T}}\Gamma)^{-1/2}\\ \Gamma(I+\Gamma^{\operatorname{T}}\Gamma)^{-1/2}\end{bmatrix}=\begin{bmatrix}W&\\ &Z\end{bmatrix}\begin{bmatrix}I-(I+\Xi^{\operatorname{T}}\Xi)^{-1/2}\\ \Xi(I+\Xi^{\operatorname{T}}\Xi)^{-1/2}\end{bmatrix}W^{\operatorname{T}},

where for the middle matrix on the right, $I-(I+\Xi^{\operatorname{T}}\Xi)^{-1/2}$ is diagonal and $\Xi(I+\Xi^{\operatorname{T}}\Xi)^{-1/2}$ is leading diagonal. Hence the singular values of the middle matrix are given by: for each singular value $\gamma$ of $\Gamma$ ,

	$\displaystyle\sqrt{\left(1-\frac{1}{\sqrt{1+\gamma^{2}}}\right)^{2}+\left(\frac{\gamma}{\sqrt{1+\gamma^{2}}}\right)^{2}}$	$\displaystyle=\sqrt{2\left(1-\frac{1}{\sqrt{1+\gamma^{2}}}\right)}$
		$\displaystyle=\frac{\sqrt{2}\,\gamma}{\big{[}\sqrt{1+\gamma^{2}}(\sqrt{1+\gamma^{2}}+1)\big{]}^{1/2}}$
		$\displaystyle\leq\gamma\leq\\|\Gamma\\|_{2}.$

Therefore, we get

\|\check{U}_{1}-U_{1}\|_{2}=\left\|\begin{bmatrix}I-(I+\Xi^{\operatorname{T}}\Xi)^{-1/2}\\ \Xi(I+\Xi^{\operatorname{T}}\Xi)^{-1/2}\end{bmatrix}\right\|_{2}=\frac{\sqrt{2}\,\|\Gamma\|_{2}}{\big{[}\sqrt{1+\|\Gamma\|_{2}^{2}}(\sqrt{1+\|\Gamma\|_{2}^{2}}+1)\big{]}^{1/2}}\leq\|\Gamma\|_{2},

yielding (4.17a) in light of (4.11). Similarly, we have (4.17b).

Finally, we prove (4.18). We have

	$\displaystyle\check{\Sigma}_{1}-\Sigma_{1}$	$\displaystyle=\check{U}_{1}^{\operatorname{T}}\widetilde{G}\check{V}_{1}-U_{1}^{\operatorname{T}}\widetilde{G}\check{V}_{1}+U_{1}^{\operatorname{T}}\widetilde{G}\check{V}_{1}-U_{1}^{\operatorname{T}}G\check{V}_{1}+U_{1}^{\operatorname{T}}G\check{V}_{1}-U_{1}^{\operatorname{T}}GV_{1}$
		$\displaystyle=(\check{U}_{1}-U_{1})^{\operatorname{T}}\widetilde{G}\check{V}_{1}+U_{1}^{\operatorname{T}}(\widetilde{G}-G)\check{V}_{1}+U_{1}^{\operatorname{T}}G(\check{V}_{1}-V_{1}).$		(4.19)

In light of (4.6) and (4.17), we get

	$\displaystyle\\|\check{\Sigma}_{1}-\Sigma_{1}\\|_{2}$	$\displaystyle\leq\\|\check{U}_{1}-U_{1}\\|_{2}\\|\widetilde{G}\\|_{2}+\\|\widetilde{G}-G\\|_{2}+\\|G\\|_{2}\\|\check{V}_{1}-V_{1}\\|$
		$\displaystyle\leq\left(1+\frac{4\sigma_{1}}{\underaccent{\bar}{\delta}}\right)\varepsilon,$
and
	$\displaystyle\\|\check{\Sigma}_{1}^{-1}-\Sigma_{1}^{-1}\\|_{2}$	$\displaystyle=\\|\check{\Sigma}_{1}^{-1}(\Sigma_{1}-\check{\Sigma}_{1})\Sigma_{1}^{-1}\\|_{2}$
		$\displaystyle\leq\\|\check{\Sigma}_{1}^{-1}\\|_{2}\\|\Sigma_{1}-\check{\Sigma}_{1}\\|_{2}\\|\Sigma_{1}^{-1}\\|_{2}$
		$\displaystyle\leq\frac{1}{\sigma_{r}\underaccent{\bar}{\sigma}_{r}}\left(1+\frac{4\sigma_{1}}{\underaccent{\bar}{\delta}}\right)\varepsilon,$

completing the proof of (4.18). ∎

Remark 4.2.

Another upper bound on $\|\check{\Sigma}_{1}-\Sigma_{1}\|_{2}$ can be obtained as follows. Alternatively to (4.19), we have

	$\displaystyle\check{\Sigma}_{1}-\Sigma_{1}$	$\displaystyle=(I+\Gamma^{\operatorname{T}}\Gamma)^{1/2}(\Sigma_{1}+E_{11}+E_{12}\Omega)(I+\Omega^{\operatorname{T}}\Omega)^{-1/2}-\Sigma_{1}$
		$\displaystyle=(I+\Gamma^{\operatorname{T}}\Gamma)^{1/2}\Sigma_{1}(I+\Omega^{\operatorname{T}}\Omega)^{-1/2}-\Sigma_{1}$
		$\displaystyle\quad+(I+\Gamma^{\operatorname{T}}\Gamma)^{1/2}(E_{11}+E_{12}\Omega)(I+\Omega^{\operatorname{T}}\Omega)^{-1/2}$
		$\displaystyle=(I+\Gamma^{\operatorname{T}}\Gamma)^{1/2}\Sigma_{1}(I+\Omega^{\operatorname{T}}\Omega)^{-1/2}-\Sigma_{1}(I+\Omega^{\operatorname{T}}\Omega)^{-1/2}+\Sigma_{1}(I+\Omega^{\operatorname{T}}\Omega)^{-1/2}-\Sigma_{1}$
		$\displaystyle\quad+(I+\Gamma^{\operatorname{T}}\Gamma)^{1/2}(E_{11}+E_{12}\Omega)(I+\Omega^{\operatorname{T}}\Omega)^{-1/2}$
		$\displaystyle=\Big{[}(I+\Gamma^{\operatorname{T}}\Gamma)^{1/2}-I\big{]}\Sigma_{1}(I+\Omega^{\operatorname{T}}\Omega)^{-1/2}+\Sigma_{1}\big{[}(I+\Omega^{\operatorname{T}}\Omega)^{-1/2}-I\big{]}$
		$\displaystyle\quad+(I+\Gamma^{\operatorname{T}}\Gamma)^{1/2}(E_{11}+E_{12}\Omega)(I+\Omega^{\operatorname{T}}\Omega)^{-1/2},$

and therefore

$\displaystyle\\|\check{\Sigma}_{1}-\Sigma_{1}\\|_{2}$	$\displaystyle\leq\Big{\\|}(I+\Gamma^{\operatorname{T}}\Gamma)^{1/2}-I\\|_{2}\sigma_{1}+\sigma_{1}\big{\\|}(I+\Omega^{\operatorname{T}}\Omega)^{-1/2}-I\big{\\|}_{2}$
	$\displaystyle\quad+\\|(I+\Gamma^{\operatorname{T}}\Gamma)^{1/2}\\|_{2}(1+\\|\Omega\\|_{2})\varepsilon$
	$\displaystyle\leq\frac{\\|\Gamma\\|_{2}^{2}}{\sqrt{1+\\|\Gamma\\|_{2}^{2}}+1}\sigma_{1}+\sigma_{1}\frac{\\|\Omega\\|_{2}^{2}}{\sqrt{1+\\|\Omega\\|_{2}^{2}}(\sqrt{1+\\|\Omega\\|_{2}^{2}}+1)}$
	$\displaystyle\quad+\sqrt{1+\\|\Gamma\\|_{2}^{2}}(1+\\|\Omega\\|_{2})\varepsilon$
	$\displaystyle\leq\frac{(2\varepsilon/\underaccent{\bar}{\delta})^{2}}{\sqrt{1+(2\varepsilon/\underaccent{\bar}{\delta})^{2}}+1}\left(1+\frac{1}{\sqrt{1+(2\varepsilon/\underaccent{\bar}{\delta})^{2}}}\right)\sigma_{1}+\sqrt{1+\left(\frac{2\varepsilon}{\underaccent{\bar}{\delta}}\right)^{2}}\left(1+\frac{2\varepsilon}{\underaccent{\bar}{\delta}}\right)\varepsilon$
	$\displaystyle\leq\frac{1}{2}\left(1+\frac{1}{\sqrt{2}}\right)\sigma_{1}\left(\frac{2\varepsilon}{\underaccent{\bar}{\delta}}\right)^{2}+2\sqrt{2}\,\varepsilon.$	(4.20)

Comparing (4.18a) with (4.20), we find that both contain a term that depends only on $\varepsilon$ : $\varepsilon$ in the former whereas $2\sqrt{2}\,\varepsilon$ in the latter, and clearly the edge goes to (4.18a) for the term, and that both contain a term proportional to $\sigma_{1}$ , and the edge goes to (4.20) because it is $O(\sigma_{1}\varepsilon)$ v.s. $O(\sigma_{1}\varepsilon^{2})$ . In the same way as how (4.18b) is created, we can create an upper bound on $\|\check{\Sigma}_{1}^{-1}-\Sigma_{1}^{-1}\|_{2}$ , using (4.20), instead. Detail is omitted.

As we commented on [19, Theorem 3.1], (4.15) improves the first inequality in (4.14), but it relies on the latter to first establish the fact that the top $r$ singular values of $\widetilde{G}$ are exactly the $r$ singular values of $\check{\Sigma}_{1}$ .

The decomposition (4.9) we built for $\widetilde{G}$ has an SVD look, but it is not an SVD because $\check{\Sigma}_{i}$ for $i=1,2$ are not diagonal. One idea is to perform an SVD on $\check{\Sigma}_{1}$ and update $\check{U}_{1}$ , $\check{V}_{1}$ accordingly to get $\widetilde{U}_{1}$ and $\widetilde{V}_{1}$ for the dominant left and right singular vector matrices, but it is hard, if not impossible, to relate the resulting $\widetilde{U}_{1}$ and $\widetilde{V}_{1}$ to $U_{1}$ and $V_{1}$ , and in return, difficult to relate $\widetilde{X}_{1}$ and $\widetilde{Y}_{1}$ defined in in (3.3) to $X_{1},\,Y_{1}$ defined in (2.16). This is precisely the reason behind our previous comment at the end of Sections 2 and 3 that $X_{1},\,Y_{1}$ defined in (2.16) and $\widetilde{X}_{1}$ and $\widetilde{Y}_{1}$ in (3.3) are difficult to use. Fortunately these concrete forms for $X_{1},\,Y_{1}$ and $\widetilde{X}_{1}$ and $\widetilde{Y}_{1}$ are not essential as far the transfer functions are concerned because of Theorem 2.1. On the other hand, it is rather easy to relate $\check{U}_{1}$ , $\check{V}_{1}$ , and $\check{\Sigma}_{1}$ there to $U_{1}$ , $V_{1}$ , and $\Sigma_{1}$ , respectively, from the SVD of $G=S^{\operatorname{T}}R$ .

In the rest of this paper, we will assume the following setup without explicit mentioning it:

Setup. Approximate Gramians

\widetilde{P}

and

\widetilde{Q}

satisfy (4.1) such that the conditions of Theorem 4.1, including

\varepsilon/\omega<1/2

, hold. True balanced truncation is carried with

X_{1},\,Y_{1}

in (2.20), while

\widetilde{X}_{1}=[\widetilde{S},0_{n\times(m_{1}-\widetilde{r}_{1})}]\check{U}_{1},\,\,\widetilde{Y}_{1}=[\widetilde{R},0_{n\times(m_{2}-\widetilde{r}_{2})}]\check{V}_{1}\check{\Sigma}_{1}^{-1}

are used for approximate balanced truncation. Accordingly,

\widehat{A}_{11}

\widehat{B}_{1}

, and

\widehat{C}_{1}

in the reduced system (2.6) from the true balanced truncation are defined by (2.5), and

\widetilde{A}_{11}

\widetilde{B}_{1}

, and

\widetilde{C}_{1}

in the reduced system (3.4) from approximate balanced truncation by (3.2).

(4.21)

$X_{1},\,Y_{1}$ in (2.20) and $\widetilde{X}_{1},\,\widetilde{Y}_{1}$ just intriduced, produce different reduced models from the usual reduced models by balanced truncation in the literature, but keep the associated transfer functions intact, nonetheless. In particular, $\widetilde{X}_{1}$ and $\widetilde{Y}_{1}$ are introduced for our analysis only. In practice, they cannot be computed because given $\widetilde{S}$ and $\widetilde{R}$ , knowledge on what $m_{1}$ and $m_{2}$ are is not available, a priori.

4.2 Bounds on differences between reduced systems

In this subsection we will bound the differences of the coefficient matrices and transfer functions between the reduced system (2.6) from the true balanced truncation and (3.4) from an approximate balanced truncation.

First we will establish bounds on $\|X_{1}-\widetilde{X}_{1}\|_{2}$ and $\|Y_{1}-\widetilde{Y}_{1}\|_{2}$ .

Lemma 4.3.

We have


$\displaystyle\\|\widetilde{X}_{1}-X_{1}\\|_{2}$	$\displaystyle\leq\sqrt{\epsilon_{1}}+\sqrt{\\|P\\|_{2}}\,\frac{2\varepsilon}{\underaccent{\bar}{\delta}}=:\epsilon_{x},$	(4.22a)
$\displaystyle\\|\widetilde{Y}_{1}-Y_{1}\\|_{2}$	$\displaystyle\leq\frac{\sqrt{\epsilon_{2}}}{\sigma_{r}}+\frac{\sqrt{\\|Q\\|_{2}}}{\sigma_{r}}\,\left(1+\frac{\underaccent{\bar}{\delta}}{2\underaccent{\bar}{\sigma}_{r}}+\frac{2\sigma_{1}}{\underaccent{\bar}{\sigma}_{r}}\right)\,\frac{2\varepsilon}{\underaccent{\bar}{\delta}}=:\epsilon_{y}.$	(4.22b)

Proof.

Recall (4.2c). We have

	$\displaystyle\widetilde{X}_{1}-X_{1}$	$\displaystyle=[\widetilde{S},0_{n\times(m_{1}-\widetilde{r}_{1})}]\check{U}_{1}-S\check{U}_{1}+S\check{U}_{1}-SU_{1}$
		$\displaystyle=[0_{n\times\widetilde{r}_{1}},-E]\check{U}_{1}+S(\check{U}_{1}-U_{1}),$

and hence, upon using (4.5) and (4.17a), and noticing $\|S\|_{2}=\sqrt{\|P\|_{2}}$ , we arrive at (4.22a). For (4.22b), we have

	$\displaystyle\widetilde{Y}_{1}-Y_{1}$	$\displaystyle=[\widetilde{R},0_{n\times(m_{2}-\widetilde{r}_{2})}]\check{V}_{1}\check{\Sigma}_{1}^{-1}-[\widetilde{R},0_{n\times(m_{2}-\widetilde{r}_{2})}]\check{V}_{1}\Sigma_{1}^{-1}$
		$\displaystyle\quad+[\widetilde{R},0_{n\times(m_{2}-\widetilde{r}_{2})}]\check{V}_{1}\Sigma_{1}^{-1}-[\widetilde{R},0_{n\times(m_{2}-\widetilde{r}_{2})}]V_{1}\Sigma_{1}^{-1}$
		$\displaystyle\quad+[\widetilde{R},0_{n\times(m_{2}-\widetilde{r}_{2})}]V_{1}\Sigma_{1}^{-1}-RV_{1}\Sigma_{1}^{-1}$
		$\displaystyle=[\widetilde{R},0_{n\times(m_{2}-\widetilde{r}_{2})}]\check{V}_{1}\left(\check{\Sigma}_{1}^{-1}-\Sigma_{1}^{-1}\right)$
		$\displaystyle\quad+[\widetilde{R},0_{n\times(m_{2}-\widetilde{r}_{2})}]\left(\check{V}_{1}-V_{1}\right)\Sigma_{1}^{-1}+[0,-F]V_{1}\Sigma_{1}^{-1},$

and, therefore, by Lemma 4.1, and (4.17) and (4.18), we get

	$\displaystyle\big{\\|}\widetilde{Y}_{1}-Y_{1}\big{\\|}_{2}$	$\displaystyle\leq\big{\\|}\widetilde{R}\big{\\|}_{2}\big{\\|}\check{\Sigma}_{1}^{-1}-\Sigma_{1}^{-1}\big{\\|}_{2}+\big{\\|}\widetilde{R}\big{\\|}_{2}\big{\\|}\check{V}_{1}-V_{1}\big{\\|}_{2}\big{\\|}\Sigma_{1}^{-1}\big{\\|}_{2}+\\|F\\|_{2}\big{\\|}\Sigma_{1}^{-1}\big{\\|}_{2}$
		$\displaystyle\leq\sqrt{\\|Q\\|_{2}}\,\frac{1}{\sigma_{r}\underaccent{\bar}{\sigma}_{r}}\left(1+\frac{4\sigma_{1}}{\underaccent{\bar}{\delta}}\right)\varepsilon+\sqrt{\\|Q\\|_{2}}\,\frac{2\varepsilon}{\underaccent{\bar}{\delta}}\frac{1}{\sigma_{r}}+\frac{\sqrt{\epsilon_{2}}}{\sigma_{r}},$

yielding (4.22b). ∎

The differences between the coefficient matrices of the two reduced systems are bounded in Theorem 4.2 below, where the use of any unitarily invariant norm does not require additional care for proofs, and yet may be of independent interest.

Theorem 4.2.

For any unitarily invariant norm $\|\cdot\|_{\operatorname{ui}}$ , we have


$\displaystyle\frac{\big{\\|}\widetilde{A}_{11}-\widehat{A}_{11}\big{\\|}_{\operatorname{ui}}}{\\|A\\|_{\operatorname{ui}}}$	$\displaystyle\leq\sqrt{\\|P\\|_{2}}\epsilon_{y}+\frac{\sqrt{\\|Q\\|_{2}}}{\sigma_{r}}\,\epsilon_{x}=:\epsilon_{a},$	(4.23a)
$\displaystyle\frac{\big{\\|}\widetilde{B}_{1}-\widehat{B}_{1}\big{\\|}_{\operatorname{ui}}}{\\|B\\|_{\operatorname{ui}}}$	$\displaystyle\leq\epsilon_{y}=:\epsilon_{b},$	(4.23b)
$\displaystyle\frac{\big{\\|}\widetilde{C}_{1}-\widehat{C}_{1}\big{\\|}_{\operatorname{ui}}}{\\|C\\|_{\operatorname{ui}}}$	$\displaystyle\leq\epsilon_{x}=:\epsilon_{c},$	(4.23c)

and


$\displaystyle\frac{\big{\\|}\widetilde{B}_{1}\widetilde{B}_{1}^{\operatorname{T}}-\widehat{B}_{1}\widehat{B}_{1}^{\operatorname{T}}\big{\\|}_{\operatorname{ui}}}{\\|BB^{\operatorname{T}}\\|_{\operatorname{ui}}}$	$\displaystyle\leq\sqrt{\\|Q\\|_{2}}\left(\frac{1}{\underaccent{\bar}{\sigma}_{r}}+\frac{1}{\sigma_{r}}\right)\,\epsilon_{y}=:\epsilon_{b2},$	(4.24a)
$\displaystyle\frac{\big{\\|}\widetilde{C}_{1}\widetilde{C}_{1}^{\operatorname{T}}-\widehat{C}_{1}\widehat{C}_{1}^{\operatorname{T}}\big{\\|}_{\operatorname{ui}}}{\\|CC^{\operatorname{T}}\\|_{\operatorname{ui}}}$	$\displaystyle\leq 2\sqrt{\\|P\\|_{2}}\,\epsilon_{x}=:\epsilon_{c2},$	(4.24b)

Proof.

In light of (4.4), it is not difficult to show that

\|X_{1}\|_{2}\leq\sqrt{\|P\|_{2}},\quad\big{\|}\widetilde{X}_{1}\big{\|}_{2}\leq\sqrt{\|P\|_{2}},\quad\|Y_{1}\|_{2}\leq\frac{\sqrt{\|Q\|_{2}}}{\sigma_{r}},\quad\big{\|}\widetilde{Y}_{1}\big{\|}_{2}\leq\frac{\sqrt{\|Q\|_{2}}}{\underaccent{\bar}{\sigma}_{r}},

(4.25)

except for the last one, for which we have

\big{\|}\widetilde{Y}_{1}\big{\|}_{2}\leq\big{\|}\widetilde{R}\big{\|}_{2}\big{\|}\check{\Sigma}_{1}^{-1}\big{\|}_{2}\\\leq\sqrt{\|Q\|_{2}}\,\frac{1}{\underaccent{\bar}{\sigma}_{r}},

which gives the last inequality in (4.25). Next we have

$\displaystyle\widetilde{A}_{11}-\widehat{A}_{11}$	$\displaystyle=\widetilde{Y}_{1}^{\operatorname{T}}A\widetilde{X}_{1}-Y_{1}^{\operatorname{T}}A\widetilde{X}_{1}+Y_{1}^{\operatorname{T}}A\widetilde{X}_{1}-Y_{1}^{\operatorname{T}}AX_{1}$
	$\displaystyle=(\widetilde{Y}_{1}-Y_{1})^{\operatorname{T}}A\widetilde{X}_{1}+Y_{1}^{\operatorname{T}}A(\widetilde{X}_{1}-X_{1}),$	(4.26)
$\displaystyle\widetilde{B}_{1}-\widetilde{B}_{1}$	$\displaystyle=\widetilde{Y}_{1}^{\operatorname{T}}B-Y_{1}^{\operatorname{T}}B=(\widetilde{Y}_{1}-Y_{1})^{\operatorname{T}}B,$
$\displaystyle\widetilde{C}_{1}-\widehat{C}_{1}$	$\displaystyle=\widetilde{X}_{1}^{\operatorname{T}}C-X_{1}^{\operatorname{T}}C=(\widetilde{X}_{1}-X_{1})^{\operatorname{T}}C,$

and

	$\displaystyle\widetilde{B}_{1}\widetilde{B}_{1}^{\operatorname{T}}-\widehat{B}_{1}\widehat{B}_{1}^{\operatorname{T}}$	$\displaystyle=\widetilde{Y}_{1}^{\operatorname{T}}BB^{\operatorname{T}}\widetilde{Y}_{1}-Y_{1}^{\operatorname{T}}BB^{\operatorname{T}}Y_{1}$
		$\displaystyle=\widetilde{Y}_{1}^{\operatorname{T}}BB^{\operatorname{T}}\widetilde{Y}_{1}-\widetilde{Y}_{1}^{\operatorname{T}}BB^{\operatorname{T}}Y_{1}+\widetilde{Y}_{1}^{\operatorname{T}}BB^{\operatorname{T}}Y_{1}-Y_{1}^{\operatorname{T}}BB^{\operatorname{T}}Y_{1}$
		$\displaystyle=\widetilde{Y}_{1}^{\operatorname{T}}BB^{\operatorname{T}}(\widetilde{Y}_{1}-Y_{1})+(\widetilde{Y}_{1}-Y_{1})^{\operatorname{T}}BB^{\operatorname{T}}Y_{1},$
	$\displaystyle\widetilde{C}_{1}\widetilde{C}_{1}^{\operatorname{T}}-\widehat{C}_{1}\widehat{C}_{1}^{\operatorname{T}}$	$\displaystyle=\widetilde{X}_{1}^{\operatorname{T}}CC^{\operatorname{T}}\widetilde{X}_{1}-X_{1}^{\operatorname{T}}CC^{\operatorname{T}}X_{1}$
		$\displaystyle=\widetilde{X}_{1}^{\operatorname{T}}CC^{\operatorname{T}}\widetilde{X}_{1}-\widetilde{X}_{1}^{\operatorname{T}}CC^{\operatorname{T}}X_{1}+\widetilde{X}_{1}^{\operatorname{T}}CC^{\operatorname{T}}X_{1}-X_{1}^{\operatorname{T}}CC^{\operatorname{T}}X_{1}$
		$\displaystyle=\widetilde{X}_{1}^{\operatorname{T}}CC^{\operatorname{T}}(\widetilde{X}_{1}-X_{1})+(\widetilde{X}_{1}-X_{1})^{\operatorname{T}}CC^{\operatorname{T}}X_{1}.$

Take any unitarily invariant norm, e.g., on (4.26), to get

\big{\|}\widetilde{A}_{11}-\widehat{A}_{11}\big{\|}_{\operatorname{ui}}\leq\big{\|}(\widetilde{Y}_{1}-Y_{1})^{\operatorname{T}}\big{\|}_{2}\|A\|_{\operatorname{ui}}\big{\|}\widetilde{X}_{1}\big{\|}_{2}+\big{\|}Y_{1}^{\operatorname{T}}\big{\|}_{2}\|A\|_{\operatorname{ui}}\big{\|}\widetilde{X}_{1}-X_{1}\big{\|}_{2},

and use (4.22) and (4.25) to conclude (4.23). ∎

Remark 4.3.

With the last inequality in (4.25), alternatively to (4.26), we may use

	$\displaystyle\widetilde{A}_{11}-\widehat{A}_{11}$	$\displaystyle=\widetilde{Y}_{1}^{\operatorname{T}}A\widetilde{X}_{1}-\widetilde{Y}_{1}^{\operatorname{T}}AX_{1}+\widetilde{Y}_{1}^{\operatorname{T}}AX_{1}-Y_{1}^{\operatorname{T}}AX_{1}$
		$\displaystyle=\widetilde{Y}_{1}^{\operatorname{T}}A(\widetilde{X}_{1}-X_{1})+(\widetilde{Y}_{1}-Y_{1})^{\operatorname{T}}AX_{1},$

and get

\frac{\big{\|}\widetilde{A}_{11}-\widehat{A}_{11}\big{\|}_{\operatorname{ui}}}{\|A\|_{\operatorname{ui}}}\leq\frac{\sqrt{\|Q\|_{2}}}{\underaccent{\bar}{\sigma}_{r}}\,\epsilon_{x}+\sqrt{\|P\|_{2}}\,\epsilon_{y},

which is slightly worse than (4.23a) because $0<\underaccent{\bar}{\delta}<\underaccent{\bar}{\sigma}_{r}=\sigma_{r}-\varepsilon\leq\sigma_{r}$ .

Previously, we have introduced $H_{\operatorname{bt}}(\cdot)$ in (2.17) and $\widetilde{H}_{\operatorname{bt}}(\cdot)$ in (3.5) for the transfer functions for the reduced systems by the true and approximate balanced truncation, respectively. Let

H_{\operatorname{d}}(s)=H_{\operatorname{bt}}(s)-\widetilde{H}_{\operatorname{bt}}(s),

the difference between the transfer functions, where subscript ‘d’ is used here and in what follows to stand for ‘difference’ between the related things from the true balanced truncation and its approximation.

We are interested in established bounds for $\left\|H_{\operatorname{d}}(\cdot)\right\|_{{\cal H}_{\infty}}$ and $\left\|H_{\operatorname{d}}(\cdot)\right\|_{{\cal H}_{2}}$ . To this end, we introduce

K_{1}=\int_{0}^{\infty}e^{\widehat{A}_{11}t}e^{\widehat{A}_{11}^{\operatorname{T}}t}dt,\quad K_{2}=\int_{0}^{\infty}e^{\widehat{A}_{11}^{\operatorname{T}}t}e^{\widehat{A}_{11}t}dt,

(4.27)

the solutions to $\widehat{A}_{11}K_{1}+K_{1}\widehat{A}_{11}^{\operatorname{T}}+I_{r}=0$ and $\widehat{A}_{11}^{\operatorname{T}}K_{2}+K_{2}\widehat{A}_{11}+I_{r}=0$ , respectively, and let

\eta_{1}=\|A_{11}\|_{2}\|K_{1}\|_{2}\epsilon_{a},\quad\eta_{2}=\|A_{11}\|_{2}\|K_{2}\|_{2}\epsilon_{a}.

(4.28)

Both $K_{1}$ and $K_{2}$ are well-defined because $\widehat{A}_{11}$ is from the exact balanced truncation and hence inherits its stability from the original state matrix $A$ .

$H_{\operatorname{d}}(s)$ is the transfer function of the system

\left\{\begin{array}[]{lll}{\hat{\boldsymbol{x}}}^{\prime}_{r}(t)&=\widehat{A}_{11}\hat{\boldsymbol{x}}_{r}(t)+\widehat{B}_{1}\boldsymbol{u}(t),\quad\mbox{given $\hat{\boldsymbol{x}}_{r}(0)=\widetilde{\boldsymbol{x}}_{r}(0)$},\\ {\widetilde{\boldsymbol{x}}}^{\prime}_{r}(t)&=\widetilde{A}_{11}\widetilde{\boldsymbol{x}}_{r}(t)+\widetilde{B}_{1}\boldsymbol{u}(t),\\ {\boldsymbol{z}(t)}&=\widehat{C}_{1}^{\operatorname{T}}\hat{\boldsymbol{x}}_{r}(t)-\widetilde{C}_{1}^{\operatorname{T}}\widetilde{\boldsymbol{x}}_{r}(t),\end{array}\right.

(4.29)

or in short,

\mathscr{S}_{\operatorname{d}}=\left(\begin{array}[]{cc|c}\widehat{A}_{11}&&\widehat{B}_{1}\\ &\widetilde{A}_{11}&\widetilde{B}_{1}\\ \hline\cr\vphantom{\vrule height=12.80365pt,width=0.9pt,depth=2.84544pt}\widehat{C}_{1}^{\operatorname{T}}&-\widetilde{C}_{1}^{\operatorname{T}}&\end{array}\right)=:\left(\begin{array}[]{c|c}A_{\operatorname{d}}&B_{\operatorname{d}}\\ \hline\cr\vphantom{\vrule height=12.80365pt,width=0.9pt,depth=2.84544pt}C_{\operatorname{d}}^{\operatorname{T}}&\end{array}\right),

Denoted by $P_{\operatorname{d}},\,Q_{\operatorname{d}}\in\mathbb{R}^{2r\times r}$ , the controllability and observability Gramians of (4.29), respectively. They are the solutions to


$\displaystyle A_{\operatorname{d}}P_{\operatorname{d}}+P_{\operatorname{d}}A_{\operatorname{d}}^{\operatorname{T}}+B_{\operatorname{d}}B_{\operatorname{d}}^{\operatorname{T}}$	$\displaystyle=0,$	(4.30a)
$\displaystyle A_{\operatorname{d}}^{\operatorname{T}}Q_{\operatorname{d}}+Q_{\operatorname{d}}A_{\operatorname{d}}+C_{\operatorname{d}}C_{\operatorname{d}}^{\operatorname{T}}$	$\displaystyle=0,$	(4.30b)

respectively. It is well-known that

\left\|H_{\operatorname{d}}(\cdot)\right\|_{{\cal H}_{\infty}}=\sqrt{\lambda_{\max}(P_{\operatorname{d}}Q_{\operatorname{d}})},\quad\left\|H_{\operatorname{d}}(\cdot)\right\|_{{\cal H}_{2}}=\sqrt{\operatorname{tr}(B_{\operatorname{d}}^{\operatorname{T}}Q_{\operatorname{d}}B_{\operatorname{d}})}=\sqrt{\operatorname{tr}(C_{\operatorname{d}}^{\operatorname{T}}P_{\operatorname{d}}C_{\operatorname{d}})}.

(4.31)

The ${\cal H}_{2}$ -norm of continuous system (1.1) is the energy of the associated impulse response in the time domain [2]. Our goals are then turned into estimating the largest eigenvalues of $P_{\operatorname{d}}Q_{\operatorname{d}}$ and the traces.

Lemma 4.4.

If $\eta_{i}<1/2$ for $i=1,2$ , then

P_{\operatorname{d}}=\underbrace{\begin{bmatrix}I_{r}&I_{r}\\ I_{r}&I_{r}\end{bmatrix}}_{=:P_{0}}+\underbrace{\begin{bmatrix}0&\Delta P_{12}\\ (\Delta P_{12})^{\operatorname{T}}&\Delta P_{22}\end{bmatrix}}_{=:\Delta P_{0}},\quad Q_{\operatorname{d}}=\underbrace{\begin{bmatrix}\Sigma_{1}^{2}&-\Sigma_{1}^{2}\\ -\Sigma_{1}^{2}&\Sigma_{1}^{2}\end{bmatrix}}_{=:Q_{0}}+\underbrace{\begin{bmatrix}0&\Delta Q_{12}\\ (\Delta Q_{12})^{\operatorname{T}}&\Delta Q_{22}\end{bmatrix}}_{=:\Delta Q_{0}},

(4.32)

where $\Delta P_{ij},\,\Delta Q_{ij}\in\mathbb{R}^{r\times r}$ and satisfy


$\displaystyle\\|\Delta P_{12}\\|_{2}$	$\displaystyle\leq\frac{\\|K_{1}\\|_{2}}{1-\eta_{1}}\left(\big{\\|}\widehat{B}_{1}\big{\\|}_{2}\\|B\\|_{2}\epsilon_{b}+\\|A\\|_{2}\epsilon_{a}\right)=:\xi_{1},$	(4.33a)
$\displaystyle\\|\Delta P_{22}\\|_{2}$	$\displaystyle\leq\frac{\\|K_{1}\\|_{2}}{1-2\eta_{1}}\left(\\|BB^{\operatorname{T}}\\|_{2}\epsilon_{b2}+2\\|A\\|_{2}\epsilon_{a}\right)=:\xi_{2},$	(4.33b)

and


$\displaystyle\\|\Delta Q_{12}\\|_{2}$	$\displaystyle\leq\frac{\\|K_{2}\\|_{2}}{1-\eta_{2}}\left(\big{\\|}\widehat{C}_{1}\big{\\|}_{2}\\|C\\|_{2}\epsilon_{c}+\\|A\\|_{2}\epsilon_{a}\right)=:\zeta_{1},$	(4.34a)
$\displaystyle\\|\Delta Q_{22}\\|_{2}$	$\displaystyle\leq\frac{\\|K_{2}\\|_{2}}{1-2\eta_{2}}\left(\\|CC^{\operatorname{T}}\\|_{2}\epsilon_{c2}+2\\|A\\|_{2}\epsilon_{a}\right)=:\zeta_{2}.$	(4.34b)

Proof.

Partition both $P_{\operatorname{d}},\,Q_{\operatorname{d}}$ as

P_{\operatorname{d}}=\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle r\\\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle P_{11}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle P_{12}\\\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle P_{12}^{\operatorname{T}}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle P_{22}$\hfil\kern 5.0pt\crcr}}}}\right]$}},\quad Q_{\operatorname{d}}=\hbox{}\vbox{\kern 0.86108pt\hbox{$\kern 0.0pt\kern 2.5pt\kern-5.0pt\left[\kern 0.0pt\kern-2.5pt\kern-5.55557pt\vbox{\kern-0.86108pt\vbox{\vbox{ \halign{\kern\arraycolsep\hfil\@arstrut$\kbcolstyle#$\hfil\kern\arraycolsep& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep&& \kern\arraycolsep\hfil$\@kbrowstyle#$\ifkbalignright\relax\else\hfil\fi\kern\arraycolsep\cr 5.0pt\hfil\@arstrut$\scriptstyle$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle\scriptstyle r\\\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle Q_{11}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle Q_{12}\\\scriptstyle r$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle Q_{12}^{\operatorname{T}}$\hfil\kern 5.0pt&5.0pt\hfil$\scriptstyle Q_{22}$\hfil\kern 5.0pt\crcr}}}}\right]$}}.

We start by investigating $P_{\operatorname{d}}$ first. Blockwise, (4.30a) is equivaent to the following three equations:


$\displaystyle\widehat{A}_{11}P_{11}+P_{11}\widehat{A}_{11}^{\operatorname{T}}+\widehat{B}_{1}\widehat{B}_{1}^{\operatorname{T}}$	$\displaystyle=0,$	(4.35a)
$\displaystyle\widehat{A}_{11}P_{12}+P_{12}\widetilde{A}_{11}^{\operatorname{T}}+\widehat{B}_{1}\widetilde{B}_{1}^{\operatorname{T}}$	$\displaystyle=0,$	(4.35b)
$\displaystyle\widetilde{A}_{11}P_{22}+P_{22}\widetilde{A}_{11}^{\operatorname{T}}+\widetilde{B}_{1}\widetilde{B}_{1}^{\operatorname{T}}$	$\displaystyle=0.$	(4.35c)

It follows from Section 2.3 that $P_{11}=I_{r}$ , and from Lemma B.1 that both $P_{12}$ and $P_{22}$ are near $I_{r}$ , and therefore the form of $P_{\operatorname{d}}$ as in (4.32). Specifically, by (4.23) and Lemma B.1, we have

P_{12}=I_{r}+\Delta P_{12},\quad P_{22}=I_{r}+\Delta P_{22}

with

	$\displaystyle\\|\Delta P_{12}\\|_{2}$	$\displaystyle\leq\frac{\\|K_{1}\\|_{2}}{1-\eta_{1}}\left({\big{\\|}\widehat{B}_{1}\big{\\|}_{2}\\|\widetilde{B}_{1}-\widehat{B}_{1}\\|_{2}}+\\|\widetilde{A}_{11}-\widehat{A}_{11}\\|_{2}\right),$
	$\displaystyle\\|\Delta P_{22}\\|_{2}$	$\displaystyle\leq\frac{\\|K_{1}\\|_{2}}{1-2\eta_{1}}\left(\\|\widetilde{B}_{1}\widetilde{B}_{1}^{\operatorname{T}}-\widehat{B}_{1}\widehat{B}_{1}\\|_{2}+2\\|\widetilde{A}_{11}-\widehat{A}_{11}\\|_{2}\right).$

They, together with Theorem 4.2, yield (4.33).

We now turn our attention to $Q_{\operatorname{d}}$ . Blockwise, (4.30b) is equivalent to the following three equations:


$\displaystyle\widehat{A}_{11}^{\operatorname{T}}Q_{11}+Q_{11}\widehat{A}_{11}+\widehat{C}_{1}\widehat{C}_{1}^{\operatorname{T}}$	$\displaystyle=0,$	(4.36a)
$\displaystyle\widehat{A}_{11}^{\operatorname{T}}Q_{12}+Q_{12}\widetilde{A}_{11}-\widehat{C}_{1}\widetilde{C}_{1}^{\operatorname{T}}$	$\displaystyle=0,$	(4.36b)
$\displaystyle\widetilde{A}_{11}^{\operatorname{T}}Q_{22}+Q_{22}\widetilde{A}_{11}+\widetilde{C}_{1}\widetilde{C}_{1}^{\operatorname{T}}$	$\displaystyle=0.$	(4.36c)

It follows from Section 2.3 that $Q_{11}=\Sigma_{1}^{2}$ , and from Lemma B.1 that both $-Q_{12}$ and $Q_{22}$ are near $\Sigma_{1}^{2}$ , and therefore the form of $Q_{\operatorname{d}}$ as in (4.32). Specifically, by (4.23) and Lemma B.1, we have

Q_{12}=-\Sigma_{1}^{2}+\Delta Q_{12},\quad Q_{22}=\Sigma_{1}^{2}+\Delta Q_{22}

with

	$\displaystyle\\|\Delta Q_{12}\\|_{2}$	$\displaystyle\leq\frac{\\|K_{2}\\|_{2}}{1-\eta_{2}}\left({\big{\\|}\widehat{C}_{1}\big{\\|}_{2}\\|\widetilde{C}_{1}-\widehat{C}_{1}\\|_{2}}+\\|\widetilde{A}_{11}-\widehat{A}_{11}\\|_{2}\right),$
	$\displaystyle\\|\Delta Q_{22}\\|_{2}$	$\displaystyle\leq\frac{\\|K_{2}\\|_{2}}{1-2\eta_{2}}\left(\\|\widetilde{C}_{1}\widetilde{C}_{1}^{\operatorname{T}}-\widehat{C}_{1}\widehat{C}_{1}\\|_{2}+2\\|\widetilde{A}_{11}-\widehat{A}_{11}\\|_{2}\right).$

They, together with Theorem 4.2, yield (4.34). ∎

Remark 4.4.

Bounds on $\|\Delta P_{ij}\|_{\operatorname{F}}$ and $\|\Delta Q_{ij}\|_{\operatorname{F}}$ can also be established, only a little more complicated than (4.33) and (4.34), upon using Lemma B.1 with the Frobenius norm and noticing

\|I_{r}\|_{\operatorname{F}}=\sqrt{r},\quad\|\Sigma_{1}^{2}\|_{\operatorname{F}}=\left(\sum_{i=1}^{r}\sigma_{i}^{4}\right)^{1/2}\leq\sqrt{r}\,\sigma_{1}^{2}.

In fact, we will have

	$\displaystyle\\|\Delta P_{12}\\|_{\operatorname{F}}$	$\displaystyle\leq\frac{\\|K_{1}\\|_{2}}{1-\eta_{1}}\left({\big{\\|}\widehat{B}_{1}\big{\\|}_{2}\\|\widetilde{B}_{1}-\widehat{B}_{1}\\|_{\operatorname{F}}}+\sqrt{r}\,\\|\widetilde{A}_{11}-\widehat{A}_{11}\\|_{2}\right)$
		$\displaystyle\leq\frac{\\|K_{1}\\|_{2}}{1-\eta_{1}}\left(\big{\\|}\widehat{B}_{1}\big{\\|}_{2}\\|B\\|_{\operatorname{F}}\epsilon_{b}+\sqrt{r}\,\\|A\\|_{2}\epsilon_{a}\right),$
	$\displaystyle\\|\Delta P_{22}\\|_{\operatorname{F}}$	$\displaystyle\leq\frac{\\|K_{1}\\|_{2}}{1-2\eta_{1}}\left(\\|\widetilde{B}_{1}\widetilde{B}_{1}^{\operatorname{T}}-\widehat{B}_{1}\widehat{B}_{1}\\|_{\operatorname{F}}+2\sqrt{r}\,\\|\widetilde{A}_{11}-\widehat{A}_{11}\\|_{2}\right)$
		$\displaystyle\leq\frac{\\|K_{1}\\|_{2}}{1-2\eta_{1}}\left(\\|BB^{\operatorname{T}}\\|_{\operatorname{F}}\epsilon_{b2}+2\sqrt{r}\,\\|A\\|_{2}\epsilon_{a}\right),$
	$\displaystyle\\|\Delta Q_{12}\\|_{\operatorname{F}}$	$\displaystyle\leq\frac{\\|K_{2}\\|_{2}}{1-\eta_{2}}\left({\big{\\|}\widehat{C}_{1}\big{\\|}_{2}\\|\widetilde{C}_{1}-\widehat{C}_{1}\\|_{\operatorname{F}}}+\\|\Sigma_{1}^{2}\\|_{\operatorname{F}}\\|\widetilde{A}_{11}-\widehat{A}_{11}\\|_{2}\right)$
		$\displaystyle\leq\frac{\\|K_{2}\\|_{2}}{1-\eta_{2}}\left(\big{\\|}\widehat{C}_{1}\big{\\|}_{2}\\|C\\|_{\operatorname{F}}\epsilon_{c}+\sqrt{r}\,\sigma_{1}^{2}\\|A\\|_{2}\epsilon_{a}\right),$
	$\displaystyle\\|\Delta Q_{22}\\|_{\operatorname{F}}$	$\displaystyle\leq\frac{\\|K_{2}\\|_{2}}{1-2\eta_{2}}\left(\\|\widetilde{C}_{1}\widetilde{C}_{1}^{\operatorname{T}}-\widehat{C}_{1}\widehat{C}_{1}\\|_{\operatorname{F}}+2\\|\Sigma_{1}^{2}\\|_{\operatorname{F}}\\|\widetilde{A}_{11}-\widehat{A}_{11}\\|_{2}\right)$
		$\displaystyle\leq\frac{\\|K_{2}\\|_{2}}{1-2\eta_{2}}\left(\\|CC^{\operatorname{T}}\\|_{\operatorname{F}}\epsilon_{c2}+2\sqrt{r}\,\sigma_{1}^{2}\\|A\\|_{2}\epsilon_{a}\right).$

But these bounds are not materially better than these straightforwardly obtained from (4.33) and (4.34), together with $\|M\|_{\operatorname{F}}\leq\sqrt{r}\|M\|_{2}$ for any $M\in\mathbb{R}^{r\times r}$ .

Theorem 4.3.

If $\eta_{i}<1/2$ for $i=1,2$ , then

$\displaystyle\left\\|H_{\operatorname{d}}(\cdot)\right\\|_{{\cal H}_{\infty}}$	$\displaystyle\leq\sqrt{2\sigma_{1}^{2}(\xi_{1}+\xi_{2})+2(\zeta_{1}+\zeta_{2})+(\xi_{1}+\xi_{2})\zeta_{1}+\zeta_{2})}=:\epsilon_{{\operatorname{d}},\infty},$	(4.37)
$\displaystyle\left\\|H_{\operatorname{d}}(\cdot)\right\\|_{{\cal H}_{2}}$	$\displaystyle\leq\sqrt{\min\{r,m\}}\left[\sigma_{1}^{2}\big{(}\big{\\|}\widehat{B}_{1}\big{\\|}_{2}+\big{\\|}\widetilde{B}_{1}\big{\\|}_{2}\big{)}\\|B\\|_{2}\epsilon_{b}\right.$
	$\displaystyle\hskip 56.9055pt\left.+\left(2\\|\widehat{B}_{1}^{\operatorname{T}}\\|_{2}\big{\\|}\widetilde{B}_{1}\big{\\|}_{2}\,\zeta_{1}+\\|\widetilde{B}_{1}^{\operatorname{T}}\\|_{2}^{2}\,\zeta_{2}\right)\right]^{1/2}=:\epsilon_{{\operatorname{d}},2},$	(4.38)

where $\xi_{i}$ and $\zeta_{i}$ for $i=1,2$ are defined in Lemma 4.4, and $m$ and $p$ is the numbers of columns of $B$ and $C$ , respectively.

Proof.

Recall (4.32). Noticing that

	$\displaystyle\\|P_{0}\\|_{2}=2,\quad\\|Q_{0}\\|_{2}=2\sigma_{1}^{2},\quad P_{0}Q_{0}=0,$
	$\displaystyle\\|\Delta P_{0}\\|_{2}\leq\\|\Delta P_{12}\\|_{2}+\\|\Delta P_{22}\\|_{2},\quad\\|\Delta Q_{0}\\|_{2}\leq\\|\Delta Q_{12}\\|_{2}+\\|\Delta Q_{22}\\|_{2},$

we get

	$\displaystyle\lambda_{\max}(P_{\operatorname{d}}Q_{\operatorname{d}})$	$\displaystyle\leq\\|P_{\operatorname{d}}Q_{\operatorname{d}}\\|_{2}$
		$\displaystyle\leq\\|P_{0}Q_{0}+P_{0}\Delta Q_{0}+(\Delta P_{0})Q_{0}+(\Delta P_{0})(\Delta Q_{0})\\|_{2}$
		$\displaystyle\leq 2\left(\\|\Delta Q_{12}\\|_{2}+\\|\Delta Q_{22}\\|_{2}\right)+2\sigma_{1}^{2}\left(\\|\Delta P_{12}\\|_{2}+\\|\Delta P_{22}\\|_{2}\right)$
		$\displaystyle\quad+\left(\\|\Delta Q_{12}\\|_{2}+\\|\Delta Q_{22}\\|_{2}\right)\left(\\|\Delta P_{12}\\|_{2}+\\|\Delta P_{22}\\|_{2}\right),$

which together with (4.33) and (4.34) lead to (4.37), upon noticing (4.31).

Next we prove (4.38). We claim that

	$\displaystyle\operatorname{tr}\big{(}B_{\operatorname{d}}^{\operatorname{T}}Q_{0}B_{\operatorname{d}}\big{)}$	$\displaystyle\leq\min\{r,m\}\,\sigma_{1}^{2}\big{(}\big{\\|}\widehat{B}_{1}\big{\\|}_{2}+\big{\\|}\widetilde{B}_{1}\big{\\|}_{2}\big{)}\\|B\\|_{2}\epsilon_{b},$		(4.39)
	$\displaystyle\|\operatorname{tr}(B_{\operatorname{d}}^{\operatorname{T}}[\Delta Q_{0}]B_{\operatorname{d}})\|$	$\displaystyle\leq\min\{r,m\}\left(2\\|\widehat{B}_{1}^{\operatorname{T}}\\|_{2}\big{\\|}\widetilde{B}_{1}\big{\\|}_{2}\,\zeta_{1}+\\|\widetilde{B}_{1}^{\operatorname{T}}\\|_{2}^{2}\,\zeta_{2}\right).$		(4.40)

Note that, for any square matrix $M$ ,

\operatorname{tr}(M)\leq\operatorname{rank}(M)\,\|M\|_{2}.

Using Theorem 4.2, we have

	$\displaystyle\operatorname{tr}\big{(}B_{\operatorname{d}}^{\operatorname{T}}Q_{0}B_{\operatorname{d}}\big{)}$	$\displaystyle=\operatorname{tr}\big{(}\widehat{B}_{1}^{\operatorname{T}}\Sigma_{1}^{2}\widehat{B}_{1}-2\widehat{B}_{1}^{\operatorname{T}}\Sigma_{1}^{2}\widetilde{B}_{1}+\widetilde{B}_{1}^{\operatorname{T}}\Sigma_{1}^{2}\widetilde{B}_{1}\big{)}$
		$\displaystyle=\operatorname{tr}\big{(}\widehat{B}_{1}^{\operatorname{T}}\Sigma_{1}^{2}[\widehat{B}_{1}-\widetilde{B}_{1}]\big{)}+\operatorname{tr}([\widetilde{B}_{1}-\widehat{B}_{1}]^{\operatorname{T}}\Sigma_{1}^{2}\widetilde{B}_{1}\big{)}$
		$\displaystyle\leq\min\{r,m\}\big{\\|}\widehat{B}_{1}^{\operatorname{T}}\Sigma_{1}^{2}[\widehat{B}_{1}-\widetilde{B}_{1}]\big{\\|}_{2}+\big{\\|}[\widetilde{B}_{1}-\widehat{B}_{1}]^{\operatorname{T}}\Sigma_{1}^{2}\widetilde{B}_{1}\big{\\|}_{2}$
		$\displaystyle\leq\min\{r,m\}\big{\\|}\Sigma_{1}^{2}\big{\\|}_{2}\big{(}\big{\\|}\widehat{B}_{1}\big{\\|}_{2}+\big{\\|}\widetilde{B}_{1}\big{\\|}_{2}\big{)}\\|\widehat{B}_{1}-\widetilde{B}_{1}\\|_{2}$
		$\displaystyle\leq\min\{r,m\}\,\sigma_{1}^{2}\big{(}\big{\\|}\widehat{B}_{1}\big{\\|}_{2}+\big{\\|}\widetilde{B}_{1}\big{\\|}_{2}\big{)}\\|B\\|_{2}\epsilon_{b},$

proving (4.39), and

	$\displaystyle\|\operatorname{tr}\big{(}B_{\operatorname{d}}^{\operatorname{T}}[\Delta Q_{0}]B_{\operatorname{d}}\big{)}\|$	$\displaystyle=\|\operatorname{tr}\big{(}2\widehat{B}_{1}^{\operatorname{T}}[\Delta Q_{12}]\widetilde{B}_{1}\big{)}+\widetilde{B}_{1}^{\operatorname{T}}[\Delta Q_{22}]\widetilde{B}_{1}\big{)}\|$
		$\displaystyle\leq\min\{r,m\}\left(2\big{\\|}\widehat{B}_{1}^{\operatorname{T}}[\Delta Q_{12}]\widetilde{B}_{1}\big{\\|}_{2}+\big{\\|}\widetilde{B}_{1}^{\operatorname{T}}[\Delta Q_{22}]\widetilde{B}_{1}\big{\\|}_{2}\right)$
		$\displaystyle\leq\min\{r,m\}\left(2\big{\\|}\widehat{B}_{1}^{\operatorname{T}}\big{\\|}_{2}\big{\\|}\widetilde{B}_{1}\big{\\|}_{2}\\|\Delta Q_{12}\\|_{2}+\big{\\|}\widetilde{B}_{1}^{\operatorname{T}}\big{\\|}_{2}^{2}\\|\Delta Q_{22}\\|_{2}\right),$

yielding (4.40). With (4.39) and (4.40), we are ready to show (4.38). We have

	$\displaystyle\\|H_{\operatorname{d}}(\cdot)\\|_{{\cal H}_{2}}^{2}$	$\displaystyle=\operatorname{tr}\big{(}B_{\operatorname{d}}^{\operatorname{T}}Q_{\operatorname{d}}B_{\operatorname{d}}\big{)}=\operatorname{tr}\big{(}B_{\operatorname{d}}^{\operatorname{T}}Q_{0}B_{\operatorname{d}}\big{)}+\operatorname{tr}\big{(}B_{\operatorname{d}}^{\operatorname{T}}[\Delta Q_{0}]B_{\operatorname{d}}\big{)}$
		$\displaystyle\leq\min\{r,m\}\,\sigma_{1}^{2}\big{(}\big{\\|}\widehat{B}_{1}\big{\\|}_{2}+\big{\\|}\widetilde{B}_{1}\big{\\|}_{2}\big{)}\\|B\\|_{2}\epsilon_{b}$
		$\displaystyle\quad+\min\{r,m\}\left(2\big{\\|}\widehat{B}_{1}^{\operatorname{T}}\big{\\|}_{2}\big{\\|}\widetilde{B}_{1}\big{\\|}_{2}\,\zeta_{1}+\big{\\|}\widetilde{B}_{1}^{\operatorname{T}}\big{\\|}_{2}^{2}\,\zeta_{2}\right),$

as expected. ∎

Remark 4.5.

Alternatively, basing on the second expression in (4.31) for $\|H_{\operatorname{d}}(\cdot)\|_{{\cal H}_{2}}$ , we can derive a different bound. Similarly to (4.39) and (4.40), we claim that

	$\displaystyle\operatorname{tr}\big{(}C_{\operatorname{d}}^{\operatorname{T}}P_{0}C_{\operatorname{d}}\big{)}$	$\displaystyle\leq\min\{r,p\}\big{(}\big{\\|}\widehat{C}_{1}\big{\\|}_{2}+\big{\\|}\widetilde{C}_{1}\big{\\|}_{2}\big{)}\\|C\\|_{2}\epsilon_{c},$		(4.41)
	$\displaystyle\|\operatorname{tr}\big{(}C_{\operatorname{d}}^{\operatorname{T}}[\Delta P_{0}]C_{\operatorname{d}}\big{)}\|$	$\displaystyle\leq\min\{r,p\}\left(2\big{\\|}\widehat{C}_{1}^{\operatorname{T}}\big{\\|}_{2}\big{\\|}\widetilde{C}_{1}\big{\\|}_{2}\,\xi_{1}+\big{\\|}\widetilde{C}_{1}^{\operatorname{T}}\big{\\|}_{2}^{2}\,\xi_{2}\right).$		(4.42)

They can be proven, analogously along the line we proved (4.39) and (4.40), as follows:

	$\displaystyle\operatorname{tr}\big{(}C_{\operatorname{d}}^{\operatorname{T}}P_{0}C_{\operatorname{d}}\big{)}$	$\displaystyle\leq\min\{r,p\}\big{(}\big{\\|}\widehat{C}_{1}\big{\\|}_{2}+\big{\\|}\widetilde{C}_{1}\big{\\|}_{2}\big{)}\\|\widehat{C}_{1}-\widetilde{C}_{1}\\|_{2}$
		$\displaystyle\leq\min\{r,p\}\big{(}\big{\\|}\widehat{C}_{1}\big{\\|}_{2}+\big{\\|}\widetilde{C}_{1}\big{\\|}_{2}\big{)}\\|C\\|_{2}\epsilon_{c},$
	$\displaystyle\|\operatorname{tr}\big{(}C_{\operatorname{d}}^{\operatorname{T}}[\Delta P_{0}]C_{\operatorname{d}}\big{)}\|$	$\displaystyle=\|\operatorname{tr}\big{(}2\widehat{C}_{1}^{\operatorname{T}}[\Delta P_{12}]\widetilde{C}_{1}\big{)}+\widetilde{C}_{1}^{\operatorname{T}}[\Delta P_{22}]\widetilde{C}_{1}\big{)}\|$
		$\displaystyle\leq\min\{r,m\}\left(2\big{\\|}\widehat{C}_{1}^{\operatorname{T}}[\Delta P_{12}]\widetilde{C}_{1}\big{\\|}_{2}+\big{\\|}\widetilde{C}_{1}^{\operatorname{T}}[\Delta P_{22}]\widetilde{C}_{1}\big{\\|}_{2}\right)$
		$\displaystyle\leq\min\{r,m\}\left(2\big{\\|}\widehat{C}_{1}^{\operatorname{T}}\big{\\|}_{2}\big{\\|}\widetilde{C}_{1}\big{\\|}_{2}\\|\Delta P_{12}\\|_{2}+\big{\\|}\widetilde{C}_{1}^{\operatorname{T}}\big{\\|}_{2}^{2}\\|\Delta P_{22}\\|_{2}\right).$

Finally,

	$\displaystyle\\|H_{\operatorname{d}}(\cdot)\\|_{{\cal H}_{2}}^{2}$	$\displaystyle=\operatorname{tr}\big{(}C_{\operatorname{d}}^{\operatorname{T}}P_{\operatorname{d}}C_{\operatorname{d}}\big{)}=\operatorname{tr}\big{(}C_{\operatorname{d}}^{\operatorname{T}}P_{0}C_{\operatorname{d}}\big{)}+\operatorname{tr}\big{(}C_{\operatorname{d}}^{\operatorname{T}}[\Delta P_{0}]C_{\operatorname{d}}\big{)}$
		$\displaystyle\leq\min\{r,p\}\big{(}\big{\\|}\widehat{C}_{1}\big{\\|}_{2}+\big{\\|}\widetilde{C}_{1}\big{\\|}_{2}\big{)}\\|C\\|_{2}\epsilon_{c}$
		$\displaystyle\quad+\min\{r,p\}\left(2\big{\\|}\widehat{C}_{1}^{\operatorname{T}}\big{\\|}_{2}\big{\\|}\widetilde{C}_{1}\big{\\|}_{2}\,\xi_{1}+\big{\\|}\widetilde{C}_{1}^{\operatorname{T}}\big{\\|}_{2}^{2}\,\xi_{2}\right),$

yielding a different $\epsilon_{{\operatorname{d}},2}$ from the one in (4.38). It is not clear which one is smaller.

Norms of the coefficient matrices for the reduced systems appear in the bounds in Theorem 4.3. They can be replaced by the norms of the corresponding coefficient matrices for the original system with the help of the next lemma.

Lemma 4.5.

For any unitarily invariant norm $\|\cdot\|_{\operatorname{ui}}$ , we have


$\displaystyle\frac{\\|\widehat{A}_{11}\\|_{\operatorname{ui}}}{\\|A\\|_{\operatorname{ui}}}$	$\displaystyle\leq\frac{\sqrt{\\|P\\|_{2}\\|Q\\|_{2}}}{\sigma_{r}},\quad$	$\displaystyle\frac{\\|\widetilde{A}_{11}\\|_{\operatorname{ui}}}{\\|A\\|_{\operatorname{ui}}}$	$\displaystyle\leq\frac{\sqrt{\\|P\\|_{2}\\|Q\\|_{2}}}{\sigma_{r}}+\epsilon_{a},$	(4.43a)
$\displaystyle\frac{\big{\\|}\widehat{B}_{1}\big{\\|}_{\operatorname{ui}}}{\\|B\\|_{\operatorname{ui}}}$	$\displaystyle\leq\frac{\sqrt{\\|Q\\|_{2}}}{\sigma_{r}},\quad$	$\displaystyle\frac{\big{\\|}\widetilde{B}_{1}\big{\\|}_{\operatorname{ui}}}{\\|B\\|_{\operatorname{ui}}}$	$\displaystyle\leq\frac{\sqrt{\\|Q\\|_{2}}}{\sigma_{r}}+\epsilon_{b},$	(4.43b)
$\displaystyle\frac{\big{\\|}\widehat{C}_{1}\big{\\|}_{\operatorname{ui}}}{\\|C\\|_{\operatorname{ui}}}$	$\displaystyle\leq\sqrt{\\|P\\|_{2}},\quad$	$\displaystyle\frac{\big{\\|}\widehat{C}_{1}\big{\\|}_{\operatorname{ui}}}{\\|C\\|_{\operatorname{ui}}}$	$\displaystyle\leq\sqrt{\\|P\\|_{2}}+\epsilon_{c},$	(4.43c)

where $\epsilon_{a}$ , $\epsilon_{b}$ , and $\epsilon_{c}$ are as in (4.23).

Proof.

We have by (2.20)

	$\displaystyle\big{\\|}\widehat{A}_{11}\big{\\|}_{\operatorname{ui}}=\big{\\|}Y_{1}^{\operatorname{T}}AX_{1}\big{\\|}_{\operatorname{ui}}$	$\displaystyle\leq\big{\\|}Y_{1}^{\operatorname{T}}\big{\\|}_{2}\\|A\\|_{\operatorname{ui}}\\|X_{1}\\|_{2}$
		$\displaystyle\leq\\|R\\|_{2}\big{\\|}\Sigma_{1}^{-1}\big{\\|}_{2}\\|A\\|_{\operatorname{ui}}\\|S\\|_{2}$
		$\displaystyle=\frac{\sqrt{\\|P\\|_{2}\\|Q\\|_{2}}}{\sigma_{r}}\,\\|A\\|_{\operatorname{ui}},$
	$\displaystyle\big{\\|}\widehat{B}_{1}\big{\\|}_{\operatorname{ui}}=\\|Y_{1}^{\operatorname{T}}B\\|_{\operatorname{ui}}$	$\displaystyle\leq\frac{\sqrt{\\|Q\\|_{2}}}{\sigma_{r}}\,\\|B\\|_{\operatorname{ui}},$
	$\displaystyle\big{\\|}\widehat{C}_{1}\big{\\|}_{\operatorname{ui}}=\\|Y_{1}^{\operatorname{T}}C\\|_{\operatorname{ui}}$	$\displaystyle\leq\sqrt{\\|P\\|_{2}}\,\\|C\\|_{\operatorname{ui}},$

Therefore

\big{\|}\widetilde{A}_{11}\big{\|}_{\operatorname{ui}}\leq\big{\|}\widehat{A}_{11}\big{\|}_{\operatorname{ui}}+\big{\|}\widetilde{A}_{11}-\widehat{A}_{11}\big{\|}_{\operatorname{ui}}\leq\left(\frac{\sqrt{\|P\|_{2}\|Q\|_{2}}}{\sigma_{r}}+\epsilon_{a}\right)\|A\|_{\operatorname{ui}},

and similarly for $\big{\|}\widetilde{B}_{1}\big{\|}_{\operatorname{ui}}$ and $\big{\|}\widetilde{C}_{1}\big{\|}_{\operatorname{ui}}$ . The proofs of the other two inequalities are similar. ∎

4.3 Transfer function for approximate balanced truncation

In this subsection, we establish bounds to measure the quality of the reduced system (3.4) from approximate balanced truncation as an approximation to the original system (1.1). Even though the projection matrices $\widetilde{X}_{1},\,\widetilde{Y}_{1}$ we used for approximate balanced truncation are different from the ones in practice, the transfer function as a result remains the same, nonetheless. Therefore our bounds are applicable in real applications. These bounds are the immediate consequences of Theorem 4.3 upon using

\big{\|}H(\cdot)-\widetilde{H}_{\operatorname{bt}}(\cdot)\big{\|}\leq\big{\|}H(\cdot)-H_{\operatorname{bt}}(\cdot)\big{\|}+\big{\|}H_{\operatorname{bt}}(\cdot)-\widetilde{H}_{\operatorname{bt}}(\cdot)\big{\|}

for $\|\cdot\|=\|\cdot\|_{{\cal H}_{\infty}}$ and $\|\cdot\|_{{\cal H}_{2}}$ .

Theorem 4.4.

Under the conditions of Theorem 4.3, we have

	$\displaystyle\big{\\|}H(\cdot)-\widetilde{H}_{\operatorname{bt}}(\cdot)\big{\\|}_{{\cal H}_{\infty}}$	$\displaystyle\leq 2\sum_{k=r+1}^{n}\sigma_{k}+\epsilon_{{\operatorname{d}},\infty},$		(4.44)
	$\displaystyle\big{\\|}H(\cdot)-\widetilde{H}_{\operatorname{bt}}(\cdot)\big{\\|}_{{\cal H}_{2}}$	$\displaystyle\leq\big{\\|}H(\cdot)-H_{\operatorname{bt}}(\cdot)\big{\\|}_{{\cal H}_{2}}+\epsilon_{{\operatorname{d}},2},$		(4.45)

where $\epsilon_{{\operatorname{d}},\infty}$ and $\epsilon_{{\operatorname{d}},2}$ are as in Theorem 4.3.

An immediate explanation to both inequalities (4.44) and (4.45) is that the reduced system (3.4) from the approximate balanced reduction as an approximation to the original system (1.1) is worse than the one from the true balanced reduction by no more than $\epsilon_{{\operatorname{d}},\infty}$ and $\epsilon_{{\operatorname{d}},2}$ in terms of the ${\cal H}_{\infty}$ - and ${\cal H}_{2}$ -norm, respectively. Both $\epsilon_{{\operatorname{d}},\infty}$ and $\epsilon_{{\operatorname{d}},2}$ can be traced back to the initial approximation errors $\epsilon_{1}$ and $\epsilon_{2}$ in the computed Gramians as specified in (4.1) albeit complicatedly. To better understand what $\epsilon_{{\operatorname{d}},\infty}$ and $\epsilon_{{\operatorname{d}},2}$ are in terms of $\epsilon_{1}$ and $\epsilon_{2}$ , we summarize all quantities that lead to them, up to the first order in

\epsilon_{\operatorname{}}:=\max\{\epsilon_{1},\epsilon_{2}\}.

Then $\varepsilon\leq\rho\sqrt{\epsilon_{\operatorname{}}}+\epsilon_{\operatorname{}}$ in (4.6). Let $\rho=\max\big{\{}\sqrt{\|P\|_{2}},\sqrt{\|Q\|_{2}}\big{\}}$ . We have

$\displaystyle\epsilon_{x}$	$\displaystyle\leq\left(1+\frac{2\rho^{2}}{\delta}\right)\sqrt{\epsilon_{\operatorname{}}}+O(\epsilon_{\operatorname{}}),$	$\displaystyle\quad(\mbox{see \eqref{eq:diff(X1Y1):BTaBT}})$
$\displaystyle\epsilon_{y}$	$\displaystyle\leq\frac{1}{\sigma_{r}}\left[1+\left(1+\frac{\delta}{2\sigma_{r}}+\frac{2\sigma_{1}}{\sigma_{r}}\right)\,\frac{2\rho^{2}}{\delta}\right]\sqrt{\epsilon_{\operatorname{}}}+O(\epsilon_{\operatorname{}}),$	$\displaystyle\quad(\mbox{see \eqref{eq:diff(X1Y1):BTaBT}})$
$\displaystyle\epsilon_{a}$	$\displaystyle\leq\frac{\rho}{\sigma_{r}}\epsilon_{x}+\rho\epsilon_{y},$	$\displaystyle\quad(\mbox{see \eqref{eq:diff(ABC):BTaBT}})$
$\displaystyle\epsilon_{b}$	$\displaystyle=\epsilon_{y},$	$\displaystyle\quad(\mbox{see \eqref{eq:diff(ABC):BTaBT}})$
$\displaystyle\epsilon_{c}$	$\displaystyle=\epsilon_{x},$	$\displaystyle\quad(\mbox{see \eqref{eq:diff(ABC):BTaBT}})$
$\displaystyle\epsilon_{b2}$	$\displaystyle=\frac{2\rho}{\sigma_{r}}\,\epsilon_{y}+O(\epsilon_{\operatorname{}}),$	$\displaystyle\quad(\mbox{see \eqref{eq:diff(BBCC):BTaBT}})$
$\displaystyle\epsilon_{c2}$	$\displaystyle=2\rho\epsilon_{x},$	$\displaystyle\quad(\mbox{see \eqref{eq:diff(BBCC):BTaBT}})$
$\displaystyle\xi_{1}$	$\displaystyle\leq\\|K_{1}\\|_{2}\left(\\|A\\|_{2}\epsilon_{a}+\frac{\rho}{\sigma_{r}}\\|B\\|_{2}^{2}\epsilon_{b}\right)+O(\epsilon_{\operatorname{}}),$	$\displaystyle\quad(\mbox{see \eqref{eq:bd4Pd}})$
$\displaystyle\xi_{2}$	$\displaystyle\leq\\|K_{1}\\|_{2}\left(2\\|A\\|_{2}\epsilon_{a}+\\|B\\|_{2}^{2}\epsilon_{b2}\right)+O(\epsilon_{\operatorname{}}),$	$\displaystyle\quad(\mbox{see \eqref{eq:bd4Pd}})$
$\displaystyle\zeta_{1}$	$\displaystyle\leq\\|K_{2}\\|_{2}\left(\\|A\\|_{2}\epsilon_{a}+{\rho}\\|C\\|_{2}^{2}\epsilon_{c}\right)+O(\epsilon_{\operatorname{}}),$	$\displaystyle\quad(\mbox{see \eqref{eq:bd4Qd}})$
$\displaystyle\zeta_{2}$	$\displaystyle\leq\\|K_{2}\\|_{2}\left(2\\|A\\|_{2}\epsilon_{a}+\\|C\\|_{2}^{2}\epsilon_{c2}\right)+O(\epsilon_{\operatorname{}}),$	$\displaystyle\quad(\mbox{see \eqref{eq:bd4Qd}})$
$\displaystyle\epsilon_{{\operatorname{d}},\infty}$	$\displaystyle=\sqrt{2\sigma_{1}^{2}(\xi_{1}+\xi_{2})+2(\zeta_{1}+\zeta_{2})}+O(\sqrt{\epsilon_{\operatorname{}}}),$	$\displaystyle\quad(\mbox{see \eqref{eq:bd4Hd}})$
$\displaystyle\epsilon_{{\operatorname{d}},2}$	$\displaystyle\leq\sqrt{\min\{r,m\}}\\|B\\|_{2}\left[2\sigma_{1}^{2}\frac{\rho}{\sigma_{r}}\epsilon_{b}+\left(\frac{\rho}{\sigma_{r}}\right)^{2}\big{(}2\zeta_{1}+\zeta_{2}\big{)}\right]^{1/2}+O(\sqrt{\epsilon_{\operatorname{}}}).$	$\displaystyle\quad(\mbox{see \eqref{eq:bd4Hd:nrmH2}})$

Alternatively, for $\epsilon_{{\operatorname{d}},2}$ , also by Remark 4.5

\epsilon_{{\operatorname{d}},2}\leq\sqrt{\min\{r,p\}}\|C\|_{2}\left[2\rho\epsilon_{b}+\rho^{2}\big{(}2\xi_{1}+\xi_{2}\big{)}\right]^{1/2}+O(\sqrt{\epsilon_{\operatorname{}}}).

It can be seen that both $\epsilon_{{\operatorname{d}},\infty}$ and $\epsilon_{{\operatorname{d}},2}$ are of $O(\epsilon_{\operatorname{}}^{1/4})$ , pretty disappointing.

5 Concluding Remarks

For a continuous linear time-invariant dynamic system, the existing global error bound that bounds the error between a reduced model via balanced truncation and the original dynamic system assumes that the reduced model is constructed from two exact controllability and observability Gramians. But in practice, the Gramians are usually approximated by some computed low-rank approximations, especially when the original dynamic system is large scale. Thus, rigorously speaking, the existing global error bound, although indicative about the accuracy in the reduced system, is not really applicable. In this paper, we perform an error analysis, assuming the reduced model is constructed from two low-rank approximations of the Gramians, making up the deficiency in the current theory for measuring the quality of the reduced model obtained by approximate balanced truncation. Error bounds have been obtained for the purpose.

So far, we have been focused on continuous linear time-invariant dynamic systems, but our techniques should be extendable to discrete time-invariant dynamic systems without much difficulty.

Throughout this paper, our presentation is restricted to the real number field $\mathbb{R}$ . This restriction is more for simplicity and clarity than the capability of our techniques. In fact, our approach can be straightforwardly modified to cover the complex number case: replace all transposes $(\cdot)^{\operatorname{T}}$ of vectors/matrices with complex conjugate transposes $(\cdot)^{\operatorname{H}}$ .

Appendix

Appendix A Some results on subspaces

Consider two subspaces ${\cal U}$ and $\widetilde{\cal U}$ with dimension $r$ of ${\mathbb{R}}^{n}$ and let $U\in{\mathbb{R}}^{n\times r}$ and $\widetilde{U}\in{\mathbb{R}}^{n\times r}$ be orthonormal basis matrices of ${\cal U}$ and ${\cal U}$ , respectively, i.e.,

U^{\operatorname{T}}U=I_{r},\,{\cal U}={\cal R}(U),\quad\mbox{and}\quad\widetilde{U}^{\operatorname{T}}\widetilde{U}=I_{r},\,\widetilde{\cal U}={\cal R}(\widetilde{U}),

and denote by $\tau_{j}$ for $1\leq j\leq r$ in the descending order, i.e., $\tau_{1}\geq\cdots\geq\tau_{r}$ , the singular values of $\widetilde{U}^{\operatorname{T}}U$ . The $r$ canonical angles $\theta_{j}({\cal U},{\cal\widetilde{U}})$ between ${\cal U}$ to ${\cal\widetilde{U}}$ are defined by

0\leq\theta_{j}({\cal U},{\cal\widetilde{U}}):=\arccos\tau_{j}\leq\frac{\pi}{2}\quad\mbox{for $1\leq j\leq r$}.

They are in the ascending order, i.e., $\theta_{1}({\cal U},{\cal\widetilde{U}})\leq\cdots\leq\theta_{r}({\cal U},{\cal\widetilde{U}})$ . Set

\Theta({\cal U},{\cal\widetilde{U}})=\operatorname{diag}(\theta_{1}({\cal U},{\cal\widetilde{U}}),\ldots,\theta_{r}({\cal U},{\cal\widetilde{U}})).

It can be seen that these angles are independent of the orthonormal basis matrices $U$ and $\widetilde{U}$ which are not unique.

We sometimes place a matrix in one of or both arguments of $\theta_{j}(\,\cdot\,,\,\cdot\,)$ and $\Theta(\,\cdot\,,\,\cdot\,)$ with an understanding that it is about the subspace spanned by the columns of the matrix argument.

It is known that $\|\sin\Theta({\cal U},{\cal\widetilde{U}})\|_{2}$ defines a distance metric between ${\cal U}$ and ${\cal\widetilde{U}}$ [24, p.95].

The next two lemmas and their proofs are about how to pick up two bi-orthogonal basis matrices of two subspaces with acute canonical angles. The results provide a foundation to some of our argument in the paper.

Lemma A.1.

Let ${\cal X}_{1}$ and ${\cal Y}_{1}$ be two subspaces with dimension $r$ of ${\mathbb{R}}^{n}$ . Then

\|\sin\Theta({\cal X}_{1},{\cal Y}_{1})\|_{2}<1

if and only if $Y_{1}^{\operatorname{T}}X_{1}$ is nonsingular for any two basis matrices $X_{1},\,Y_{1}\in\mathbb{R}^{n\times r}$ of ${\cal X}_{1}$ and ${\cal Y}_{1}$ , respectively.

Proof.

Suppose that $\|\sin\Theta({\cal X}_{1},{\cal Y}_{1})\|_{2}<1$ , and let $X_{1},\,Y_{1}\in\mathbb{R}^{n\times r}$ be basis matrices of ${\cal X}_{1}$ and ${\cal Y}_{1}$ , respectively. Then

U=X_{1}(X_{1}^{\operatorname{T}}X_{1})^{-1/2},\quad\widetilde{U}=Y_{1}(Y_{1}^{\operatorname{T}}Y_{1})^{-1/2}

(A.1)

are two orthonormal basis matrices of ${\cal X}_{1}$ and ${\cal Y}_{1}$ , respectively. The singular values of $\widetilde{U}^{\operatorname{T}}U$ are $\cos\theta_{j}({\cal X}_{1},{\cal Y}_{1})$ for $1\leq j\leq r$ which are positive because $\|\sin\Theta({\cal X}_{1},{\cal Y}_{1})\|_{2}<1$ , and hence $\widetilde{U}^{\operatorname{T}}U$ is nonsingular, and since

\widetilde{U}^{\operatorname{T}}U=(Y_{1}^{\operatorname{T}}Y_{1})^{-1/2}Y_{1}^{\operatorname{T}}X_{1}(X_{1}^{\operatorname{T}}X_{1})^{-1/2},

(A.2)

$Y_{1}^{\operatorname{T}}X_{1}$ is nonsingular.

Conversely, let $X_{1},\,Y_{1}\in\mathbb{R}^{n\times r}$ be basis matrices of ${\cal X}_{1}$ and ${\cal Y}_{1}$ , respectively, and suppose that $Y_{1}^{\operatorname{T}}X_{1}$ is nonsingular. Set $U$ and $\widetilde{U}$ as in (A.1). Then $\widetilde{U}^{\operatorname{T}}U$ is nonsingular by (A.2), which means $\cos\theta_{j}({\cal X}_{1},{\cal Y}_{1})>0$ for $1\leq j\leq r$ , implying

\|\sin\Theta({\cal X}_{1},{\cal Y}_{1})\|_{2}=\max_{j}\sqrt{1-\cos^{2}\theta_{j}({\cal X}_{1},{\cal Y}_{1})}<1,

as was to be shown. ∎

Lemma A.2.

Let ${\cal X}_{1}$ and ${\cal Y}_{1}$ be two subspaces with dimension $r$ of ${\mathbb{R}}^{n}$ and suppose that $\|\sin\Theta({\cal X}_{1},{\cal Y}_{1})\|_{2}<1$ , i.e., the canonical angles between the two subspaces are acute.

(a)

There exist basis matrices $X_{1},\,Y_{1}\in\mathbb{R}^{n\times r}$ of ${\cal X}_{1}$ and ${\cal Y}_{1}$ , respectively, such that $Y_{1}^{\operatorname{T}}X_{1}=I_{r}$ ;
(b)

Given a basis matrix $X_{1}\in\mathbb{R}^{n\times r}$ of ${\cal X}_{1}$ , there exists a basis matrix $Y_{1}\in\mathbb{R}^{n\times r}$ of ${\cal Y}_{1}$ such that $Y_{1}^{\operatorname{T}}X_{1}=I_{r}$ ;

(c)

Given basis matrices $X_{1},\,Y_{1}\in\mathbb{R}^{n\times r}$ of ${\cal X}_{1}$ and ${\cal Y}_{1}$ , respectively, such that $Y_{1}^{\operatorname{T}}X_{1}=I_{r}$ , there exist matrices $X_{2},\,Y_{2}\in\mathbb{R}^{n\times(n-r)}$ such that

[Y_{1},Y_{2}]^{\operatorname{T}}[X_{1},X_{2}]=\begin{bmatrix}Y_{1}^{\operatorname{T}}X_{1}&Y_{1}^{\operatorname{T}}X_{2}\\ Y_{2}^{\operatorname{T}}X_{1}&Y_{2}^{\operatorname{T}}X_{2}\end{bmatrix}=I_{n}.

Proof.

For item (a), first we pick two orthonormal basis matrices $U,\,\widetilde{U}\in\mathbb{R}^{n\times r}$ of ${\cal X}_{1}$ and ${\cal Y}_{1}$ , respectively. The assumption $\|\sin\Theta({\cal X}_{1},{\cal Y}_{1})\|_{2}<1$ implies that the singular values of $\widetilde{U}^{\operatorname{T}}U$ are $\cos\theta_{j}({\cal X}_{1},{\cal Y}_{2})$ for $1\leq j\leq r$ are positive, and hence $\widetilde{U}^{\operatorname{T}}U$ is nonsingular. Now take $X_{1}=U(\widetilde{U}^{\operatorname{T}}U)^{-1}$ and $Y_{1}=\widetilde{U}$ .

For item (b), we note $U=X_{1}(X_{1}^{\operatorname{T}}X_{1})^{-1/2}$ is an orthonormal basis matrix of ${\cal X}_{1}$ . Let $\widetilde{U}$ be an orthonormal basis matrix of ${\cal Y}_{1}$ . As we just argued,

\widetilde{U}^{\operatorname{T}}U=(\widetilde{U}^{\operatorname{T}}X_{1})(X_{1}^{\operatorname{T}}X_{1})^{-1/2}

is nonsingular, implying $\widetilde{U}^{\operatorname{T}}X_{1}$ is nonsingular. Now take $Y_{1}=\widetilde{U}(\widetilde{U}^{\operatorname{T}}X_{1})^{-\operatorname{T}}$ .

Finally for item (c), let $\widetilde{V},\,V\in\mathbb{R}^{n\times(n-r)}$ be any orthonormal basis matrices of ${\cal X}_{1}^{\bot}$ and ${\cal Y}_{1}^{\bot}$ , the orthogonal complements of ${\cal X}_{1}$ and ${\cal Y}_{1}$ , respectively, i.e.,

\widetilde{V}^{\operatorname{T}}\widetilde{V}=V^{\operatorname{T}}V=I_{n-r},\quad\widetilde{V}^{\operatorname{T}}X_{1}=V^{\operatorname{T}}Y_{1}=0.

We claim that $X:=[X_{1},V]\in\mathbb{R}^{n\times n}$ is nonsingular; otherwise there exists

0\neq\boldsymbol{x}=\begin{bmatrix}\boldsymbol{z}\\ \boldsymbol{y}\end{bmatrix},\quad\boldsymbol{z}\in\mathbb{R}^{r},\,\,\boldsymbol{y}\in\mathbb{R}^{n-r}

such that $X\boldsymbol{x}=0$ , i.e., $X_{1}\boldsymbol{z}+V\boldsymbol{y}=0$ , pre-multiplying which by $Y_{1}^{\operatorname{T}}$ leads to $\boldsymbol{z}=0$ , which implies $V\boldsymbol{y}=0$ , which implies $\boldsymbol{y}=0$ because $V$ is an orthonormal basis matrix of ${\cal Y}_{1}^{\bot}$ , which says $\boldsymbol{x}=0$ , a contradiction. Similarly, we know $Y:=[Y_{1},\widetilde{V}]\in\mathbb{R}^{n\times n}$ is nonsingular, and so is

Y^{\operatorname{T}}X=[Y_{1},\widetilde{V}]^{\operatorname{T}}[X_{1},V]=\begin{bmatrix}Y_{1}^{\operatorname{T}}X_{1}&Y_{1}^{\operatorname{T}}V\\ \widetilde{V}^{\operatorname{T}}X_{1}&\widetilde{V}^{\operatorname{T}}V\end{bmatrix}=\begin{bmatrix}I_{r}&0\\ 0&\widetilde{V}^{\operatorname{T}}V\end{bmatrix},

implying $\widetilde{V}^{\operatorname{T}}V$ is nonsingular. Now take $X_{2}=V(\widetilde{V}^{\operatorname{T}}V)^{-1}$ and $Y_{2}=\widetilde{V}$ . ∎

Appendix B Perturbation for Lyapunov equation

In this section, we will establish a lemma on the change of the solution to

A^{\operatorname{H}}X+XA+W=0

(B.1)

subject to perturbations to $A$ and $W$ , along the technical line of [13], where $W$ may not necessarily be Hermitian. It is known as the Lyapunov equation if $W$ is Hermitian, but here it may not be. The result in the lemma below is used during our intermediate estimates of transfer function. In conforming to [13], we will state the result for complex matrices: $\mathbb{C}^{n\times n}$ is the set of all $n$ -by- $n$ complex matrices and $A^{\operatorname{H}}$ denotes the complex conjugate of $A$ .

Lemma B.1.

Suppose that $A\in\mathbb{C}^{n\times n}$ is stable, i.e., all of its eigenvalues are located in the left half of the complex plane, and let

K=\int_{0}^{\infty}e^{A^{\operatorname{H}}t}e^{At}dt,

which is the unique solution to the Lyapunov equation $A^{\operatorname{H}}X+XA+I_{n}=0$ . Let $W\in\mathbb{C}^{n\times n}$ (not necessarily Hermitian) and $X\in\mathbb{C}^{n\times n}$ is the solution to the matrix equation (B.1). Perturb $A$ and $W$ to $A+\Delta A_{i}$ ( $i=1,2$ ) and $W+\Delta W$ , respectively, and suppose that the perturbed equation

(A+\Delta A_{1})^{\operatorname{H}}(X+\Delta X)+(X+\Delta X)(A+\Delta A_{2})+(W+\Delta W)=0,

(B.2)

has a solution $X+\Delta X$ , where the trivial case either $A=0$ or $W=0$ is excluded. If

\eta:=\|K\|_{2}\sum_{i=1}^{2}{\|\Delta A_{i}\|_{2}}<1,

(B.3)

then for any unitarily invariant norm $\|\cdot\|_{\operatorname{ui}}$

\|\Delta X\|_{\operatorname{ui}}\leq\frac{\|K\|_{2}}{1-\eta}\left(\|\Delta W\|_{\operatorname{ui}}+\|X\|_{\operatorname{ui}}\sum_{i=1}^{2}\|\Delta A_{i}\|_{2}\right).

(B.4)

Equation (B.1) is not necessarily a Lyapunov equation because $W$ is allowed non-Hermitian, not to mention (B.2) for which two different perturbations are allowed to $A$ at its two occurrences. Equation (B.1) has a unique solution $X$ because $A$ is assumed stable, but a solution to the perturbed equation (B.2) is assumed to exist. It is not clear if the assumption (B.3) ensures both $A+\Delta A_{i}$ for $i=1,2$ are stable and thereby guarantees that (B.2) has a unique solution, too, something worthy further investigation.

Proof of Lemma B.1.

Modifying the proof of [13, Theorem 2.1], instead of [13, Ineq. (2.11)] there, we have

\|\Delta X\|_{\operatorname{ui}}\leq\left(\|\Delta W\|_{\operatorname{ui}}+\sum_{i=1}^{2}\|\Delta A_{i}\|_{2}\big{[}\|X\|_{\operatorname{ui}}+\|\Delta X\|_{\operatorname{ui}}\big{]}\right)\|K\|_{2},

yielding (B.4). ∎

In [13] for the case $\Delta A_{i}$ for $i=1,2$ are the same and denoted by $\Delta A$ , under the condition of Lemma B.1 but without assuming (B.3), it is proved that

\frac{\|\Delta X\|_{2}}{\|X+\Delta X\|_{2}}\leq 2\|A+\Delta A\|_{2}\|K\|_{2}\left(\frac{\|\Delta A\|_{2}}{\|A+\Delta A\|_{2}}+\frac{\|\Delta W\|_{2}}{\|W+\Delta W\|_{2}}\right),

(B.5)

elegantly formulated in such a way that all changes are measured relatively. We can achieve the same thing, too. In fact, under the condition of Lemma B.1 but without assuming (B.3), it can be shown that

\frac{\|\Delta X\|_{\operatorname{ui}}}{\|X+\Delta X\|_{\operatorname{ui}}}\leq\sum_{i=1}^{2}\|A+\Delta A_{i}\|_{2}\|K\|_{2}\left(\frac{\|\Delta A_{i}\|_{2}}{\|A+\Delta A_{i}\|_{2}}+\frac{\|\Delta W\|_{\operatorname{ui}}}{\|W+\Delta W\|_{\operatorname{ui}}}\right).

(B.6)

But, as we argued at the beginning, (B.4) is more convenient for us to use in our intermediate estimations.

References

[1] B. D. O. Anderson and J. B. Moore. Linear Optimal Control. Prentice-Hall, Englewood Cliffs, 1971.
[2] A. C. Antoulas. Approximation of Large-Scale Dynamical Systems. Advances in Design and Control. SIAM, Philadelphia, PA, 2005.
[3] A. C. Antoulas, D. C. Sorensen, and S. Gugercin. A survey of model reduction methods for large-scale systems. In Vadim Olshevsky, editor, Structured Matrices in Mathematics, Computer Science, and Engineering I: Proceedings of an AMS-IMS-SIAM joint summer research conference, University of Colorado, Boulder, June 27-July 1, 1999, volume 280 of Contemporary Mathematics, pages 193–219. American Mathematical Society, Providence, Rhode Island, 2001.
[4] A. C. Antoulas, D. C. Sorensen, and Y. Zhou. On the decay rate of Hankel singular values and related issues. Sys. Contr. Lett., 46(5):323–342, August 2002.
[5] Zhaojun Bai. Krylov subspace techniques for reduced-order modeling of large-scale dynamical systems. Appl. Numer. Math., 43:9–44, 2002.
[6] Jonathan Baker, Mark Embree, and John Sabino. Fast singular value decay for Lyapunov solutions with nonnormal coefficients. SIAM J. Matrix Anal. Appl., 36(2):656–668, 2015.
[7] R. H. Bartels and G. W. Stewart. Algorithm 432: The solution of the matrix equation $AX-BX=C$ . Commun. ACM, 8:820–826, 1972.
[8] Bernhard Beckermann and Alex Townsend. Bounds on the singular values of matrices with displacement structure. SIAM Rev., 61(2):319–344, 2019.
[9] Peter Benner, R.-C. Li, and Ninoslav Truhar. On ADI method for Sylvester equations. J. Comput. Appl. Math., 233(4):1035–1045, 2009.
[10] Peter Benner, Mario Ohlberger, Albert Cohen, and Karen Willcox, editors. Model Reduction and Approximation: Theory and Algorithms. Computational Science & Engineering. SIAM, Philadelphia, 2017.
[11] Roland W. Freund. Model reduction methods based on Krylov subspaces. Acta Numer., 12:267–319, 2003.
[12] Serkan Gugercin and Athanasios C. Antoulas. A survey of model reduction by balanced truncation and some new results. Int. J. Control, 77(8):748–766, 2004.
[13] C. Hewer and C. Kenney. The sensitivity of the stable Lyapunov equation. SIAM J. Control Optim., 26:321–344, 1988.
[14] Jing-Rebecca Li and Jacob White. Low-rank solution of Lyapunov equations. SIAM J. Matrix Anal. Appl., 24(1):260–280, 2002.
[15] Jing-Rebecca Li and Jacob White. Low-rank solution of Lyapunov equations. SIAM Rev., 46(4):693–713, 2004.
[16] R.-C. Li. Matrix perturbation theory. In L. Hogben, R. Brualdi, and G. W. Stewart, editors, Handbook of Linear Algebra, page Chapter 21. CRC Press, Boca Raton, FL, 2nd edition, 2014.
[17] R.-C. Li and Z. Bai. Structure-preserving model reduction using a Krylov subspace projection formulation. Comm. Math. Sci., 3(2):179–199, 2005.
[18] R.-C. Li and L.-H. Zhang. Convergence of block Lanczos method for eigenvalue clusters. Numer. Math., 131:83–113, 2015.
[19] Ren-Cang Li, Ninoslav Truhar, and Lei-Hong Zhang. On Stewart’s perturbation theorem for svd. Ann. Math. Sci. Appl., 2024. to appear.
[20] Thilo Penzl. Eigenvalue decay bounds for solutions of Lyapunov equations: the symmetric case. Sys. Contr. Lett., 40(2):139–144, June 2000.
[21] J. Sabino. Solution of Large-Scale Lyapunov Equations via the Block Modified Smith Method. PhD thesis, Rice University, Houston, Texas, 2006.
[22] D.C. Sorensen and A.C. Antoulas. The Sylvester equation and approximate balanced reduction. Linear Algebra Appl., 351-352:671–700, 2002.
[23] G. W. Stewart and J.-G. Sun. Matrix Perturbation Theory. Academic Press, Boston, 1990.
[24] Ji-Guang Sun. Matrix Perturbation Analysis. Academic Press, Beijing, 1987. In Chinese.
[25] Kemin Zhou, John C. Doyle, and Keith Glover. Robust and Optimal Control. Prentice Hall, Upper Saddle River, New Jersey, 1995.

	$\displaystyle\\|\widetilde{S}\\|_{2}=\sqrt{\\|\widetilde{P}\\|_{2}}$	$\displaystyle\leq\\|S\\|_{2}=\sqrt{\\|P\\|_{2}},$	$\displaystyle\quad\big{\\|}\widetilde{R}\big{\\|}_{2}=\sqrt{\\|\widetilde{Q}\\|_{2}}$	$\displaystyle\leq\\|R\\|_{2}=\sqrt{\\|Q\\|_{2}},$		(4.4)
	$\displaystyle\\|E\\|_{2}$	$\displaystyle\leq\sqrt{\epsilon_{1}},$	$\displaystyle\quad\\|F\\|_{2}$	$\displaystyle\leq\sqrt{\epsilon_{2}},$		(4.5)

	$\displaystyle\\|\widetilde{G}-G\\|_{2}$	$\displaystyle=\left\\|\begin{bmatrix}0&\widetilde{S}^{\operatorname{T}}F\\ E^{\operatorname{T}}\widetilde{R}&E^{\operatorname{T}}F\end{bmatrix}\right\\|_{2}$
		$\displaystyle\leq\max\left\{\sqrt{\\|P\\|_{2}\epsilon_{2}},\sqrt{\\|Q\\|_{2}\epsilon_{1}}\right\}+\sqrt{\epsilon_{1}\epsilon_{2}}=:\varepsilon.$		(4.6)

	$\displaystyle\left\\|\begin{bmatrix}0&\widetilde{S}^{\operatorname{T}}F\\ E^{\operatorname{T}}\widetilde{R}&E^{\operatorname{T}}F\end{bmatrix}\right\\|_{2}$	$\displaystyle\leq\left\\|\begin{bmatrix}0&\widetilde{S}^{\operatorname{T}}F\\ E^{\operatorname{T}}\widetilde{R}&0\end{bmatrix}\right\\|_{2}+\left\\|\begin{bmatrix}0&0\\ 0&E^{\operatorname{T}}F\end{bmatrix}\right\\|_{2}$
		$\displaystyle=\max\left\{\\|\widetilde{S}^{\operatorname{T}}F\\|_{2},\\|E^{\operatorname{T}}\widetilde{R}\\|_{2}\right\}+\\|E^{\operatorname{T}}F\\|_{2}$
		$\displaystyle\leq\max\left\{\sqrt{\\|P\\|_{2}\epsilon_{2}},\sqrt{\\|Q\\|_{2}\epsilon_{1}}\right\}+\sqrt{\epsilon_{1}\epsilon_{2}},$


$\displaystyle\\|\check{U}_{1}-U_{1}\\|_{2}$	$\displaystyle=\frac{\sqrt{2}\,\\|\Gamma\\|_{2}}{\big{[}\sqrt{1+\\|\Gamma\\|_{2}^{2}}(\sqrt{1+\\|\Gamma\\|_{2}^{2}}+1)\big{]}^{1/2}}\leq\\|\Gamma\\|_{2}\leq\frac{2\varepsilon}{\underaccent{\bar}{\delta}},$	(4.17a)
$\displaystyle\\|\check{V}_{1}-V_{1}\\|_{2}$	$\displaystyle=\frac{\sqrt{2}\,\\|\Omega\\|_{2}}{\big{[}\sqrt{1+\\|\Omega\\|_{2}^{2}}(\sqrt{1+\\|\Omega\\|_{2}^{2}}+1)\big{]}^{1/2}}\leq\\|\Omega\\|_{2}\leq\frac{2\varepsilon}{\underaccent{\bar}{\delta}},$	(4.17b)

	$\displaystyle\\|\check{\Sigma}_{1}-\Sigma_{1}\\|_{2}$	$\displaystyle\leq\\|\check{U}_{1}-U_{1}\\|_{2}\\|\widetilde{G}\\|_{2}+\\|\widetilde{G}-G\\|_{2}+\\|G\\|_{2}\\|\check{V}_{1}-V_{1}\\|$
		$\displaystyle\leq\left(1+\frac{4\sigma_{1}}{\underaccent{\bar}{\delta}}\right)\varepsilon,$
and
	$\displaystyle\\|\check{\Sigma}_{1}^{-1}-\Sigma_{1}^{-1}\\|_{2}$	$\displaystyle=\\|\check{\Sigma}_{1}^{-1}(\Sigma_{1}-\check{\Sigma}_{1})\Sigma_{1}^{-1}\\|_{2}$
		$\displaystyle\leq\\|\check{\Sigma}_{1}^{-1}\\|_{2}\\|\Sigma_{1}-\check{\Sigma}_{1}\\|_{2}\\|\Sigma_{1}^{-1}\\|_{2}$
		$\displaystyle\leq\frac{1}{\sigma_{r}\underaccent{\bar}{\sigma}_{r}}\left(1+\frac{4\sigma_{1}}{\underaccent{\bar}{\delta}}\right)\varepsilon,$