%\setcounter{chapter}{-1}
\chapter{Mathematical Background}
\label{chp:m_review}
The aim of this chapter is to review some of the basic mathematical techniques that will be applied in later chapters.
\section{Notation}
We start with Cantor's definition of a set: ``A set is any collection of definite, distinguishable objects of our thought, to be conceived as a whole''~\cite[page 1]{maddox88}. The objects, in the definition, are the elements or members of the set. Usually, we will use capital letters to denote sets and lower case letters to denote elements of the set. Thus, if $X$ is a set, then we write $x\in X$ to mean that $x$ is an element of $X$. When an object $x$ is not an element of a set $X$, we write $x \notin X$. Often, we will use the curly bracket notation for sets. For example, $\{x:x<7\}$ denotes the set of all $x$ less than 7. The colon, :, following $x$ is read ``such that''. We will be working with the following sets:

\begin{tabular}{l}
$N=\{1,2,3,\ldots \}$, the set of all positive integers,\\
$Z=\{0,1,-1,2,-2,\ldots \}$, the set of all integers,\\
$\mb{\Bb R}$, the set of all real numbers,\\
$\mb{\Bb C}$, the set of all complex numbers.
\end{tabular}


\section{Vectors}
We define an $n$-dimensional column vector to be an $n$-tuple of numbers
\begin{equation}
\left[\begin{array}{c}
x_1\\
x_2\\
\vdots\\
x_n
\end{array}\right],
\label{eq:n_tuple}
\end{equation}
where $n$ is a positive integer. Thus a vector is an ordered set of scalars. We will use lower-case bold letters, for example $\b x$, to denote such $n$-tuples. We call the numbers $x_1, x_2,\ldots,x_n$ the {\em coordinates} or {\em components} of the vector $\b x$. The number of components of a vector is its dimension. When referring to a particular component of a vector, we use its component, for example $x_i$, or the subscripted name of the vector, for example $(\b x)_i$ to denote $x_i$. The set of all $n$-tuples, that is $n$-vectors, with real coefficients will be denoted $\mb{\Bb R}^n$ and called the real $n$-space. The transpose of the vector $\b x$, denoted $\b x^T$, is an $n$-dimensional row vector
\begin{equation}
\b x^T=\left[\begin{array}{c}
x_1\\
x_2\\
\vdots\\
x_n
\end{array}\right]^T=\left[\begin{array}{cccc}
x_1&x_2\cdots x_n
\end{array}\right].
\label{eq:tran_tuple}
\end{equation}

The sum of two $n$-vectors $\b x$ and $\b y$ is an $n$-vector, denoted $\b x+\b y$, whose $i$-th component is $x_i+y_i$. 
The following rules are satisfied:
\begin{enumerate}
\item $(\b x+\b y)+\b z=\b x+(\b y+\b z)$, that is, vector addition is associative.
\item $\b x+\b y=\b y+\b x$, which means that vector addition is commutative.
\item Let $\b 0=\left[\begin{array}{cccc}
0&0&\cdots&0
\end{array}\right]^T$ denote the zero vector. Then,
\[
\b 0+\b x=\b x+\b 0=\b x.
\]
\item Let $\b x=\left[\begin{array}{cccc}
x_1&x_2&\cdots&x_n
\end{array}\right]^T$ and let $-\b x=\left[\begin{array}{cccc}
-x_1&-x_2&\cdots&-x_n
\end{array}\right]^T$. Then,
\[
\b x+(-\b x)=\b 0.
\]
\end{enumerate}


The product of an $n$-vector $\b x$ by a scalar $c$ is an $n$-th vector, denoted $c\b x$, whose $i$-th component is $cx_i$.


For two vectors $\b x$ and $\b y$ from $\mb{\Bb R}^n$, we define their {\em scalar product}, also called the {\em inner product}, to be
\begin{equation}
\langle \b x,\b y \rangle=\b x^T\b y=x_1y_1+x_2y_2+\cdots+x_ny_n.
\label{eq:scalar_prod}
\end{equation}
The scalar product is commutative,
\[
\langle\b x,\b y\rangle=\langle\b y,\b x\rangle,
\]
and distributivity over vector addition holds,
\[
\langle\b x,\b y+\b z\rangle=\langle\b x,\b y\rangle+\langle\b x,\b z\rangle.
\]
Two vectors $\b x$ and $\b y$ are said to be orthogonal if $\langle \b x,\b y \rangle=0$.

A {\em vector norm}, denoted $\|\cdot\|$, satisfies the following properties:
\begin{itemize}
\item[{\bf VN 1.}] $\|\b x\| > 0$ for any nonzero vector.
\item[{\bf VN 2.}] For any scalar $c$, $\|c \b x\|=|c|\|\b x\|$.
\item[{\bf VN 3.}] For any two vectors $\b x$ and $\b y$, $\|\b x + \b y\| \le \|\b x\|+\|\b y\|$---{\em the triangle inequality}.
\end{itemize}

The {\em Euclidean norm}, or the {\em 2-norm}, of a vector $\b x \in \mb{\Bb R}^n$, denoted $\|\b x\|_2$, is defined as
\begin{equation}
\|\b x\|_2=\sqrt{\langle \b x,\b x\rangle}=\sqrt{x_1^2+x_2^2+\cdots+x_n^2}.
\label{eq:norm_def}
\end{equation}
Note that $\|\b x\|_2^2=\langle \b x,\b x\rangle=\b x^T\b x$. The Euclidean norm is a special class of norms called the {\em H\"{o}lder} or {\em p-norms} defined by
\begin{equation}
\|\b x\|_p=\left(|x_1|^p+\cdots+|x_n|^p \right)^{1/p},\quad p\ge 1.
\label{eq:defofpnorm}
\end{equation}
Other often used norms that are special cases of the above class of norms are the {\em 1-norm} and the {\em infinity norm} defined respectively as
\begin{equation}
\|\b x\|_1=|x_1|+\cdots+|x_n|
\label{eq:defof1norm}
\end{equation}
and
\begin{equation}
\|\b x\|_{\infty}=\max_i |x_i|.
\label{eq:defofinfnorm}
\end{equation}
Most of the time, we will be using the Euclidean norm for vectors. So for the sake of brevity, we will write $\|\b x\|$ to mean $\|\b x\|_2$.
\begin{thm}[Pythagoras]
If $\b x$ and $\b y$ are orthogonal, then
\[
\|\b x+\b y\|^2=\|\b x\|^2+\|\b y\|^2.
\]
\label{thm:pythagoras}
\end{thm}
{\bf Proof}\quad To prove the theorem, we use the above properties of the scalar product to obtain
\begin{eqnarray*}
\|\b x+\b y\|^2 &=& \left( \b x+\b y\right)^T\left(\b x+\b y \right)\\
&=& \|\b x\|^2+2\langle \b x,\b y\rangle +\|\b y\|^2\\
&=& \|\b x\|^2+\|\b y\|^2,
\end{eqnarray*}
since by assumption $\langle \b x,\b y\rangle=0$. The proof of the theorem is complete.\\
\hspace*{5.5in}$\Box$

We now introduce the notion of {\em projection}. Let $\b x$, $\b y$ be two vectors, where $\b y \ne \b 0$. Let $\b z=c\b y$, $c\in \mb{\Bb R}$, be the vector such that the vector $\b x-\b z$ is orthogonal to $\b y$---see Figure~\ref{fig:projection_il} for an illustration of the above construction.
\begin{figure}
\vskip.2in
\setlength{\unitlength}{.07in}
\begin{center}
\begin{picture}(60,40)(0,0)
\thicklines
\put(10,10){\vector(4,1){40}}
\put(51,20){$\b y$}
\put(10,10){\vector(1,1){20.1}}
\put(27,29){$\b x$}
\put(10,10){\vector(4,1){24}}
\put(34,14){$\b z=c\b y$}
\put(34,16){\vector(-1,4){3.5}}
\put(8,9){$\b 0$}
\put(32,28){$\b x-\b z$}
\end{picture}
\end{center}
\caption{Computing the projection of $\b x$ along $\b y$.}
\label{fig:projection_il}
\end{figure}
Thus,
\begin{equation}
\left( \b x-\b z\right)^T\b y=0.
\label{eq:vecs_xyzp}
\end{equation}
Substituting $\b z=c\b y$ into~(\ref{eq:vecs_xyzp}) yields
\begin{equation}
\left( \b x-\b z\right)^T\b y=\left( \b x-c\b y\right)^T\b y=0.
\label{eq:after_subssyzp}
\end{equation}
Hence,
\begin{equation}
c=\frac{\langle\b x,\b y \rangle}{\|\b y\|^2}.
\label{eq:project_coef}
\end{equation}
The number $c$ given by~(\ref{eq:project_coef}) is called the component of $\b x$ along $\b y$. The projection of $\b x$ along $\b y$ is the vector
\begin{equation}
\b z=c\b y=\frac{\langle\b x,\b y \rangle}{\|\b y\|^2}\b y.
\label{eq:pro_vec_al}
\end{equation}
\begin{thm}[Cauchy-Schwarz Inequality]
Let $\b x, \b y \in \mb{\Bb R}^n$. Then,
\begin{equation}
\left\arrowvert \langle \b x,\b y \rangle\right\arrowvert \le \|\b x\| \|\b y\|.
\label{eq:schwarz_ineq}
\end{equation}
\end{thm}
{\bf Proof}\quad If $\b y = \b 0$, then both sides of~(\ref{eq:schwarz_ineq}) are equal to 0, and the assertion holds. Suppose that $\b y\ne \b 0$, and let $c$ denote the component of $\b x$ along $\b y$. Hence, $c$ is given by~(\ref{eq:project_coef}). We then write
\begin{equation}
\b x=\b x-\b z+\b z=\b x-c\b y+c\b y.
\label{eq:aminusaz}
\end{equation}
For an illustration of~(\ref{eq:aminusaz}), consult Figure~\ref{fig:projection_il}. Applying the Pythagoras Theorem to~(\ref{eq:aminusaz}), we obtain
\begin{equation}
\|\b x\|^2=\|\b x-c\b y\|^2+\|c\b y\|^2=\|\b x-c\b y\|^2+c^2\|\b y\|^2.
\label{eq:after_pytag}
\end{equation}
Hence,
\begin{equation}
c^2\|\b y\|^2\le \|\b x\|^2.
\label{eq:pyta_less}
\end{equation}
On the other hand, by~(\ref{eq:project_coef}),
\begin{equation}
c^2\|\b y\|^2=\frac{\langle\b x,\b y \rangle^2}{\|\b y\|^4}\|\b y\|^2=\frac{\langle\b x,\b y \rangle^2}{\|\b y\|^2}.
\label{eq:on_the_other}
\end{equation}
Combining~(\ref{eq:pyta_less}) and~(\ref{eq:on_the_other}) yields
\begin{equation}
\frac{\langle\b x,\b y \rangle^2}{\|\b y\|^2}\le \|\b x\|^2.
\label{eq:com_and_done}
\end{equation}
Multiplying both sides of~(\ref{eq:com_and_done}) by $\|\b y\|^2$ and taking the square root, concludes the proof.\\
\hspace*{5.5in}$\Box$

\begin{thm}[Triangle inequality]
Let $\b x$, $\b y$ be two vectors in $\mb{\Bb R}^n$. Then,
\begin{equation}
\|\b x+\b y\|\le \|\b x\|+\|\b y\|.
\label{eq:triangle_in}
\end{equation}
\label{thm:triangle_in}
\end{thm}
{\bf Proof}\quad Both sides of~(\ref{eq:triangle_in}) are nonnegative. Hence, to prove that~(\ref{eq:triangle_in}) holds it suffices to show that the squares of the both sides satisfy the desired inequality,
\begin{equation}
\left(\b x+\b y \right)^T\left(\b x+\b y \right)\le \left(\|\b x\|+\|\b y\| \right)^2.
\label{eq:triangle_sq}
\end{equation}
Expanding $\left(\b x+\b y \right)^T\left(\b x+\b y \right)$ and applying the Schwarz inequality yields
\begin{eqnarray}
\left(\b x+\b y \right)^T\left(\b x+\b y \right)&=&\|\b x\|^2+2\langle\b x,\b y \rangle +\|\b y\|^2\nonumber\\
&\le& |\b x\|^2+2\|\b x\|\|\b y \| +\|\b y\|^2\nonumber\\
&=& \left(\|\b x\|+\|\b y\| \right)^2,
\label{eq:triangle_sql}
\end{eqnarray}
which completes the proof of the theorem.\\
\hspace*{5.5in}$\Box$

\begin{thm}
Let $\b x$, $\b y$ be two vectors in $\mb{\Bb R}^n$. Then,
\[
\left|\|\b x\|-\|\b y \|\right|\le \|\b x-\b y\|.
\]
\label{thm:abs_twonmrm}
\end{thm}
{\bf Proof}\quad We have
\[
\|\b x\|=\|\left( \b x-\b y\right)+\b y\|\le \|\b x-\b y\|+\|\b y\|.
\]
Hence,
\begin{equation}
\|\b x\|-\|\b y\| \le \|\b x-\b y\|.
\label{eq:first_abs}
\end{equation}
The above implies that
\begin{equation}
\|\b x-\b y\|=\|\b y-\b x\|\ge \|\b y\|-\|\b x\|=-\left(\|\b x\|-\|\b y\| \right).
\label{eq:sec_abs}
\end{equation}
Combining~(\ref{eq:first_abs}) and~(\ref{eq:sec_abs}) yields
\[
\left|\|\b x\|-\|\b y \|\right|\le \|\b x-\b y\|
\]
which completes the proof.\\
\hspace*{5.5in}$\Box$

\section{Matrices and Determinants}
Let $m$, $n$ be two integers such that $m, n \ge 1$. An array of numbers
\[
\left[\begin{array}{cccc}
a_{11} & a_{12} & \cdots & a_{1n}\\
a_{21} & a_{22} & \cdots & a_{2n}\\
\vdots & \vdots & \ddots & \vdots\\
a_{m1} & a_{m2} & \cdots & a_{mn}
\end{array}\right]
\]
is called a {\em matrix}. Thus a matrix is a set of scalars subject to two orderings. We say that it is an $m$ by $n$ matrix, or $m\times n$ matrix. This matrix has $m$ rows and $n$ columns. If $m=n$, then we say that it is a square matrix. An $m\times n$ matrix $\b A$ whose elements are real numbers will be denoted as $\b A \in \mb{\Bb R}^{m\times n}$. If the elements of the matrix $\b A$ are complex numbers, then we will write $\b A \in \mb{\Bb C}^{m\times n}$. We define the zero matrix, denoted $\b O$, to be the matrix such that $a_{ij}=0$ for all $i$, $j$. We define the $n\times n$ identity matrix, denoted $\b I_n$, as the matrix having diagonal components all equal to 1 and all other components all equal to 0.

Given an $m\times n$ matrix $\b A$. The $n\times m$ matrix $\b B$ such that $b_{ji}=a_{ij}$ is called the transpose of $\b A$ and is denoted $\b A^T$. Thus, the transpose of a matrix is obtained by changing rows into columns and vice versa. A matrix $\b A$ which is equal to its transpose, that is, $\b A=\b A^T$, is called symmetric.

We now define the sum of two matrices. Given two matrices $\b A$ and $\b B$ of the same size. We define $\b A+\b B$ to be the matrix whose elements are $a_{ij}+b_{ij}$, that is, we add matrices of the same size componentwise.

Let $c$ be a number and $\b A$ be a matrix. We define $c\b A$ to be the matrix whose components are $c a_{ij}$. Thus, to obtain $c\b A$, we multiply each component of $\b A$ by $c$.

We now define the product of two matrices. Let $\b A$ be an $m\times n$ matrix and $\b B$ be an $n\times s$ matrix. Let $\b a_i^T$ denote the $i$-th row of $\b A$, and $\b b_j$ the $j$-th column of $\b B$, Then, the product of $\b A\b B$ is an $m \times s$ matrix $\b C$ defined as
\begin{equation}
\b C=\b A\b B=\left[\begin{array}{cccc}
\b a_1^T\b b_1 & \b a_1^T\b b_2 & \cdots & \b a_1^T\b b_s\\
\b a_2^T\b b_1 & \b a_2^T\b b_2 & \cdots & \b a_2^T\b b_s\\
\vdots & \vdots & \ddots & \vdots\\
\b a_m^T\b b_1 & \b a_m^T\b b_2 & \cdots & \b a_m^T\b b_s
\end{array}\right].
\label{eq:mat_product}
\end{equation}
It is not always true that $\b A\b B=\b B\b A$, that is, matrix multiplication is not necessarily commutative.

Let $\b A$, $\b B$, and $\b C$ be matrices such that $\b A$, $\b B$ can be multiplied, $\b B$, $\b C$ can be multiplied, and $\b B$, $\b C$ can be added. Then,
\begin {equation}
\b A(\b B+\b C)=\b A\b B+\b A\b C.
\label{eq:mat_distribute}
\end{equation}
The above property is called the distributive law. If $x$ is a number, then
\begin{equation}
\b A(x\b B)=x(\b A\b B).
\label{eq:mat_scalar}
\end{equation}
Furthermore, matrix multiplication is associative since
\begin{equation}
(\b A\b B)\b C=\b A(\b B\b C).
\label{eq:mat_associate}
\end{equation}

Suppose now that $\b A$, $\b B$ can be multiplied. Then,
\begin{equation}
(\b A\b B)^T=\b B^T\b A^T.
\label{eq:transpose_prod}
\end{equation}
Thus the transpose of a product is the product of transposed matrices in reverse order.

Let $\b A$ be an $n\times n$ matrix. An inverse of $\b A$, if it exists, is an $n\times n$ matrix $\b B$ such that
\begin{equation}
\b A\b B=\b B\b A=\b I_n.
\label{eq:mat_inverse}
\end{equation}
If an inverse exists, then it is unique. A square matrix that is invertible is said to be {\em nonsingular}. We denote the inverse of $\b A$ by $\b A^{-1}$. Hence, we have
\[
\b A^{-1}\b A=\b A\b A^{-1}=\b I_n.
\]
The transpose of an inverse is the inverse of the transpose,
\begin{equation}
\left( \b A^{-1}\right)^T=\left( \b A^T\right)^{-1},
\label{eq:mat_trans_inverse}
\end{equation}
and we denote it by $\b A^{-T}$.

We now present the Sherman-Morrison-Woodbury formula for inverting a special kind of matrices.
\begin{lemma}
\label{lemma:smwformula}
If $\b A$, $\b D$, and $\b D^{-1}+\b C\b A^{-1}\b B$ are invertible, then
\[
\left( \b A+\b B\b D\b C\right)^{-1}=\b A^{-1}-\b A^{-1}\b B\left( \b D^{-1}+\b C\b A^{-1}\b B\right)^{-1}\b C\b A^{-1}.
\]
\end{lemma}
{\bf Proof}\quad We prove the Sherman-Morrison-Woodbury formula by verification. We have
\begin{eqnarray*}
&{}&\left( \b A+\b B\b D\b C\right)\left( \b A^{-1}-\b A^{-1}\b B\left( \b D^{-1}+\b C\b A^{-1}\b B\right)^{-1}\b C\b A^{-1}\right)\\
&{}&\quad= \left( \b A+\b B\b D\b C\right)\b A^{-1}\\
&{}&\qquad -\left( \b A+\b B\b D\b C\right)\b A^{-1}\b B\left(\b D^{-1}+\b C\b A^{-1}\b B \right)^{-1}\b C\b A^{-1}\\
&{}&\quad=\b I+\b B\b D\b C\b A^{-1}\\
&{}&\qquad -\left(\b B+\b B\b D\b C\b A^{-1}\b B \right)\left(\b D^{-1}+\b C\b A^{-1}\b B \right)^{-1}\b C\b A^{-1}\\
&{}&\quad=\b I+\b B\b D\b C\b A^{-1}\\
&{}&\qquad -\b B\b D\left(\b D^{-1}+\b C\b A^{-1}\b B \right)\left(\b D^{-1}+\b C\b A^{-1}\b B \right)^{-1}\b C\b A^{-1}\\
&{}&\quad=\b I+\b B\b D\b C\b A^{-1}-\b B\b D\b C\b A^{-1}\\
&{}&\quad=\b I,
\end{eqnarray*}
as desired.\\
\hspace*{5.5in}$\Box$

In the case when $\b D=\b I$, the Sherman-Morrison-Woodbury formula takes the form
\begin{equation}
\left( \b A+\b B\b C\right)^{-1}=\b A^{-1}-\b A^{-1}\b B\left( \b I+\b C\b A^{-1}\b B\right)^{-1}\b C\b A^{-1}.
\label{eq:specialsmw}
\end{equation}
Formula~(\ref{eq:specialsmw}) is sometimes referred to as the {\em matrix inversion formula}.


A set of nonzero vectors $\{\b q_1,\b q_2,\ldots,\b q_n\}$ is said to be orthogonal if
\[
\langle \b q_i,\b q_j \rangle = 0,\quad\mb{for}\quad i \ne j.
\]
If each vector in an orthogonal set is of unity length, that is, $\|\b q_i\|=1$ for all $i=1,2\ldots,n$, then the set of vectors is called {\em orthonormal}. A square matrix with orthonormal columns will be called a {\em unitary matrix}. Sometimes in the literature, a unitary matrix is called an {\em orthogonal matrix}. If $\b Q\in \mb{\Bb R}^{n \times n}$ is a unitary matrix, then
\begin{equation}
\b Q^T\b Q=\b Q\b Q^T=\b I_n.
\label{eq:orthogonal_mat}
\end{equation}

For a square, $n \times n$, matrix $\b A$, its trace, denoted $\mb{trace}\b A$, is the sum of its diagonal elements,
\begin{equation}
\mb{trace} \b A=\sum_{i=1}^n a_{ii}.
\label{eq:mat_trace}
\end{equation}
If $\b B$ is also an $n \times n$ matrix, then
\begin{equation}
\mb{trace}(\b A\b B)=\mb{trace}(\b B\b A).
\label{eq:two_mat_trace}
\end{equation}

The {\em determinant\/} of a square matrix $\b A$, denoted $\mb{det}\b A$, can be evaluated using the Laplace expansion as follows. Let $c_{ij}$ denote the cofactor corresponding to the element $a_{ij}$ of an $n\times n$ matrix $\b A$. The cofactor $c_{ij}$ is $(-1)^{i+j}$ times the determinant of the $(n-1)\times (n-1)$ submatrix of $\b A$ that is obtained by deleting the $i$-th row and $j$-th column of $\b A$. Then, for any fixed $i$, $1\le i \le n$, we have
\begin{equation}
\det\b A=\sum_{j=1}^n a_{ij}c_{ij}.
\label{eq:mat_det}
\end{equation}
This is the expansion of the determinant along the $i$-th row. In a similar fashion, we can expand the determinant along a column. The determinant has the following properties:
\begin{itemize}
\item[{\bf D 1.}] $\det \b A = \det \b A^T$
\item[{\bf D 2.}] $\det\left( \b A\b B\right)= \det\b A \det \b B = \det \left(\b B\b A \right)$
\item[{\bf D 3.}] $\det \b O = 0$
\item[{\bf D 4.}] $\det \b I_n = 1$.
\end{itemize}

If $\b A$ is $n \times n$ and $\b v$ is a nonzero $n \times 1$ vector such that for some scalar $\lambda$,
\begin{equation}
\b A\b v=\lambda \b v,
\label{eq:eigen_eigen}
\end{equation}
then $\b v$ is an {\em eigenvector\/} corresponding to the {\em eigenvalue\/} $\lambda$. In order that $\lambda$ be an eigenvalue of $\b A$, it is necessary and sufficient for the matrix $\lambda \b I_n-\b A$ to be singular, that is, $\det\left(\lambda \b I_n - \b A \right)=0$. This leads to an $n$-th order polynomial equation in $\lambda$ of the form
\begin{equation}
\det\left(\lambda \b I_n - \b A \right)=\lambda^n+ a_{n-1}\lambda^{n-1}+\cdots+a_1\lambda+a_0=0.
\label{eq:character_eq}
\end{equation}
This equation, called the characteristic polynomial equation, must have $n$, possibly nondistinct, complex roots which are the eigenvalues of $\b A$. The {\em spectral radius\/} of $\b A$, denoted $\rho(\b A)$ is defined as
\[
\rho(\b A)=\max_{1\le i\le n}|\lambda_i(\b A)|,
\]
where $\lambda_i(\b A)$ denotes the $i$-th eigenvalue of $\b A$.

Any $m\times n$ matrix $\b A$ can be represented as
\begin{equation}
\b A=\b U\b S\b V^T,
\label{eq:svddecompose}
\end{equation}
where $\b U$ is an $m\times m$ unitary matrix, $\b V$ is an $n\times n$ unitary matrix, and $\b S$ is an $m\times n$ diagonal matrix of the form
\[
\b S=\mb{diag}(\sigma_1,\sigma_2,\ldots,\sigma_p),\quad\mb{where}\quad p=\min(m,n)\quad\mb{and} \quad\sigma_i\ge 0,\quad i=1,2,\ldots,p.
\]
The nonnegative numbers $\{\sigma_i\}$ are called the {\em singular values\/} of $\b A$, and~(\ref{eq:svddecompose}) is called the {\em singular value decomposition\/}~(SVD) of $\b A$. Following the convention, we will write
\[
\sigma_1\ge \sigma_2\ge \cdots\ge \sigma_p,
\]
so that $\sigma_1=\sigma_1(\b A)$ denotes the largest singular value of $\b A$. The singular values of $\b A$ are the square roots of the eigenvalues of $\b A^T\b A$ if $m\ge n$ or of $\b A\b A^T$ if $m<n$. If $\b A$ is symmetric, its singular values are the absolute values of its eigenvalues.


\section{Quadratic Forms}
A real {\em quadratic form\/} is the product $\b x^T\b Q\b x$, where $\b Q$ is a specified $n\times n$ matrix with real coefficients and $\b x$ is an $n\times 1$ real, column, vector. Without loss of generality $\b Q$ can be taken to be symmetric, that is, $\b Q=\b Q^T$. Indeed, if $\b Q$ is not symmetric, we can always replace it with a symmetric matrix without changing the quadratic form values. Indeed, since $z=\b x^T\b Q\b x$ is a scalar, hence $z^T=z$, and we have
\[
\left( \b x^T\b Q\b x\right)^T=\b x^T\b Q^T\b x=\b x^T\b Q\b x.
\]
Therefore, we can write
\[
\b x^T\b Q\b x=\b x^T\left( \frac{1}{2}\b Q+\frac{1}{2}\b Q^T\right)\b x=\b x^T\left(\frac{\b Q+\b Q^T}{2} \right)\b x,
\]
where $\frac{1}{2}\left(\b Q+\b Q^T \right)$ is a symmetric matrix. Thus, the quadratic form values are unchanged if $\b Q$ is replaced by the symmetric matrix $\frac{1}{2}\left(\b Q+\b Q^T \right)$. From now on, in our study of quadratic forms, we assume that $\b Q$ is symmetric.

A quadratic form $\b x^T\b Q\b x$, or equivalently, the matrix $\b Q$ is called {\em positive semi-definite\/} if $\b x^T\b Q\b x \ge 0$ for all $\b x$. It is {\em positive definite\/} if $\b x^T\b Q\b x > 0$ for all $\b x\ne \b 0$. If $\b Q$ is positive semi-definite, then we write $\b Q\ge 0$, if it is positive definite, we write $\b Q>0$. We write $\b Q_1 \ge \b Q_2$ to mean that $\b Q_1-\b Q_2 \ge 0$. We say that $\b Q$ is {\em negative semi-definite\/} if $\b x^T\b Q\b x \le 0$ for all $\b x$, it is {\em negative definite\/} if $\b x^T\b Q\b x < 0$ for all $\b x \ne \b 0$. Note that $\b Q \le 0$ if and only if $-\b Q \ge 0$.

We now present tests for definiteness properties of symmetric matrices. We state these tests in terms of {\em principal minors\/} of a given real, symmetric, $n \times n$ matrix $\b Q$ with elements $q_{ij}$. We will also give tests for definiteness in terms of eigenvalues of $\b Q$. The principal minors of $\b Q$ are $\det \b Q$ itself and the determinants of submatrices of $\b Q$ obtained by removing successively an $i$-th row and $i$-th column. Specifically, for $p=1,2,\ldots,n$, the principal minors are
\[
\Delta_p\left(\begin{array}{cccc}
i_1, & i_2, &\ldots, & i_p\\
i_1, & i_2, &\ldots, & i_p
\end{array} \right)=\det\left[\begin{array}{cccc}
q_{i_1 i_1}&q_{i_1 i_2}&\cdots&q_{i_1 i_p}\\
q_{i_2 i_1}&q_{i_2 i_2}&\cdots&q_{i_2 i_p}\\
\vdots& \vdots & \ddots & \vdots\\
q_{i_p i_1}&q_{i_p i_2}&\cdots&q_{i_p i_p}
\end{array}\right].
\]
Since $\b Q=\b Q^T$, its eigenvalues are all real.
\begin{thm}
\label{thm:testforpsd}
A symmetric matrix is positive definite (or positive semi-definite)
if and only if its eigenvalues are all positive (nonnegative).
\end{thm}

If a given matrix $\b Q\in \mb{\Bb R}^{n\times n}$ is symmetric, then its eigenvalues are all real and there exists a unitary matrix $\b U$ such that
\begin{equation}
\b U^{T}\b Q\b U=\left[\begin{array}{cccc}
\lambda_1 & 0 & \cdots &0\\
0 & \lambda_2 &\cdots & 0\\
\vdots&{}&\ddots& \vdots\\
0 & 0& \cdots& \lambda_n
\end{array}\right]=\b{\Lambda}.
\label{eq:diaofsymauutr}
\end{equation}
Using the above, we can easily show that a symmetric positive semi-definite matrix $\b Q$ has a positive semi-definite symmetric square root, denoted $\b Q^{1/2}$ or $\sqrt{\b Q}$, that satisfies
\begin{equation}
\b Q^{1/2}\b Q^{1/2}=\b Q.
\label{eq:squmatrpro}
\end{equation}
Indeed, it follows from Theorem~\ref{thm:testforpsd} that if $\b Q=\b Q^T \ge 0$ then its eigenvalues are all nonnegative. Hence,
\begin{equation}
\b{\Lambda}^{1/2}=\left[\begin{array}{cccc}
\lambda_1^{1/2} & 0 & \cdots &0\\
0 & \lambda_2^{1/2} &\cdots & 0\\
\vdots&{}&\ddots& \vdots\\
0 & 0& \cdots& \lambda_n^{1/2}
\end{array}\right]
\label{eq:sqroooflam}
\end{equation}
is well defined. Using the fact that $\b U^T\b U=\b I_n$, we can define
\begin{equation}
\b Q^{1/2}=\b U\b{\Lambda}^{1/2}\b U^T,
\label{eq:deofsqroqpsd}
\end{equation}
which has the desired properties.

If $\b Q=\b Q^T$ is positive semi-definite, and not positive definite, then some of its eigenvalues are equal to zero. Suppose that $\b Q$ has $r$ nonzero eigenvalues. Then, there is a unitary matrix, say $\b U$, such that
\begin{equation}
\b U^{T}\b Q\b U=\left[\begin{array}{ccccccc}
\lambda_1 & 0 & \cdots&0&0 &\cdots&0\\
0 & \lambda_2 &\cdots&0&0 &\cdots& 0\\
\vdots&{}&\ddots& \vdots&{}&\ddots&\vdots\\
0 & 0& \cdots& \lambda_r&0&\cdots&0\\
0&0&\cdots&0&0&\cdots&0\\
\vdots&{}&{}&{}&{}&\ddots&\vdots\\
0&0&\cdots&0&0&\cdots&0
\end{array}\right].
\label{eq:diapsdmauum}
\end{equation}
Let $\b V\in \mb{\Bb R}^{n\times r}$ be a matrix that consists of the first $r$ columns of the matrix $\b U$, and let
\begin{equation}
\b C=\left[\begin{array}{cccc}
\lambda_1^{1/2} & 0 & \cdots &0\\
0 & \lambda_2^{1/2} &\cdots & 0\\
\vdots&{}&\ddots& \vdots\\
0 & 0& \cdots& \lambda_r^{1/2}
\end{array}\right]
\b V^T.
\label{eq:srtmofpsd}
\end{equation}
We obtain a full rank factorization of $\b Q$,
\begin{equation}
\b Q=\b C^T\b C.
\label{eq:furafaofq}
\end{equation}
We can generalize the above, using~(\ref{eq:svddecompose}), to obtain a full rank factorization of a non-square matrix. Specifically, given an $m\times n$ matrix $\b A$ of rank $r$, then we can represent it as
\[
\b A=\b B\b C,
\]
where $\b B$ is an $m\times r$ matrix and $\b C$ is a $r\times n$ matrix such that
\[
\mb{rank}(\b B)=\mb{rank}(\b C)=\mb{rank}(\b A)=r.
\]

\section{The Kronecker Product}
\label{sec:kroneckerpro}
Let $\b A$ be an $m \times n$ and $\b B$ be a $p \times q$ matrices. The {\em Kronecker product\/} of $\b A$ and $\b B$ is denoted $\b A \otimes \b B$ and is a $mp \times nq$ matrix defined by
\begin{equation}
\b A \otimes \b B=\left[\begin{array}{cccc}
a_{11}\b B & a_{12} \b B&\cdots& a_{1n}\b B\\
a_{21}\b B & a_{22} \b B&\cdots& a_{2n}\b B\\
\vdots& \vdots & \ddots & \vdots\\
a_{m1}\b B & a_{m2} \b B&\cdots& a_{mn}\b B
\end{array}\right].
\label{eq:kronecker_def}
\end{equation}
Thus, the matrix $\b A \otimes \b B$ consists of $mn$ blocks. Using the definition of the Kronecker product, we can verify the following two properties:
\begin{eqnarray}
\left(\b A \otimes \b B\right)\left(\b C \otimes \b D \right) &=& \b A \b C \otimes \b B \b D,\\
\left(\b A \otimes \b B\right)^T &=& \b A^T \otimes \b B^T.
\end{eqnarray}

Let
\[
\b X=\left[\begin{array}{cccc}
\b x_1,&\b x_2,&\ldots,&\b x_m
\end{array}\right]
\]
be a $n \times m$ matrix, where $\b x_i$, $i=1,2,\ldots,m$ are the columns of $\b X$. Then, the {\em stacking operator\/} is defined as
\begin{eqnarray}
s(\b X) &=& \left[\begin{array}{cccc}
\b x^T_1,&\b x^T_2,&\ldots,&\b x^T_m
\end{array}\right]^T\nonumber\\
&=&\left[\begin{array}{cccccccccccc}
x_{11}&x_{21}&\cdots&x_{n1}&x_{12}&x_{22}&\cdots&x_{1m}&\cdots&x_{1n}&\cdots &x_{nm}
\end{array}\right]^T
\end{eqnarray}
is the column $nm$-vector formed from the columns of $\b X$ taken in order. Let now $\b A$ be an $n\times n$, and $\b C$ and $\b X$ be $n\times m$ matrices. Then, the matrix equation
\begin{equation}
\b A \b X = \b C
\label{eq:axc_eq}
\end{equation}
can be written as
\begin{equation}
\left(\b I_m \otimes \b A \right) s(\b X)=s(\b C).
\end{equation}
Let $\b B$ be an $m \times m$ matrix, then the matrix equation
\begin{equation}
\b X \b B =\b C
\label{eq:xbc_eq}
\end{equation}
can be written as
\begin{equation}
\left(\b B^T \otimes \b I_n \right) s(\b X)=s(\b C).
\end{equation}
Using the above two facts, we can verify that the matrix equation
\begin{equation}
\b A \b X + \b X \b B = \b C,
\label{eq:abc_eq}
\end{equation}
where $\b A$ is $n \times n$, $\b B$ is $m\times m$, and $\b C$ is $n\times m$, can be written as
\begin{equation}
\left( \b I_m \otimes \b A + \b B^T \otimes \b I_n \right) s(\b X)=s(\b C).
\label{eq:abc_equieq}
\end{equation}

Let $\lambda_i$, $\b w_i$ be the eigenvalues and eigenvectors, respectively of $\b A$, and $\mu_j$ and $\b v_j$ the eigenvalues and eigenvectors of $m\times m$ matrix $\b B$. Then,
\begin{eqnarray}
\left( \b A \otimes \b B \right) \left(\b w_i \otimes \b v_j \right) &=& \b A \b w_i \otimes \b B \b v_j\nonumber\\
&=& \lambda _i \b w_i \otimes \mu _j \b v_j\nonumber\\
&=& \lambda _i \mu _j \left( \b w_i \otimes \b v_j \right).
\end{eqnarray}
Thus, the eigenvalues of $\b A \otimes \b B$ are $\lambda _i \mu_j$, and their respective eigenvectors are $\b w_i \otimes \b v_j$ for $i=1,2,\ldots ,n$, $j=1,2,\ldots ,m$.

We use the above result to study the matrix equation given by~(\ref{eq:abc_eq}). We represent~(\ref{eq:abc_eq}) as
\begin{equation}
\b D \,s(\b X)=s(\b C),
\label{eq:abc_asd}
\end{equation}
where
\begin{equation}
\b D = \b I_m \otimes \b A + \b B^T \otimes \b I_n.
\end{equation}
The solution to~(\ref{eq:abc_asd}) is unique if, and only if, the $mn \times mn$ matrix $\b D$ is nonsingular. To find the condition for this to hold, consider the following matrix
\begin{equation}
\left(\b I_m + \varepsilon \b B^T \right) \otimes \left( \b I_n + \varepsilon \b A \right)=\b I_m \otimes \b I_n + \varepsilon \b D + \varepsilon^2 \b B^T \otimes \b A
\label{eq:long_matrix}
\end{equation}
whose eigenvalues are
\begin{equation}
\left(1+\varepsilon \mu_j\right) \left(1+\varepsilon \lambda_i \right) = 1+\varepsilon \left( \mu_j+\lambda_i \right) + \varepsilon^2 \mu_j \lambda_i
\end{equation}
because
\[
\lambda_i \left(\b I_n + \b A \right)=1+\lambda_i \left(\b A \right).
\]
Comparing terms in $\varepsilon$, we conclude that the eigenvalues of $\b D$ are $\lambda_i+\mu_j$, $i=1,2,\ldots ,n$, $j=1,2,\ldots ,m$. Hence, $\b D$ is nonsingular if and only if
\begin{equation}
\lambda_i+\mu_j \ne 0.
\end{equation}
The above is the necessary and sufficient condition for the solution $\b X$ of the matrix equation $\b A \b X + \b X \b B = \b C$ to be unique.

Another useful identity is
\begin{equation}
s(\b A\b B\b C)=\left(\b C^T\otimes \b A \right)s(\b B)
\label{eq:stack_three}
\end{equation}

For further information on the subject of the Kronecker product, we refer the reader to Brewer~\cite{brewer78}.




\section{Linear Maps}
A function
\[
\b F:\mb{\Bb R}^n \to \mb{\Bb R}^m
\]
is said to be {\em linear\/} if it satisfies the following two properties:
\begin{itemize}
\item[{\bf LM 1.}] For any $\b x$, $\b y$ in $\mb{\Bb R}^n$, we have
\[
\b F(\b x+\b y)=\b F(\b x)+\b F(\b y),
\]
\item[{\bf LM 2.}] If $c$ is a number, then
\[
\b F(c\,\b x)=c\b F(\b x).
\]
\end{itemize}
Linear functions are sometimes referred to as {\em linear maps\/} or {\em linear operators\/}.

Linear functions from $\mb{\Bb R}^n$ into $\mb{\Bb R}^m$ can be identified with the $m\times n$ real matrices by
\[
\b F:\mb{\Bb R}^n \to \mb{\Bb R}^m\quad \mb{is linear if and only if}\quad \b F(\b x)=\b A\b x,
\]
where $\b A\in \mb{\Bb R}^{m\times n}$.

If $\b F:\mb{\Bb R}^n \to \mb{\Bb R}^m$ is linear, then there exists a constant $L \ge 0$ such that
\begin{equation}
\|\b F(\b x)\| \le L\|\b x\|
\label{eq:lin_bounded}
\end{equation}
for all $\b x \in \mb{\Bb R}^n$. The smallest such constant $L$ for which this is true is denoted by $\|\b F\|$ and is called the norm of the linear map $\b F$. If $\b F$ is the linear map induced by the matrix $\b A$, then we define $\|\b A\|=\|\b F\|$. For any matrix $\b A$ and vector norm $\|\cdot\|$, we can compute the matrix norm using the formula
\begin{equation}
\|\b A\|=\sup_{\b x\ne \b 0}\frac{\|\b A\b x\|}{\|\b x\|}.
\label{eq:mat_norm_def}
\end{equation} 
The matrix norm, that is, the norm of the linear map $\b F$ given by~(\ref{eq:mat_norm_def}) is equivalent to
\begin{equation}
\fbox{$\displaystyle \|\b A\|=\max_{\|\b x\|=1} \|\b A\b x\|$}
\label{eq:linnorm_gen}
\end{equation}
A matrix norm produced using~(\ref{eq:mat_norm_def}) or~(\ref{eq:linnorm_gen}) is called the {\em induced norm}. Note that
\[
\|\b A\b x-\b A\b y\|=\|\b A(\b x-\b y)\|\le \|\b A\|\|\b x-\b y\|,
\]
that is, the linear map, $\b F$, represented by the matrix $\b A$ is uniformly continuous.

We will now comment on the notation used in~(\ref{eq:mat_norm_def}) and~(\ref{eq:linnorm_gen}). First, the same symbol is used for the matrix norm as for the norm of a vector. Second, the norms appearing on the right side of~(\ref{eq:mat_norm_def}) and~(\ref{eq:linnorm_gen}) are the norms of the vectors $\b x \in \mb{\Bb R}^n$, and $\b A\b x \in \mb{\Bb R}^m$. Thus, if $\mb{\Bb R}^n$ is equipped with the vector norm $\|\cdot\|_s$ while $\mb{\Bb R}^m$ is equipped with the vector norm $\|\cdot\|_t$, then~(\ref{eq:linnorm_gen}) should be interpreted as
\[
\|\b A\|=\max_{\|\b x\|_s=1} \|\b A\b x\|_t.
\]
Vector norms and a matrix norm are said to be {\em compatible\/} or {\em consistent\/} if they satisfy,
\begin{equation}
\|\b A\b x\|_t \le \|\b A\|\|\b x\|_s.
\label{eq:compconfn}
\end{equation}

We will now show that~(\ref{eq:linnorm_gen}) is well defined. For this, we need the following lemma.
\begin{lemma}
A vector norm $\| \cdot\|:\mb{\Bb R}^p \to \mb{\Bb R}$ is a uniformly continuous function.
\label{lemma:norm_cont}
\end{lemma}
{\bf Proof}\quad
By Theorem~\ref{thm:abs_twonmrm}, on page~\pageref{thm:abs_twonmrm},
\[
\left|\|\b x\|-\|\b y \|\right|\le \|\b x-\b y\|.
\]
Hence, for any $\varepsilon > 0$, we can take $\delta=\varepsilon$ such that for any $\b x, \b y \in \mb{\Bb R}^p$, if $\|\b x-\b y\|<\delta$ then $\left|\|\b x-\b y \|\right|<\varepsilon$, as was to be shown.\\
\hspace*{5.5in} $\Box$

By the continuity of the vector norm and the Theorem of Weierstrass, a vector $\b x_0$ can be found such that $\|\b x_0\|=1$ and $\|\b A\|=\|\b A\b x_0\|=\max_{\|\b x\|=1} \|\b A\b x\|$. Thus, the expression~(\ref{eq:linnorm_gen}) is well defined.


The matrix norm satisfies the following conditions:
\begin{itemize}
\item[{\bf MN 1.}] $\|\b A\| >0$ if $\b A \ne \b O$ and $\|\b O\| = 0$.
\item[{\bf MN 2.}] $\|c\b A\| = |c|\|\b A\|$, where $c$ is a constant.
\item[{\bf MN 3.}] $\|\b A+\b B\| \le \|\b A\|+\|\b B\|$.
\item[{\bf MN 4.}] $\|\b A\b B\| \le \|\b A\|\|\b B\|$.
\end{itemize}
The most frequently used matrix norms are the $p$-norms
\[
\|\b A\|_p=\max_{\|\b x\|_p=1} \|\b A\b x\|_p.
\]
The induced matrix norms corresponding to the one-, two-, and infinity-norms for vectors are:
\begin{itemize}
\item $\|\b A\|_1=\max_j\|\b a_j\|_1$, that is, $\|\b A\|_1$ is the maximum of absolute column sum;
\item $\|\b A\|_2=\sigma_1(\b A)$, which is the largest singular value of $\b A$;
\item $\|\b A\|_{\infty}=\max_i\|\b a_i^T\|_1$, which is the maximum of absolute row sum.
\end{itemize}
A matrix norm that is not induced by a vector norm is the {\em Frobenius norm}, defined for an $m\times n$ matrix $\b A$ as
\begin{equation}
\|\b A\|_F=\left(\sum_{i=1}^m\sum_{j=1}^n a_{ij}^2 \right)^{1/2}.
\label{eq:frobenius_norm}
\end{equation}
The Frobenius norm satisfies the relation
\begin{equation}
\|\b A\|_F^2=\mb{trace}\left(\b A^T\b A \right).
\label{eq:frob_trace}
\end{equation}
The Frobenius norm is compatible with the Euclidean vector norm, that is,
\begin{equation}
\|\b A\b x\|_2 \le \|\b A\|_F\|\b x\|_2.
\label{eq:frobe_euclide}
\end{equation}



\section{Vector Spaces}
A given set of vectors of the same dimension $\{\b x_1,\b x_2,\ldots,\b x_n \}$ is {\em linearly independent} if the zero vector, $\b 0$, can be written as a non-trivial {\em linear combination} of the given vectors,
\[
\alpha_1\b x_1+\alpha_2\b x_2+\cdots+\alpha_n\b x_n=\b 0,\quad\mb{with}\quad \alpha_i\ne 0\quad\mb{for at least one } i.
\]
In other words, the set of vectors $\{\b x_1,\b x_2,\ldots,\b x_n \}$ is linearly independent if
\[
\alpha_1\b x_1+\alpha_2\b x_2+\cdots+\alpha_n\b x_n=\b 0,\quad\mb{implies}\quad \alpha_1=\alpha_2=\cdots=\alpha_n=0.
\]
A set ${\cal S}$ containing vectors of dimension $n$ is a {\em subspace} of $\mb{\Bb R}^n$ if for any scalars $\alpha$ and $\beta$,
\begin{equation}
\b x, \b y \in {\cal S}\quad\mb{implies}\quad \alpha \b x+\beta \b y \in {\cal S}.
\label{eq:def_subspace}
\end{equation}
It follows from the above that a subspace must contain the zero vector, $\b 0$. The set containing only the zero vector is a subspace. This set, $\{\b 0\}$, is called a {\em trivial} subspace.

A matrix-vector product
\begin{equation}
\b A\b x=\left[\begin{array}{cccc}
\b a_1&\b a_2&\cdots&\b a_n
\end{array}\right]\left[\begin{array}{c}
x_1\\
x_2\\
\vdots\\
x_n
\end{array}\right]=x_1\b a_1+x_2\b a_2+\cdots +x_n\b a_n
\label{eq:mat_vec_prod}
\end{equation}
is a linear combination of the columns of the matrix. Thus, linear dependence of the columns of $\b A$ is equivalent to the condition
\begin{equation}
\b A\b z=\b 0\quad\mb{for some vector}\quad \b z \ne \b 0.
\label{eq:lin_dep_condm}
\end{equation}
Linear independence of the columns of $\b A$ is equivalent to the condition
\begin{equation}
\b A\b z=\b 0\quad\mb{implies}\quad \b z=\b 0.
\label{eq:lin_independm}
\end{equation}
The set of all vectors that are linear combinations of the columns of $\b A$ is called the {\em range\/}, or {\em image\/}, of $\b A$ and will be denoted $\mb{range}(\b A)$. Thus,
\begin{equation}
\mb{range}(\b A)=\{\b v:\b v=\b A\b x \quad\mb{for some}\quad \b x\}.
\label{eq:range_defn}
\end{equation}
For a given $m\times n$ matrix $\b A$ and an $m$-vector $\b b$, we have $\b b\in \mb{range}(\b A)$ if and only if there exists an $n$-vector $\b x$ such that $\b A\b x=\b b$. We say that a system of linear equations represented by $\b A\b x=\b b$ is compatible if $\b b \in \mb{range}(\b A)$.

The range of $\b A$ is a subspace of $\mb{\Bb R}^m$. Indeed, for any $\b b_1, \b b_2 \in \mb{range}(\b A)$, we have $\b b_1=\b A \b x_1$ and $\b b_2=\b A\b x_2$ for some $\b x_1$ and $\b x_2$. Using the fact that $\b A$ represents a linear map, for any scalars $\alpha_1$ and $\alpha_2$,
\begin{equation}
\alpha_1\b b_1+\alpha_2\b b_2=\alpha_1\b A\b x_1+\alpha_2\b A\b x_2=\b A\left( \alpha_1\b x_1+\alpha_2 \b x_2 \right).
\label{eq:range_sub}
\end{equation}
Hence, if $\b b_1, \b b_2 \in \mb{range}(\b A)$, then $\alpha_1\b b_1+\alpha_2 \b b_2 \in \mb{range}(\b A)$, which means that the range of $\b A$ is a subspace of $\mb{\Bb R}^m$.

Similarly, we define the range of $\b A^T$ to be the set
\begin{equation}
\mb{range}(\b A^T)=\left\{\b y : \b y=\b A^T\b z \quad\mb{for some}\quad \b z\right\}.
\label{eq:rangeatr}
\end{equation}

For a non-trivial subspace ${\cal S}$ there is a unique smallest positive integer called the {\em dimension\/} of the subspace ${\cal S}$ and will be denoted $\dim {\cal S}$. The dimension of the subspace is equal to the smallest number of linearly independent vectors of the subspace such that every vector of the subspace can be expressed as a linear combination of these vectors. Any such set of vectors is called a {\em basis\/} of the subspace. The vectors in a basis are not unique.

The number of linearly independent columns of a matrix $\b A$ is called the {\em rank\/} of $\b A$ and will be denoted $\mb{rank}(\b A)$. The rank of $\b A$ is also equal to the number of linearly independent rows of $\b A$. For the product of a $m\times n$ matrix $\b A$ and a $n \times p$ matrix $\b B$, we have {\em Sylvester's inequalities}:
\begin{equation}
\mb{rank}(\b A)+\mb{rank}(\b B)-n\le \mb{rank}(\b A\b B)\le \min \{\mb{rank}(\b A),\mb{rank}(\b B) \}.
\label{sylvesterineq}
\end{equation}

For an $m\times n$ matrix $\b A$, the set of all $n$-vectors orthogonal to the rows of $\b A$ is called the {\em kernel\/}, or {\em null space\/}, of $\b A$ and will be denoted $\ker (\b A)$. Hence,
\begin{equation}
\ker (\b A)=\left\{ \b x: \b A\b x=\b 0 \right\}.
\label{eq:kernel_def}
\end{equation}
The kernel of $\b A$ is a subspace of $\mb{\Bb R}^n$. Indeed, if $\b x_1, \b x_2 \in \ker (\b A)$, then $\b A\b x_1 =\b 0$ and $\b A\b x_2 = \b 0$. Hence, for any scalars $\alpha_1$ and $\alpha_2$, we have
\begin{equation}
\b A\left( \alpha_1\b x_1+\alpha_2\b x_2\right)=\alpha_1\b A\b x_1+\alpha_2 \b A\b x_2=\ 0
\label{eq:kernelsubs}
\end{equation}
because $\b A$ represents a linear map. Thus, if $\b x_1, \b x_2 \in \ker (\b A)$, then $\alpha_1\b x_1+\alpha_2\b x_2 \in \ker (\b A)$, which means that $\ker (\b A)$ is a subspace of $\mb{\Bb R}^n$. Similarly, we define the kernel of $\b A^T$ to be the set
\begin{equation}
\ker \left(\b A^T\right)=\left\{ \b w: \b A^T\b w = \b 0\right\}.
\label{eq:kerenlatrans}
\end{equation}
The intersection of $\mb{range}(\b A)$ and $\ker \left( \b A^T \right)$ contains only the zero vector, that is,
\begin{equation}
\mb{range}(\b A) \cap \ker \left( \b A^T\right)=\{\b 0\}.
\label{eq:interseckerandr}
\end{equation}

For a given subspace ${\cal S}$ of $\mb{\Bb R}^n$, the orthogonal complement ${\cal S}^{\perp}$ of ${\cal S}$ is defined by
\begin{equation}
{\cal S}^{\perp}=\left\{\b x\in \mb{\Bb R}^n: \b x^T\b y=0\quad\mb{for all}\quad \b y \in {\cal S} \right\}.
\label{eq:gencomplement}
\end{equation}

The vectors in $\ker \left( \b A^T\right)$ are orthogonal to the vectors in $\mb{range}(\b A)$, and vice versa. We say that $\mb{range}(\b A)$ and $\ker \left(\b A^T \right)$ are {\em orthogonal complements\/} of one another and write
\begin{equation}
\mb{range}(\b A)^{\perp}=\ker\left( \b A^T\right),
\label{eq:complaatr}
\end{equation}
and
\begin{equation}
\ker \left(\b A^T \right)^{\perp}=\mb{range}(\b A).
\label{eq:complatra}
\end{equation}
Similarly,
\begin{equation}
\mb{range}\left(\b A^T\right)^{\perp}=\ker\left( \b A\right),
\label{eq:complatrpa}
\end{equation}
and
\begin{equation}
\ker \left(\b A\right)^{\perp}=\mb{range}\left(\b A^T\right).
\label{eq:complkeraptra}
\end{equation}

If $\b A$ is an $m\times n$ matrix with linearly independent columns, which implies that $m \ge n$, then
\[
\mb{rank}(\b A)=n=\dim \mb{range}(\b A),
\]
and
\[
\dim \ker \left( \b A^T\right)=m-n.
\]

We now define a normed space. The concept of a norm in a vector space is a generalization of the length of a vector in $\mb{\Bb R}^3$.
\begin{defn}
A vector space with a norm is called a {\em normed space}.
\end{defn}
By the normed vector space $\mb{\Bb R}^n$, we usually mean the Euclidean space, which is $\mb{\Bb R}^n$ with the norm $\|\b x\|_2=\sqrt{x_1^2+\cdots +x_n^2}$. The vector space $\mb{\Bb R}^n$ equipped with any of the $p$-norm forms a normed vector space. The absolute value is a norm in $\mb{\Bb R}$.

Let $\Omega$ be a compact subset of $\mb{\Bb R}^n$ and let ${\cal C}(\Omega)$ denote the set of continuous functions $f:\mb{\Bb R}^n \to \mb{\Bb R}$ on $\Omega$. Thus, we view functions as vectors in the function space ${\cal C}(\Omega)$. Then, the function $\|\cdot\|:{\cal C}(\Omega) \to \mb{\Bb R}$ defined as 
\begin{equation}
\|f\|=\max_{\b x \in \Omega}|f(\b x)|,
\label{eq:funnormma}
\end{equation}
satisfies the norm axioms. The function space ${\cal C}(\Omega)$, with scalars from $\mb{\Bb R}$ and with the norm~(\ref{eq:funnormma}) is a linear normed space. In particular, the space of all real-valued functions, $f(x)$, continuous on some interval $[a,b]$ of the real line with the norm defined by
\begin{equation}
\|f\|=\max_{x\in[a,b]}|f(x)|,
\label{eq:rvfofone}
\end{equation}
is an example of a linear normed space. This space is denoted ${\cal C}([a,b])$.

We will also be interested in sequences of elements of ${\cal C}([a,b])$. Let $f_1,f_2,\ldots$ be a sequence of elements of ${\cal C}([a,b])$. Denote this sequence by $\{f_k\}$. In the following, $f_k(x)$, means the value of the function $f_k$ at $x$, while $f_k$ means the function itself. We say that the sequence $\{f_k\}$ converges (pointwise) on the interval $[a,b]$ if there exists a function $f$ such that for every $x\in[a,b]$, the sequence $\{f_k(x)\}$ converges to $f(x)$. In other words, the sequence $\{f_k\}$ converges (pointwise) if there exists a function $f$ such that for a given $x\in[a,b]$ and a given $\varepsilon>0$, we can find $N=N(x,\varepsilon)$ such that
\begin{equation}
|f(x)-f_k(x)|<\varepsilon\quad\mb{for}\quad k>N.
\label{eq:ptconv}
\end{equation}
The sequence $\{f_k\}$ converges uniformly, on the interval $[a,b]$, if there is a function $f$ such that for any $\varepsilon >0$ there exists $N=N(\varepsilon)$, depending on $\varepsilon$ but not on $x$, such that~(\ref{eq:ptconv}) holds for every $x\in[a,b]$. The distinction between convergence and uniform convergence is that in the later case, the same $N$ can be used for any value of $x$ to show the convergence of the sequence $\{f_k(x)\}$. One can show that the limit of an uniformly convergent sequence of continuous functions is continuous. It can be verified that the norm~(\ref{eq:rvfofone}) defines the uniform convergence since $\{f_k\}$ converges uniformly to $f$ if and only if $\max_{x\in[a,b]}|f(x)-f_k(x)|\to 0$. This the reason the norm~(\ref{eq:rvfofone}) is sometimes called the norm of uniform convergence. It can be shown that there is no norm on ${\cal C}([a,b])$ that defines pointwise convergence to an element in ${\cal C}([a,b])$---see Debnath and Mikusi\'{n}ski~\cite[page 12]{debnath90}.

\section{First Order Differential Equation}
Let $\mb{\Bb R}$ denote the set of real numbers and let $D\subseteq \mb{\Bb R}^2$ be a domain, that is, an open connected nonempty subset of $\mb{\Bb R}^2$. Let $f$ be a real valued function defined and continuous on $D$. Let $\dot{x}=dx/dt$ denote the derivative of $x$ with respect to $t$. We call
\begin{equation}
\dot{x}=f(t,x)
\label{eq:gefoode}
\end{equation}
an {\em ordinary differential equation of the first order}. By a {\em solution\/} of the differential equation~(\ref{eq:gefoode}) on an open interval $I=\{t\in \mb{\Bb R}:t_1<t<t_2\}$, we mean a real valued continuously differentiable function $\phi$ defined on $I$ such that $(t,\phi(t))\in D$ for all $t \in I$ and
\[
\dot{\phi}(t)=f(t,\phi(t))\quad\mb{for all}\quad t\in I.
\]
\begin{defn}
Given $(t_0,x_0)\in D$, the {\em initial value problem\/} for~(\ref{eq:gefoode}) is
\begin{equation}
\dot{x}=f(t,x),\quad x(t_0)=x_0.
\label{eq:invaprofsto}
\end{equation}
\end{defn}
A function $\phi$ is a solution of~(\ref{eq:invaprofsto}) if $\phi$ is a solution of the differential equation~(\ref{eq:gefoode}) on some interval $I$ containing $t_0$ and $x(t_0)=x_0$.

We can represent the initial value problem~(\ref{eq:invaprofsto}) equivalently by an integral equation,
\begin{equation}
\phi(t)=x_0+\int_{t_0}^t f(s,\phi(s)) ds.
\label{eq:eqreofinvpr}
\end{equation}
To prove this equivalence, let $\phi$ be a solution of the initial value problem~(\ref{eq:invaprofsto}). Then $x(t_0)=x_0$ and
\[
\dot{\phi}(t)=f(t,\phi(t))\quad\mb{for all}\quad t\in I.
\]
Integrating the above from $t_0$ to $t$, we obtain
\[
\int_{t_0}^t \dot{\phi}(s)ds = \int_{t_0}^t f(s,\phi(s)) ds.
\]
Hence,
\[
\phi(t)-x_0=\int_{t_0}^t f(s,\phi(s)) ds.
\]
Therefore, $\phi$ is a solution of the integral equation~(\ref{eq:eqreofinvpr}).

Conversly, let $\phi$ be a solution of the integral equation~(\ref{eq:eqreofinvpr}). Then $\phi(t_0)=x_0$, and differentiating both sides of~(\ref{eq:eqreofinvpr}) with respect to $t$ gives $\dot{\phi}(t)=f(t,\phi(t))$. Therefore, $\phi$ is also a solution of the initial value problem~(\ref{eq:invaprofsto}).


The general {\em ordinary differential linear equation of the first order\/} has the form
\begin{equation}
\dot{x}(t)=a(t)x(t)+b(t), \qquad t \in (t_1,\; t_2),
\label{eq:first_ode}
\end{equation}
where $a(t)$ and $b(t)$ are assumed to be continuous on the interval $(t_1,\; t_2)$. We associate with equation~(\ref{eq:first_ode}) the initial condition
\[
x(t_0)=x_0, \qquad x_0 \in \mb{\Bb R}.
\]
To find the solution to~(\ref{eq:first_ode}), we write it in the form
\begin{equation}
\dot{x}(t)-a(t)x(t)=b(t).
\label{eq:first_next}
\end{equation}
We then multiply~(\ref{eq:first_next}) by the ``integrating factor''
\begin{equation}
e^{-\int_{t_0}^{t}a(s)ds}
\label{eq:integrating_f}
\end{equation}
to get
\begin{equation}
\dot{x}(t)e^{-\int_{t_0}^{t}a(s)ds}-a(t)x(t)e^{-\int_{t_0}^{t}a(s)ds}=e^{-\int_{t_0}^{t}a(s)ds}b(t).
\label{eq:integrate_multi}
\end{equation}
The above equation can be represented in an equivalent form as
\begin{equation}
\frac{d}{dt}\left \{ x(t)e^{-\int_{t_0}^{t}a(s)ds}\right \}=e^{-\int_{t_0}^{t}a(s)ds}b(t).
\label{eq:first_equi}
\end{equation}
Integrating both sides of~(\ref{eq:first_equi}) yields
\begin{equation}
x(t)e^{-\int_{t_0}^{t}a(s)ds}-x_0=\int_{t_0}^{t}e^{-\int_{t_0}^{\tau}a(s)ds}b(\tau )d\tau .
\end{equation}
Hence,
\begin{equation}
x(t)=e^{\int_{t_0}^{t}a(s)ds}\left [x_0+\int_{t_0}^{t}e^{-\int_{t_0}^{\tau}a(s)ds}b(\tau )d\tau \right ],
\end{equation}
or in an equivalent form
\begin{equation}
\fbox{$\displaystyle x(t)=e^{\int_{t_0}^{t}a(s)ds}x_0+\int_{t_0}^{t}e^{\int_{\tau }^{t}a(s)ds}b(\tau )d\tau $}
\label{eq:first_soln}
\end{equation}

\subsection{The Bellman-Gronwall Lemma}
\begin{lemma}
\label{lem:bellgron}
Let $v(t)$ and $w(t)$ be continuous real functions of $t$. Let $w(t) \ge 0$, and $c$ be a real constant. If
\begin{equation}
v(t) \le c+\int_{0}^{t}w(s)v(s)ds,
\label{eq:v_func_l}
\end{equation}
then
\begin{equation}
v(t) \le ce^{\int_{0}^{t}w(s)ds}.
\end{equation}
\end{lemma}
{\bf Proof}\quad Denote the right-hand side of~(\ref{eq:v_func_l}) by $x(t)$, that is, let
\begin{equation}
x(t)=c+\int_{0}^{t}w(s)v(s)ds,
\label{eq:x_defnit}
\end{equation}
and let
\begin{equation}
z(t)=x(t)-v(t).
\label{eq:z_defn}
\end{equation}
Note that $z(t) \ge 0$, and $x(t)$ is differentiable. Differentiating both sides of~(\ref{eq:x_defnit}) and using~(\ref{eq:z_defn}) yields
\begin{eqnarray}
\dot{x}(t)&=&w(t)v(t)\nonumber\\
&=&w(t)x(t)-w(t)z(t).
\label{eq:xwvz}
\end{eqnarray}
We use~(\ref{eq:first_soln}) to represent~(\ref{eq:xwvz}) in the equivalent integral form as
\begin{equation}
x(t)=e^{\int_{0}^{t}w(s)ds}x_0-\int_{0}^{t}e^{\int_{\tau }^{t}w(s)ds}w(\tau )z(\tau )d\tau .
\label{eq:xwvz_equi}
\end{equation}
Since $z(t) \ge 0$ and $w(t) \ge 0$, we have
\begin{equation}
x(t) \le e^{\int_{0}^{t}w(s)ds}x_0.
\end{equation}
From the definition of $x(t)$ given by~(\ref{eq:x_defnit}), we have $x(0)=c$. Since
\begin{equation}
v(t)=x(t)-z(t) \le x(t),
\end{equation}
because $z(t) \ge 0$, hence,
\[
v(t) \le ce^{\int_{0}^{t}w(s)ds},
\]
and the proof is completed.\\
\hspace*{5.5in} $\Box$

The generalized Bellman-Gronwall lemma allows us to convert the investigation of integral equations into integral inequalities.

\section{Linearizing Differential Equations}
We now consider a dynamic system modeled by a set of nonlinear differential equations
\begin{eqnarray*}
\frac{dx_1}{dt} &=& f_1(x_1,x_2,\ldots,x_n,u_1,u_2,\ldots,u_m)\\
\frac{dx_2}{dt} &=& f_2(x_1,x_2,\ldots,x_n,u_1,u_2,\ldots,u_m)\\
{}&\vdots&{}\\
\frac{dx_n}{dt} &=& f_n(x_1,x_2,\ldots,x_n,u_1,u_2,\ldots,u_m).
\end{eqnarray*}
We assume that the functions $f_i$, $i=1,2,\ldots,n$ are continuously differentiable. The above set of differential equations can be represented in vector form as
\[
\dot{\b x}=\b f(\b x,\b u).
\]
Let $\b u_0=\left[\begin{array}{cccc}
u_{10} & u_{20}&\cdots&u_{m0}
\end{array}\right]^T$ be a constant input that forces the system to settle into a constant, equilibrium, state $\b x_0=\left[\begin{array}{cccc}
x_{10} & x_{20}&\cdots&x_{n0}
\end{array}\right]^T$, that is,
\[
\b f(\b x_0,\b u_0)=\b 0.
\]
We now perturb the equilibrium state by allowing
\[
\b x=\b x_0+\Delta \b x,\quad \b u=\b u_0+\Delta \b u.
\]
Taylor's expansion yields
\begin{eqnarray*}
\frac{d}{dt}\b x &=& \frac{d}{dt}\b x_0+\frac{d}{dt}\Delta \b x\\
&=& \b f(\b x_0+\Delta x,\b u_0+\Delta u)\\
&=& \b f(\b x_0,\b u_0)+\frac{\partial \b f}{\partial \b x}(\b x_0, \b u_0)\,\Delta \b x+\frac{\partial \b f}{\partial \b u}(\b x_0, \b u_0)\,\Delta \b u+\mb{higher order terms},
\end{eqnarray*}
where
\[
\frac{\partial \b f}{\partial \b x}(\b x_0,\b u_0)=\left.\left[\begin{array}{ccc}
\frac{\partial f_1}{\partial x_1}&\cdots&\frac{\partial f_1}{\partial x_n}\\
\vdots & {}&\vdots\\
\frac{\partial f_n}{\partial x_1}&\cdots&\frac{\partial f_n}{\partial x_n}
\end{array}\right]\right\arrowvert_{\b x=\b x_0 \atop
\b u=\b u_0}
\!\!\mb{and}\quad\frac{\partial \b f}{\partial \b u}(\b x_0,\b u_0)=\left.\left[\begin{array}{ccc}
\frac{\partial f_1}{\partial u_1}&\cdots&\frac{\partial f_1}{\partial u_m}\\
\vdots & {}&\vdots\\
\frac{\partial f_n}{\partial u_1}&\cdots&\frac{\partial f_n}{\partial u_m}
\end{array}\right]\right\arrowvert_{\b x=\b x_0 \atop
\b u=\b u_0}
\]
are the Jacobian matrices of $\b f$ with respect $\b x$ and $\b u$ respectively, evaluated at the equilibrium point. Note that
\[
\frac{d}{dt}\b x = \frac{d}{dt}\b x_0+\frac{d}{dt}\Delta \b x=\frac{d}{dt}\Delta \b x
\]
because $\b x_0$ is constant. Furthermore $\b f(\b x_0,\b u_0)=\b 0$. Let
\[
\b A=\frac{\partial \b f}{\partial \b x}(\b x_0,\b u_0)\quad\mb{and}\quad \b B= \frac{\partial \b f}{\partial \b u}(\b x_0,\b u_0).
\]
By neglecting higher order terms we arrive at the linear approximation
\[
\frac{d}{dt}\Delta \b x=\b A\Delta \b x+\b B\Delta \b u.
\]
Similarly, if the outputs of the nonlinear system model are of the form
\begin{eqnarray*}
y_1 &=& h_1(x_1,x_2,\ldots,x_n,u_1,u_2,\ldots,u_m)\\
y_2 &=& h_2(x_1,x_2,\ldots,x_n,u_1,u_2,\ldots,u_m)\\
{}&\vdots&{}\\
y_p &=& h_p(x_1,x_2,\ldots,x_n,u_1,u_2,\ldots,u_m),
\end{eqnarray*}
or in vector notation
\[
\b y=\b h(\b x,\b u),
\]
then Taylor's series expansion can again be used to yield the linear approximation of the above output equations. Indeed, let
\[
\b y=\b y_0+\Delta \b y,
\]
then we obtain
\[
\Delta \b y=\b C \Delta \b x+\b D \Delta \b u,
\]
where
\[
\b C=\frac{\partial \b h}{\partial \b x}(\b x_0,\b u_0)=\left.\left[\begin{array}{ccc}
\frac{\partial h_1}{\partial x_1}&\cdots&\frac{\partial h_1}{\partial x_n}\\
\vdots & {}&\vdots\\
\frac{\partial h_p}{\partial x_1}&\cdots&\frac{\partial h_p}{\partial x_n}
\end{array}\right]\right\arrowvert_{\b x=\b x_0 \atop
\b u=\b u_0}\!\in \mb{\Bb R}^{p\times n},
\]
and
\[
\b D=\frac{\partial \b h}{\partial \b u}(\b x_0,\b u_0)=\left.\left[\begin{array}{ccc}
\frac{\partial h_1}{\partial u_1}&\cdots&\frac{\partial h_1}{\partial u_m}\\
\vdots & {}&\vdots\\
\frac{\partial h_p}{\partial u_1}&\cdots&\frac{\partial h_p}{\partial u_m}
\end{array}\right]\right\arrowvert_{\b x=\b x_0 \atop
\b u=\b u_0}\!\in \mb{\Bb R}^{p\times m}
\]
are the Jacobian matrices of $\b h$ with respect $\b x$ and $\b u$ respectively, evaluated at the equilibrium point $[\b x_0^T\,\b u_0^T]^T$.


\section{Solving the State Equations}
\subsection{Solution of Uncontrolled System}
Given a dynamical system's time-invariant model of the form
\begin{equation}
\dot{\b x}(t)=\b A \b x(t)
\label{eq:xdobaxeq}
\end{equation}
subject to an initial condition
\begin{equation}
\b x(0)=\b x_0,
\label{eq:inicoax}
\end{equation}
where $\b A \in \mb{\Bb R}^{n \times n}$. The solution to~(\ref{eq:xdobaxeq}) subject to~(\ref{eq:inicoax}) is
\[
\fbox{$\displaystyle \b x(t)=e^{\b A t}\b x_0 $}
\]
The above is a generalization of the situation when $\b A$ is a scalar. The exponential matrix $e^{\b A t}$ is defined as
\[
e^{\b At}=\b I_n + t\b A + \frac{t^2}{2!}\b A^2+\frac{t^3}{3!}\b A^3+\cdots
\]
Observe that
\begin{eqnarray*}
\frac{d}{dt}e^{\b At} &=& \b Ae^{\b At}\\
&=& e^{\b At} \b A.
\end{eqnarray*}
Suppose now that the initial condition $\b x(0)=\b x_0$ is replaced by a more general condition
\[
\b x(t_0)=\b x_0,
\]
where $t_0$ is a given initial time instant. Then, the solution to $\dot{\b x}(t)=\b A \b x(t)$ subject to $\b x(t_0)=\b x_0$ is
\[
\fbox{$\displaystyle\b x(t)=e^{\b A\left(t-t_0\right)}\b x_0 $}
\]
The matrix $e^{\b A\left(t-t_0\right)}$ is often written as
\[
e^{\b A\left(t-t_0\right)}=\b \Phi \left(t,t_0\right)
\]
and is called the {\em state transition matrix\/}, since it relates the state at any instant of time $t_0$ to the state at any other time $t$. There are many methods for calculating the state transition matrix. For an overview of some of them, we refer the reader to the classic paper on the subject by Moler and Van Loan~\cite{moler78}. In the case when the matrix $\b A$ is of low dimensions, we can use the formula that results from applying the Laplace transform to the time-invariant system $\dot{\b x}=\b A \b x,\quad \b x(0)=\b x_0$. The Laplace transform of the above matrix differential equation is
\[
s\b X(s)-\b x_0= \b A\b X(s).
\]
Hence
\[
\b X(s)=\left(s\b I_n-\b A\right)^{-1}\b x_0,
\]
and the inverse transform yields
\begin{equation}
e^{\b At}={\cal L}^{-1}\left\{\left(s\b I_n-\b A\right)^{-1}\right\},
\label{eq:inlatrofunco}
\end{equation}
where ${\cal L}^{-1}$ denotes the inverse Laplace transform operator.

\subsection{Solution of Controlled System}
Consider a model of a controlled dynamical system,
\begin{equation}
\dot{\b x}(t)=\b A\b x(t) + \b B\b u(t)
\label{eq:cosymoli}
\end{equation}
subject to a given initial condition $\b x(0)=\b x_0$, where $\b B$ is an $n \times m$ matrix, and $m<n$ as is usually in practice. Premultiplying both sides of~(\ref{eq:cosymoli}) by the integrating factor, $e^{-\b At}$, yields
\[
e^{-\b At}\dot{\b x}(t)=e^{-\b At}\b A \b x(t) + e^{-\b At} \b B\b u(t).
\]
We then rearrange the above equation as
\begin{equation}
e^{-\b At}\dot{\b x}(t)-e^{-\b At}\b A \b x(t) = e^{-\b At} \b B\b u(t).
\label{eq:xdaxburear}
\end{equation}
Noting that
\[
\frac{d}{dt}\left( e^{-\b At}\b x(t) \right)=-\b Ae^{-\b At}\b x(t)+e^{-\b At}\dot{\b x}(t),
\]
we can represent~(\ref{eq:xdaxburear}) in the form
\[
\frac{d}{dt}\left(e^{-\b At}\b x(t)\right)=e^{-\b At}\b B\b u(t).
\]
Integrating both sides of the resulting equation gives
\[
e^{-\b At}\b x(t)-\b x(0)=\int_{0}^{t}e^{-\b A \tau}\b B\b u(\tau)d\tau.
\]
Hence
\[
\b x(t)=e^{\b At}\b x(0)+\int_{0}^{t}e^{\b A(t-\tau)}\b B\b u(\tau)d\tau.
\]
If the initial condition is $\b x(t_0)=\b x_0$, then integration from $t_0$ to $t$ yields
\begin{equation}
\fbox{$\displaystyle \b x(t)=e^{\b A(t-t_0)}\b x(t_0)+\int_{t_0}^{t}e^{\b A(t-\tau)}\b B\b u(\tau)d\tau $}
\label{eq:coreofcotrls}
\end{equation}

As in the case of the uncontrolled system, we can also use the Laplace transform technique to obtain the solution to~(\ref{eq:cosymoli}) subject to the initial condition $\b x(0)=\b x_0$. Indeed, taking the Laplace transform of the above equation and performing some manipulations yields
\[
\b X(s)=\left(sI_n-\b A\right)^{-1}\b x(0)+\left(sI_n-\b A\right)^{-1}\b B\b U(s).
\]
Hence
\[
\b x(t)={\cal L}^{-1}\left\{\left(sI_n-\b A\right)^{-1}\b x(0)\right\}+{\cal L}^{-1}\left\{\left(sI_n-\b A\right)^{-1}\b B\b U(s)\right\}.
\]
If $\b x_0=\b x(t_0)$, then
\[
\b x(t)={\cal L}^{-1}\left\{e^{-t_0s}\left(sI_n-\b A\right)^{-1}\b x(t_0)\right\}+{\cal L}^{-1}\left\{e^{-t_0s}\left(sI_n-\b A\right)^{-1}\b B\b U(s)\right\}.
\]





\section{Matrix Calculus Formulas}
Below useful formulas from matrix calculus are listed that involve computing the gradient (with respect to $\mb{\bm$x$}$). Recall that in our convention the gradient is a column vector. We have:
\begin{equation}
\nabla \left( \mb{\bm$x$}^T\mb{\bm$y$}\right)=\nabla \left( \mb{\bm$y$}^T\mb{\bm$x$}\right)=\mb{\bm$y$}
\label{eq:gralifo}
\end{equation}
\begin{equation}
\nabla \left( \mb{\bm$x$}^T\mb{\bm$x$}\right)=2\mb{\bm$x$}
\label{eq:xtrax}
\end{equation}
\begin{equation}
\nabla \left( \mb{\bm$x$}^T\mb{\bm$A$}\mb{\bm$y$}\right)=\mb{\bm$A$}\mb{\bm$y$}
\label{eq:xtay}
\end{equation}
\begin{equation}
\nabla \left( \mb{\bm$y$}^T\mb{\bm$A$}\mb{\bm$x$}\right)=\mb{\bm$A$}^T\mb{\bm$y$}
\label{eq:ytax}
\end{equation}
\begin{equation}
\nabla \left( \mb{\bm$x$}^T\mb{\bm$A$}\mb{\bm$x$}\right)=\left(\mb{\bm$A$}+\mb{\bm$A$}^T\right)\mb{\bm$x$}
\label{eq:xtaxgrad}
\end{equation}
Let $\mb{\bm$f$}:\mb{\Bb R}^n \to \mb{\Bb R}^n$ and $\mb{\bm$g$}:\mb{\Bb R}^n \to \mb{\Bb R}^n$ be two maps, and let $\frac{\partial \mb{\bm$f$}}{\partial \mb{\bm$x$}}$ and $\frac{\partial \mb{\bm$g$}}{\partial \mb{\bm$x$}}$ be the Jacobian matrices of $\mb{\bm$f$}$ and $\mb{\bm$g$}$, respectively. Then,
\begin{equation}
\nabla \left( \mb{\bm$f$}^T(\mb{\bm$x$})\mb{\bm$g$}(\mb{\bm$x$})\right)=\nabla \left( \mb{\bm$g$}^T(\mb{\bm$x$})\mb{\bm$f$}(\mb{\bm$x$})\right)=\left( \frac{\partial \mb{\bm$f$}}{\partial \mb{\bm$x$}}\right)^T\mb{\bm$g$}+\left( \frac{\partial \mb{\bm$g$}}{\partial \mb{\bm$x$}}\right)^T\mb{\bm$f$}
\label{eq:ftgxpr}
\end{equation}
For an $n\times n$ symmetric matrix $\mb{\bm$Q$}$,
\begin{equation}
\nabla \left( \mb{\bm$f$}^T(\mb{\bm$x$})\mb{\bm$Q$}\mb{\bm$f$}(\mb{\bm$x$})\right)=2\left( \frac{\partial \mb{\bm$f$}}{\partial \mb{\bm$x$}}\right)^T\mb{\bm$Q$}\mb{\bm$f$}
\label{eq:ftqgprgr}
\end{equation}
\section{Notes}
For more information on linear algebra, and in particular matrix theory, the reader is referred Gelfand~\cite{gelfand61} and Gantmacher~\cite{gantmacher59}. Analyses of numerical methods for matrix computations can be found in Faddeeva~\cite{faddeeva59}, Wilkinson~\cite{wilkinson65}, Householder~\cite{householder75}, Golub and Van~Loan~\cite{golub83}, as well as in Gill, Murray, and Wright~\cite{gill91}. Lang~\cite{lang87} and Seeley~\cite{seeley70} are recommended for a review of calculus of several variables. A straightforward approach to mathematical analysis is presented in Binmore~\cite{binmore82}. Principles of real analysis can be found in Natanson~\cite{natanson55}, Bartle~\cite{bartle76}, and Rudin~\cite{rudin76}. Some facts from basic functional analysis  that we touched upon in this chapter, are discussed, for example, by Maddox~\cite{maddox88}, and Debnath and Mikusi\'{n}ski~\cite{debnath90}. Lucid expositions of differential equations are in Driver~\cite{driver77}, Miller and Michel~\cite{miller82}, and in Boyce and DiPrima~\cite{boyce97}.

The English mathematician, G.~H.~Hardy~(1877--1947) argues in~\cite[page~113]{hardy92} that in truly great theorems and their proofs, ``there is a very high degree of {\em unexpectedness\/}, combined with {\em inevitability\/} and {\em economy\/}. The arguments take so odd and surprising a form; the weapons used seem so childishly simple when compared with the far-reaching results; but there is no escape from the conclusions.'' He further states that ``A mathematical proof should resemble a simple and clear-cut constellation, not a scattered cluster in the Milky Way.''
