# Co-design for Security and Performance: LMI Tools

Navid Hashemi<EMAIL_ADDRESS>Justin Ruths<EMAIL_ADDRESS>The
University of Texas at Dallas, Department of Mechanical Engineering  The
University of Texas at Dallas, Departments of Mechanical and Systems
Engineering

###### Abstract

We present a convex optimization to reduce the impact of sensor falsification
attacks in linear time invariant systems controlled by observer-based
feedback. We accomplish this by finding optimal observer and controller gain
matrices that minimize the size of the reachable set of attack-induced states.
To avoid trivial solutions, we integrate a covariance-based $\|H\|_{2}$
closed-loop performance constraint, for which we develop a novel linearization
for this typically nonlinear, non-convex problem. We demonstrate the
effectiveness of this linear matrix inequality framework through a numerical
case study.

###### keywords:

LMI, reachable set, security, robust control.

††thanks: The authors are with the Departments of Mechanical and Systems
Engineering at the University of Texas at Dallas.

, ,

##  1 Introduction

A growing awareness of security concerns in automated physical processes has
increased interest in our ability to quantify the impact of would-be
attackers. Work along these lines imposes a detector to raise alerts when
sensor measurements do not fall in line with model-based predictions, thus
constraining what an attacker can do without being discovered [8, 7, 5]. The
tuning of such detectors is a balancing act between increasing the sensitivity
to attacks while reducing the number of false alarms (alerts raised during
normal operation) [9]. The next clear step in this direction of research is to
then minimize this impact through careful control system design. Initially
work used the distance (norm) that an attacker could drive the state as a
proxy for impact [9, 13], however, ultimately the reachable set provides
important information about which components of a system are effected more
than others and can inform whether the attacked state might reach dangerous
regions of state space [10]. The quantification, or analysis, studies have
often used ellipsoidal bounds on the actual attack-induced reachable set
achieved either through iterative methods [8] or the satisfaction of linear
matrix inequalities (LMIs) [5, 10]. The latter extends gracefully to
optimization to address the design question, although doing so often requires
re-linearizing the inequalities with respect to the new design variables.

In this paper, we leverage past work on quantification of ellipsoidal bounds
on the attack-induced reachable set to design the observer and controller gain
matrices to minimize the ellipsoidal bound (and thus the reachable set) when a
linear time invariant system is controlled with estimate-based feedback. As
has been pointed out (see, e.g., [11, 10]), it is important to pair a security
minimization with a constraint on closed-loop performance, otherwise a trivial
solution exists to disconnect the feedback loop and thus cut off the effect of
the attack on the system state. Here we specify an output covariance
constrained (OCC) $\|H\|_{2}$ performance [12], which takes into account the
covariance of the noise, unlike other distribution-agnostic approaches to
robust control [3, 4]. Because we are able to know the distribution of the
noises in our system, using them in the selection of optimal gains allows us
to exploit this known structure for a tighter and more tailored result. The
problem of covariance-based $\|H\|_{2}$ design is, however, in general a non-
convex problem. To make it compatible with our LMI framework, along the way we
develop a novel convexification of the OCC $\|H\|_{2}$ design problem, by
dividing the problem into two parts: a convex optimization and a generalized
algebraic Ricatti equation.

In [11] we presented the notion of observer gain design to minimize the
ellipsoidal bound on the estimation error. Here we complete the design problem
by extending the framework to quantify and minimize the reachable set of
attack-induced states, importantly integrating the design of the controller
gain matrix as well. The more recent work [10] provides a general and
expansive framework for using an LMI approach to solve the gain design problem
for security. Our work here is distinguished by (a) using estimate feedback as
opposed to a dynamic controller, (b) using a covariance-informed $\|H\|_{2}$
performance metric as opposed to a distributionally insensitive version, and
(c) introducing a magnification factor to scale the shape matrices associated
with the $\|H\|_{2}$ and reachable set decision variables. While the dynamic
controller - characterized by a linear time invariant dynamical system with
matrices $A^{c}$, $B^{c}$, $C^{c}$, and $D^{c}$ \- used in [10] is a more
general approach, the nonlinearities caused by the estimate-based feedback
adds additional complications to the steps to linearize the constraints. In
particular, if we define $A^{c}$, $B^{c}$, $C^{c}$, and $D^{c}$ to encode
estimate-based feedback, the resulting inequalities in [10] are nonlinear and
non-convex. The magnification factor we introduce here is effectively assumed
to be unit value in [10]; in our numerical studies the value of this parameter
tends to be large, indicating that including this factor greatly improves the
quality of the optimization solution. A conference paper [6] presented the
results for the iterative approach we take here, however, due to space the
full detail of the proofs were not included. The fully convexified approach is
not present in that paper.

##  2 Background

We consider a discrete-time linear time invariant (LTI) system of the form

$\displaystyle x_{k+1}$ $\displaystyle=Fx_{k}+Gu_{k}+\nu_{k},$ (1)
$\displaystyle y_{k}$ $\displaystyle=Cx_{k}+\eta_{k},$ (2)

in which the state $x_{k}\in\mathbb{R}^{n}$, $k\in\mathbb{N}$, evolves due to
the state update provided by the state matrix $F\in\mathbb{R}^{n\times n}$,
the control input $u_{k}\in\mathbb{R}^{m}$ shaped by the input matrix
$G\in\mathbb{R}^{n\times m}$, and the i.i.d Gaussian system noise
$\nu_{k}\sim\mathcal{N}(0,R_{1})$, $R_{1}\in\mathbb{V}^{n\times n}$
($\mathbb{V}$ is the set of positive definite matrices). The output
$y_{k}\in\mathbb{R}^{p}$ aggregates a linear combination of the states, given
by the observation matrix $C\in\mathbb{R}^{p\times n}$, and the i.i.d Gaussian
measurement noise $\eta_{k}\sim\mathcal{N}(0,R_{2})$,
$R_{2}\in\mathbb{V}^{p\times p}$. For the simplicity of the exposition, we
have considered Gaussian noises, however, the approach we present is
applicable for general noise distributions. In addition we assume that $F$ is
stable, the pair $(F,C)$ is detectable, $(F,G)$ is stabilizable and system and
measurement noises are mutually independent.

In this work, we consider the scenario that the actual measurement $y_{k}$ can
be corrupted by an attack, $\delta_{k}\in\mathbb{R}^{p}$. The attack is
injected at some point between the measurement and reception of the output by
the controller,

$\bar{y}_{k}=y_{k}+\delta_{k}=Cx_{k}+\eta_{k}+\delta_{k}.$ (3)

If the attacker has access to the measurements, then it is possible for the
attack $\delta_{k}$ to cancel some or all of the original measurement $y_{k}$
\- so an additive attack can achieve arbitrary control over the “effective”
output of the system.

Because our system is stochastic, we require an estimator to produce a
prediction of the system behavior

$\hat{x}_{k+1}=F\hat{x}_{k}+Gu_{k}+L(\bar{y}_{k}-C\hat{x}_{k}),$ (4)

where $\hat{x}_{k}\in\mathbb{R}^{n}$ is the estimated state and the observer
gain $L$ is designed to force the estimate to track the system states.

We consider observer-based feedback controllers

$u_{k}=K\hat{x}_{k},$ (5)

where $K\in\mathbb{R}^{m\times n}$ is the controller gain matrix. Next, we
define the residual sequence

$r_{k}=\bar{y}_{k}-C\hat{x}_{k},$ (6)

as the difference between what we actually receive ($\bar{y}_{k}$) and expect
to receive ($C\hat{x}_{k}$), which evolves according to

$\displaystyle x_{k+1}=(F+GK)x_{k}-GKe_{k}+\nu_{k}$ (7) $\displaystyle
e_{k+1}=\big{(}F-LC\big{)}e_{k}-L\eta_{k}+\nu_{k}-L\delta_{k},$ $\displaystyle
r_{k}=Ce_{k}+\eta_{k}+\delta_{k},$

where $e_{k}=x_{k}-\hat{x}_{k}$ is the estimation error. In the absence of
attacks (i.e., $\delta_{k}=0$), we can show that the steady-state distribution
of $r_{k}$ is Gaussian with covariance,

$\displaystyle\Sigma$
$\displaystyle=\mathbf{E}[r_{k}r_{k}^{T}]=C\mathbf{E}[e_{k}e_{k}^{T}]C^{T}+\mathbf{E}[\eta_{k}\eta_{k}^{T}],$
(8) $\displaystyle=C\mathbf{P}_{e}C^{T}+R_{2},$

where the steady state covariance of the estimation error
$\mathbf{P}_{e}=\lim_{k\rightarrow\infty}P_{k}=\lim_{k\rightarrow}\mathbf{E}[e_{k}e_{k}^{T}]$
is the solution of

$\displaystyle\mathbf{P}_{e}=(F-LC)\mathbf{P}_{e}(F-LC)^{T}+LR_{2}L^{T}+R_{1}.$
(9)

In this work, we consider the chi-squared detector, although similar analysis
can be done with other detector choices [9, 13]. The chi-squared detector
constructs a quadratic distance measure $z_{k}$ to be sensitive to changes in
the variance of the distribution as well as the expected value,

$z_{k}=r_{k}^{T}\Sigma^{-1}r_{k}.$ (10)

The chi-squared detector generates alarms when the distance measure exceeds a
threshold $\alpha\in\mathbb{R}_{>0}$

$\left\\{\begin{aligned} z_{k}\leq\alpha&\quad\longrightarrow\quad\text{no
alarm},\\\ z_{k}>\alpha&\quad\longrightarrow\quad\text{alarm:
}k^{\prime}=k,\end{aligned}\right.$ (11)

such that alarm time(s) $k^{\prime}$ are produced. The $\Sigma^{-1}$ factor in
the definition of $z_{k}$ re-scales the distribution ($\mathbf{E}[z_{k}]=p$,
$\mathbf{E}[z_{k}z_{k}^{T}]=2p$) so that the threshold $\alpha$ can be
designed independent of the specific statistics (mean and covariance) of the
noises $\nu_{k}$ and $\eta_{k}$; instead, it can be selected simply based on
the number of sensors, $p$ [9].

###  2.1 Definition of Attack

Detectors are designed to identify anomalies in system behavior. If an
attacker aims to remain undetected, the choice of detector and its parameters
limit what the attacker is able to accomplish. The type of attacks we consider
here require strong knowledge of and access to system dynamics, statistics of
the noises, current estimate ($\hat{x}_{k}$), and the detector configuration.
The goal of this powerful stealthy attack is to construct the worst case
scenario to aid the design of more robust systems.

Zero-alarm attacks employ attack sequences that maintain the distance measure
at or below the threshold of detection, i.e., $z_{k}\leq\alpha$. Hence, these
attacks generate no alarms. To satisfy this condition we define the attack as

$\delta_{k}=\phi_{k}-({y}_{k}-C\hat{x}_{k})=-Ce_{k}-\eta_{k}+\phi_{k},$ (12)

where $\phi_{k}\in\mathbb{R}^{p}$ is any vector such that
$\phi_{k}^{T}\Sigma^{-1}\phi_{k}\leq\alpha$ (recall the attacker has access to
the sensor, $y_{k}$, and knowledge of the estimator, $\hat{x}_{k}$). Based on
this attack strategy,

$\displaystyle z_{k}$ $\displaystyle=r_{k}^{T}\Sigma^{-1}r_{k},$
$\displaystyle=(Ce_{k}+\eta_{k}+\delta_{k})^{T}\Sigma^{-1}(Ce_{k}+\eta_{k}+\delta_{k}),$
$\displaystyle=\phi_{k}^{T}\Sigma^{-1}\phi_{k}\leq\alpha.$ (13)

Thus $z_{k}\leq\alpha$ and no alarms are raised.

###  2.2 Reachable Set

Under a stealthy zero-alarm attack (12), the attacked system dynamics become

$\displaystyle x_{k+1}$ $\displaystyle=Fx_{k}+GK\hat{x}_{k}+\nu_{k},$ (14)
$\displaystyle\hat{x}_{k+1}$ $\displaystyle=LCx_{k}+(F+GK-
LC)\hat{x}_{k}-LCe_{k}+L\phi_{k},$ $\displaystyle e_{k+1}$
$\displaystyle=Fe_{k}-L\phi_{k}+\nu_{k}.$

We stack these into a combined state
$\xi_{k}=\left[x_{k}^{T},\,\hat{x}_{k}^{T},\,e_{k}^{T}\right]^{T}$ and
combined input $\mu_{k}=\left[\nu_{k}^{T},\,\phi_{k}^{T}\right]^{T}$,

$\xi_{k+1}=A\xi_{k}+B\mu_{k},$ (15)

with

${A=\begin{bmatrix}F&GK&0\\\ LC&F+GK-LC&-LC\\\ 0&0&F\end{bmatrix},\
B=\begin{bmatrix}I&0\\\ 0&L\\\ I&-L\end{bmatrix}.}$ (16)

######  Remark 1.

The choice of including $x_{k}$, $\hat{x}_{k}$, and $e_{k}$ seems redundant at
this point since $e_{k}=x_{k}-\hat{x}_{k}$, however, this choice is crucial as
we layer additional constraints into the design optimization.

Throughout the rest of the paper we will use a selection matrix
$E_{x}=\left[I_{n},\,0_{n\times n},\,0_{n\times n}\right]$ to pull out
quantities relevant to the state $x_{k}=E_{x}\xi_{k}$.

The reachable set of attack-induced states is then,

$\mathcal{R}=\left\\{x_{k}=E_{x}\xi_{k}\ \left|\ \begin{aligned}
&\xi_{k+1}=A\xi_{k}+B\mu_{k},\\\ &\xi_{1}={0},\
\phi_{k}^{T}\Sigma^{-1}\phi_{k}\leq\alpha,\\\
&\nu_{k}^{T}R_{1}^{-1}\nu_{k}\leq\bar{\nu},\ \forall
k\in\mathbb{N}\end{aligned}\right.\right\\},$ (17)

where the ellipsoidal bound on the attack $\phi_{k}$ is imposed by the
attacker’s desire to remain stealthy (13), and the ellipsoidal bound on the
noise is created by truncating the Gaussian system noise to a desired
probability, i.e.,
$\text{Pr}[\nu_{k}^{T}R_{1}^{-1}\nu_{k}\leq\bar{\nu}]=p_{\nu}$, where
$p_{\nu}$ is some desired (typically high) probability. In principle, the
noise has unbounded support, and hence the reachable set is unbounded. To
ensure bounded reachable sets, we apply this truncation at the desired
confidence level.

##  3 LMI Approach to Design $K$ and $L$

In the first section, we reframe an existing result more concisely, which
identifies a minimal outer ellipsoidal bound on the set of states reachable by
a stealthy (zero-alarm) attacker. We then move to consider minimizing this set
further through the design of the feedback and estimator gains $K$ and $L$. As
has been discussed in previous studies, a trivial solution exists to this
design problem - to make either $GK=0$ or $L=0$. Doing so cuts the feedback
loop and guarantees that corrupted measurements do not impact the system
state. Simultaneously, this destroys the purpose - more specifically the
performance - of the feedback loop. While many performance metrics could be
used, in Section 3.2 we impose a $\|H\|_{2}$ constraint to avoid these trivial
solutions. Unlike prior work where the performance criteria ignored the
distribution of the noise, this $\|H\|_{2}$ constraint is specific to the
covariance of the noise, thereby allowing our design optimization to leverage
this important knowledge. This output covariance constrained (OCC) $\|H\|_{2}$
constraint is non-convex; to our knowledge, this paper offers the first
convexification of the OCC $\|H\|_{2}$ criteria into an LMI framework.

###  3.1 Bounding Ellipsoid LMI (given $K$ and $L$)

Before we move on to the synthesis problem of designing the gain matrices, we
first provide a solution to the analysis problem of finding a tight outer
ellipsoidal bound of the reachable set given $K$ and $L$, when the system is
driven by the system noise and attack. A similar analysis result appears in
[5], however, there the problem is split into two optimizations - one to find
a bound on the estimation error reachable set, the result of which is used in
the second optimization to bound the state reachable set. Here, in Lemma 2, we
solve these simultaneously through the stacked states $\xi_{k}$ and inputs
$\mu_{k}$. The following lemma provides a bound on a Lyapunov-inspired
function given an elliptically bounded input.

######  Lemma 1.

[10] Let $V_{k}$ be a positive definite function with $V_{1}=0$ and
$\mu_{ik}^{T}W_{i}\mu_{ik}\leq 1$, $i=1,\dots,N$, where $W_{i}$ is positive
definite. If there exists a constant $a\in(0,1)$ and $a_{i}\in(0,1)$ such that
$\sum_{i=1}^{N}a_{i}\geq a$ and

$V_{k+1}-aV_{k}-\sum_{i=1}^{N}(1-a_{i})\mu_{ik}^{T}W_{i}\mu_{ik}\leq 0,$ (18)

then $V_{k}\leq\frac{N-a}{1-a}$.

When we select the positive definite $V_{k}$ to be a quadratic function of the
state, the result above provides an outer ellipsoidal bound on the reachable
states. We will use the notation $\mathcal{E}(\mathcal{Q})=\\{x\ |\
x^{T}\mathcal{Q}^{-1}x\leq 1$}.

######  Lemma 2.

Given the stacked system matrices $A$ and $B$ in (16), gain matrices $K$, $L$,
detector threshold $\alpha$ with steady state residual covariance $\Sigma$,
system noise truncation threshold $\bar{\nu}$ with covariance $R_{1}$, if
there exists constants $a\in[0,1)$, the solution of

$\left\\{\begin{aligned}
\min_{a_{1},a_{2},\mathcal{Q}}&\mathbf{tr}(E_{x}^{T}\mathcal{\mathcal{Q}}E_{x})\\\
\text{s.t.}&\ \ 0\leq a_{1},a_{2}<1,\quad a_{1}+a_{2}\geq a,\\\ &\
\,\begin{bmatrix}a\mathcal{Q}&\mathcal{Q}{A}^{T}&0\\\
A\mathcal{Q}&\mathcal{Q}&B\\\ 0&{B}^{T}&\frac{1-a}{2-a}W\end{bmatrix}\geq
0,\end{aligned}\right.$ (19)

provides the shape matrix $\mathcal{Q}$ of the ellipsoidal bound on the
reachable set of states, i.e.,
$\mathcal{R}\subseteq\mathcal{E}(E_{x}^{T}\mathcal{Q}E_{x})$, where

$W=\begin{bmatrix}\frac{1-a_{1}}{\bar{\nu}}R_{1}^{-1}&0\\\
0&(1-a_{2})\Pi\end{bmatrix},\ \ \Pi=\frac{\Sigma^{-1}}{\alpha}.$ (20)

Proof: The stacked dynamics (15) is driven by two inputs which are both
ellipsoidally bounded. Letting $W_{1}=R_{1}^{-1}$ and $W_{2}=\Sigma^{-1}$,
(18) becomes

$V_{k+1}-aV_{k}-\frac{1-b}{\alpha}\phi_{k}^{T}\Sigma^{-1}\phi_{k}-\frac{1-a_{1}}{\bar{\nu}}\nu_{k}^{T}R_{1}^{-1}\nu_{k}\leq
0.$ (21)

Substituting the choice
$V_{k}=\xi_{k}^{T}\big{(}\frac{2-a}{1-a}\mathcal{P}\big{)}\xi_{k}\leq\frac{2-a}{1-a}$,
$\mathcal{P}>0$, into this equation and expanding using the dynamics (15)
results in the LMI,

$\mathcal{H}=\begin{bmatrix}a\mathcal{P}&{A}^{T}\mathcal{P}&0\\\
\mathcal{P}{A}&\mathcal{P}&\mathcal{P}{B}\\\
0&{B}^{T}\mathcal{P}&\frac{1-a}{2-a}W\end{bmatrix}\geq 0,$ (22)

where $\mathcal{P}$ is the inverse of the shape matrix of the ellipsoidal
bound for the $\xi$ reachable set ($\mathcal{P}^{-1}=\mathcal{Q}$), such that
the first block $E_{x}^{T}\mathcal{Q}E_{x}$ is the shape matrix of the
ellipsoidal bound of the reachable set of system states. To make this
ellipsoidal bound tight (as small as possible), the cost is selected to
minimize the trace of the shape matrix $E_{x}^{T}\mathcal{Q}E_{x}$. To use
$\mathcal{Q}$ as the variable of the optimization instead of $\mathcal{P}$ we
apply the transformation
$T=\text{diag}\left[\mathcal{Q},\,\mathcal{Q},\,I_{n}\right]$, to (22), i.e.,
$T^{T}\mathcal{H}T$, which results in the LMI in (19). $\blacksquare$

######  Remark 2.

For any convex shape (i.e., reachable set) there are an infinite number of
tight outer ellipsoidal bounds. These different ellipsoids can be visualized
as being tangent to the reachable set at different points. The minimum trace
objective minimizes the sum of the squared principal axes, which tends to
avoid solutions with, for example, low volume but one large principal axis.

######  Remark 3.

Note that the parameter $a$ is not a decision variable of the optimization in
(19). It appears nonlinearly (multiplying $\mathcal{Q}$). Since $a$ belongs to
a compact interval, the conventional choice is to solve (19) across a grid
search in $a$ and select the minimal, feasible solution.

###  3.2 Output Covariance Constrained $\|H\|_{2}$ Constraint

The introduction of this section and past related work has identified that
trivial solutions exist for the synthesis problem unless a performance
criteria is imposed in the optimization [11, 10]. One of the distinguishing
features of this work is that we consider an output covariance constrained
$\|H\|_{2}$ constraint, which involves the covariances of the system and
sensor noises. The challenge, tackled in the next subsection, is to convexify
and linearlize this inherently nonlinear constraint. Most optimizations in the
literature either use a distributionally robust constraint that is already
convex [10, 3] or solve the OCC $\|H\|_{2}$ using iterative algorithms [15,
14]. To specify the performance, Robust Control, in general, studies the gain
observed in the signal $h_{k}=H_{1}x_{k}+H_{2}\eta_{k}+H_{3}\nu_{k}$. Here,
for the system without attack, we consider the system driven by system and
measurement noise and enforce an $\|H\|_{2}$ constraint between the output
$h_{k}=y_{k}$ and excitation
$\omega_{k}=\left[\nu_{k}^{T},\,\eta_{k}^{T}\right]^{T}$, making $H_{1}=C$,
$H_{2}=I_{p\times p}$, and $H_{3}=0_{p\times n}$.

When there is no attack the system evolves according to

$\displaystyle x_{k+1}$ $\displaystyle=Fx_{k}+GK\hat{x}_{k}+\nu_{k},$ (23)
$\displaystyle\hat{x}_{k+1}$ $\displaystyle=LCx_{k}+(F+GK-
LC)\hat{x}_{k}+L\eta_{k},$ (24) $\displaystyle y_{k}$
$\displaystyle=Cx_{k}+\eta_{k},$ (25)

which can be combined using the stacked state
$\zeta_{k}=E_{x\hat{x}}\xi_{k}=\left[x_{k}^{T},\,\hat{x}_{k}^{T}\right]^{T}$,

$\zeta_{k+1}=\hat{A}\zeta_{k}+\hat{B}\omega_{k},$ (26)

with $E_{x\hat{x}}=\left[I_{2n},\,0_{2n\times n}\right]$ making
$\hat{A}=E_{x\hat{x}}AE_{x\hat{x}}^{T}$ and
$\hat{B}=E_{x\hat{x}}BE_{x\hat{x}}^{T}$.

######  Remark 4.

It is here that we can start to appreciate the value of the seemingly
redundant definition of $\xi_{k}$ (see Remark 1). By doing so, the state
matrix without attack $\hat{A}$ can be expressed as a sub-block of the state
matrix under attack $A$. Establishing this parallel structure is key towards
being able to integrate the $\|H\|_{2}$ constraint (without attack) with the
reachable set calculation (under attack).

The OCC $||H||_{2}$ criteria specifies the gain from the noise to the output
should be less than a desired value $\bar{\gamma}$,

$\lim_{N\to\infty}\sqrt{\frac{\frac{1}{N}\sum_{k=1}^{N}y_{k}^{T}y_{k}}{\frac{1}{N}\sum_{k=1}^{N}\omega_{k}^{T}\omega_{k}}}=\sqrt{\frac{\mathbf{E}[y_{k}^{T}y_{k}]}{\mathbf{E}[\omega_{k}^{T}\omega_{k}]}}\leq\bar{\gamma}.$
(27)

######  Lemma 3.

Given the dynamics in (26), the OCC $\|H\|_{2}$ constraint in (27) is
satisfied if the steady state covariance

$\mathbf{P}=\begin{bmatrix}\mathbf{P}_{x}&\mathbf{P}_{x\hat{x}}\\\
\mathbf{P}_{x\hat{x}}^{T}&\mathbf{P}_{\hat{x}}\end{bmatrix}=\lim_{k\to\infty}\mathbf{P}_{k}=\lim_{k\to\infty}\mathbf{E}[\zeta_{k}\zeta_{k}^{T}],$
(28)

satisfies the Lyapunov equation

$\mathbf{P}=\hat{A}\mathbf{P}\hat{A}^{T}+\hat{R},\qquad\mathbf{P}\geq 0,$ (29)

and the following convex inequality holds,

$\displaystyle\mathcal{C}_{h}=\mathbf{tr}\big{(}\hat{E}_{x}^{T}C^{T}C\hat{E}_{x}$
$\displaystyle\mathbf{P}\big{)}+\mathbf{tr}(R_{2})$ (30)
$\displaystyle-\bar{\gamma}^{2}\big{(}\mathbf{tr}(R_{1})+\mathbf{tr}(R_{2})\big{)}\leq
0,$

where $\hat{E}_{x}=\left[I_{n},\,0_{n\times n}\right]$.

Proof: From (25) and the definition of $\omega_{k}$, we can calculate the
quadratic terms in (27),

$\displaystyle y_{k}^{T}y_{k}$
$\displaystyle=x_{k}^{T}C^{T}Cx_{k}+2x_{k}^{T}C^{T}\eta_{k}+\eta_{k}^{T}\eta_{k},$
(31) $\displaystyle\omega_{k}^{T}\omega_{k}$
$\displaystyle=\nu_{k}^{T}\nu_{k}+\eta_{k}^{T}\eta_{k}.$ (32)

Taking the expectation ($x_{k}$ and $\eta_{k}$ are independent),

$\displaystyle\mathbf{E}[y_{k}^{T}y_{k}]$
$\displaystyle=\mathbf{E}[x_{k}^{T}C^{T}Cx_{k}]+\mathbf{E}[\eta_{k}^{T}\eta_{k}]$
(33) $\displaystyle=\mathbf{tr}\big{(}C^{T}C\
\mathbf{E}[x_{k}x_{k}^{T}]\big{)}+\mathbf{tr}(R_{2}),$

and similarly,

$\mathbf{E}[\omega_{k}^{T}\omega_{k}]=\mathbf{E}[\nu_{k}^{T}\nu_{k}]+\mathbf{E}[\eta_{k}^{T}\eta_{k}]=\mathbf{tr}(R_{1})+\mathbf{tr}(R_{2}).$
(34)

The unknown quantity is then the covariance of the state,
$\mathbf{E}[x_{k}x_{k}^{T}]$, which is the first block of the stacked state
$\zeta_{k}$ covariance $\mathbf{P}_{k}=\mathbf{E}[\zeta_{k}\zeta_{k}^{T}]$.
This covariance follows the update, evaluating
$\mathbf{E}[\zeta_{k+1}\zeta_{k+1}^{T}]$ with (26),

$\mathbf{P}_{k+1}=\hat{A}\mathbf{P}_{k}\hat{A}^{T}+\hat{R},\qquad\hat{R}=\begin{bmatrix}R_{1}&0\\\
0&LR_{2}L^{T}\end{bmatrix}.$ (35)

Because the matrix $\hat{A}$ is stable the covariance converges to a steady
value $\lim_{k\to\infty}\mathbf{P}_{k}=\mathbf{P}$ which satisfies the
Lyapunov equation (29). Combining this $\mathbf{P}$ with (33)-(34) the
$\|H\|_{2}$ constraint becomes

$\gamma=\sqrt{\frac{\mathbf{tr}\big{(}C\mathbf{P}_{x}C^{T}\big{)}+\mathbf{tr}(R_{2})}{\mathbf{tr}(R_{2})+\mathbf{tr}(R_{1})}}\leq\bar{\gamma},$
(36)

where $\gamma$ is the actual performance and $\bar{\gamma}$ is the worst
allowable performance. Rearranging this leads to the condition (30).
$\blacksquare$

In order to use this $\|H\|_{2}$ constraint in a convex optimization we need
to linearize the Lyapunov equation constraint. We state this result as part of
a complete convex optimization problem to design the gains $K$ and $L$ to
achieve the optimal (smallest) $\|H\|_{2}$ gain.

######  Theorem 1.

Given the dynamics in (26), the smallest output covariance constrained
$\|H\|_{2}$ gain defined by (27) is

$\gamma^{*}=\sqrt{\frac{\mathbf{tr}\big{(}C\mathbf{P}_{x}^{*}C^{T}\big{)}+\mathbf{tr}(R_{2})}{\mathbf{tr}(R_{2})+\mathbf{tr}(R_{1})}},$
(37)

where $\mathbf{P}_{x}^{*}$ is the solution of

$\left\\{\begin{aligned}
\min_{\mathbf{P}_{x},\mathbf{Q}_{1},X,Y,Z}&\mathbf{tr}(C\mathbf{P}_{x}C^{T})\\\
\text{s.t.}\quad&\mathcal{C}_{L}\geq 0,\end{aligned}\right.$ (38)

with

$\mathcal{C}_{L}=\leavevmode\resizebox{352.31624pt}{}{$\begin{bmatrix}\mathbf{Q}_{1}&I&\mathbf{Q}_{1}F+XC&Z&\mathbf{Q}_{1}R_{1}&XR_{2}\\\
*&\mathbf{P}_{x}&F&F\mathbf{P}_{x}+GY&R_{1}&0\\\ *&*&\mathbf{Q}_{1}&I&0&0\\\
*&*&*&\mathbf{P}_{x}&0&0\\\ *&*&*&*&R_{1}&0\\\
*&*&*&*&*&R_{2}\end{bmatrix}$}.$ (39)

There exists at most $2n\choose n$ distinct real-valued, control gains
$L=\mathbf{Q}_{12}^{-1}X$, $K=Y\mathbf{P}_{x\hat{x}}^{-T}$ satisfying
$\gamma=\gamma^{*}$, where
$\mathbf{P}_{x\hat{x}}=(I-\mathbf{P}_{x}\mathbf{Q}_{1})\mathbf{Q}_{12}^{-T}$
and $\mathbf{Q}_{12}$ is the solution of following generalized algebraic
Ricatti equation,

$\mathbf{Q}_{12}\Gamma_{1}\mathbf{Q}_{12}+\mathbf{Q}_{12}\Gamma_{2}+\Gamma_{3}\mathbf{Q}_{12}+\Gamma_{4}=0,$
(40)

with known matrices

$\displaystyle\Gamma_{1}=GY(I-\mathbf{Q}_{1}\mathbf{P}_{x})^{-1},$ (41)
$\displaystyle\Gamma_{2}=F,$
$\displaystyle\Gamma_{3}=(\mathbf{Q}_{1}GY+XC\mathbf{P}_{x}+\mathbf{Q}_{1}F\mathbf{P}_{x}-Z)(I-\mathbf{Q}_{1}\mathbf{P}_{x})^{-1},$
$\displaystyle\Gamma_{4}=-XC.$

Proof: The formula for the optimal gain $\gamma^{*}$, (37), comes naturally
from the $\|H\|_{2}$ bound derived in (36). Since all other terms are
constant, minimizing $\mathbf{tr}(C\mathbf{P}_{x}C^{T})$ is equivalent to
minimizing the gain. This covariance is constrained by the Lyapunov equation
in (29). Here, which is a standard technique for incorporating Lyapunov
equations into convex optimizations, we replace this equality constraint with
the very similar inequality,

$\mathbf{P}-\hat{A}\mathbf{P}\hat{A}^{T}-\hat{R}\geq 0,\qquad\mathbf{P}\geq
0.$ (42)

We can now combine these two inequality constraints into one using the Schur
complement [1],

$\mathcal{C}=\begin{bmatrix}\mathbf{P}-\hat{R}&\hat{A}\mathbf{P}\\\
\mathbf{P}\hat{A}^{T}&\mathbf{P}\end{bmatrix}\geq 0.$ (43)

This relaxation is justified because the objective function
$\mathbf{tr}(C\mathbf{P}_{x}C^{T})$ minimizes the decision variable
$\mathbf{P}$ and drives the optimization to the bound of the inequality, which
would yield equality - hence driving the relaxed form (42) to the equality
(29) (see Appendix A).

We use the following transformation to linearize $\mathcal{C}$,

$\mathcal{C}_{L}=\begin{bmatrix}T_{1}&\\\
&T_{1}\end{bmatrix}^{T}\mathcal{C}\begin{bmatrix}T_{1}&\\\
&T_{1}\end{bmatrix}=\begin{bmatrix}\mathbf{P}_{L}-R_{L}&\hat{A}_{L}\\\
\hat{A}_{L}^{T}&\mathbf{P}_{L}\end{bmatrix},$ (44)

with

$T_{1}=\begin{bmatrix}\mathbf{Q}_{1}&I\\\
\mathbf{Q}_{12}^{T}&0\end{bmatrix},\qquad\mathbf{P}^{-1}=\mathbf{Q}=\begin{bmatrix}\mathbf{Q}_{1}&\mathbf{Q}_{12}\\\
\mathbf{Q}_{12}^{T}&\mathbf{Q}_{2}\end{bmatrix},$ (45)

and

$\displaystyle\mathbf{P}_{L}=T_{1}^{T}\mathbf{P}T_{1}=\begin{bmatrix}\mathbf{Q}_{1}&I\\\
I&\mathbf{P}_{x}\end{bmatrix},$ (46) $\displaystyle
R_{L}=T_{1}^{T}RT_{1}=\begin{bmatrix}\mathbf{Q}_{1}R_{1}\mathbf{Q}_{1}\ +\
\mathbf{Q}_{12}LR_{2}L^{T}\mathbf{Q}_{12}^{T}&\mathbf{Q}_{1}R_{1}\\\
R_{1}\mathbf{Q}_{1}&R_{1}\end{bmatrix},$
$\displaystyle\hat{A}_{L}=T_{1}^{T}\hat{A}\mathbf{P}T_{1}=\begin{bmatrix}\mathbf{Q}_{1}F+XC&Z\\\
F&F\mathbf{P}_{x}+GY\end{bmatrix},$

where $X$, $Y$, and $Z$ are defined as

$\displaystyle X$ $\displaystyle=\mathbf{Q}_{12}L,$ (47) $\displaystyle Y$
$\displaystyle=K\mathbf{P}_{x\hat{x}}^{T},$ (48) $\displaystyle Z$
$\displaystyle=\mathbf{Q}_{1}F\mathbf{P}_{x}+XC\mathbf{P}_{x}+\mathbf{Q}_{1}GY+\mathbf{Q}_{12}F\mathbf{P}_{x\hat{x}}^{T}$
(49) $\displaystyle\qquad\qquad+\mathbf{Q}_{12}GY-
XC\mathbf{P}_{x\hat{x}}^{T}.$

The term $\mathbf{P}_{L}-R_{L}$ can be linearized by applying a Schur
complement to recover $\mathcal{C}_{L}$ in (39). This transformation changes
the set of decision variables from ($\mathbf{P}_{x}$, $\mathbf{P}_{\hat{x}}$,
$\mathbf{P}_{x\hat{x}}$, $L$, $K$) to ($\mathbf{P}_{x}$, $\mathbf{Q}_{1}$,
$X$, $Z$, $Y$). The solution in these new decision variables is then used to
calculate $\mathbf{P}_{x\hat{x}}$ and $\mathbf{Q}_{12}$ using (49) and the
identity

$\mathbf{P}_{x}\mathbf{Q}_{1}+\mathbf{P}_{x\hat{x}}\mathbf{Q}_{12}^{T}=I,$
(50)

which comes from the first block of the definition $\mathbf{P}\mathbf{Q}=I$.
The definition of $Z$ in (49) and (50) combine to form the general algebraic
Ricatti equation (40), which, in general, has $2n\choose n$ different answers.
Finally, the gain matrices can be found by, $L=\mathbf{Q}_{12}^{-1}X$ and
$K=Y\mathbf{P}_{x\hat{x}}^{-T}$. Note that the solutions to the original
Lyapunov equality (29) are a subset of the solutions of the relaxed inequality
(42); so Theorem 1 characterizes all OCC $\|H\|_{2}$ optimal solutions.
$\blacksquare$

At the other end of the performance spectrum, security increases (size of the
reachable set decreases) as either $L$ or $GK$ approach zero. This represents
worst-case performance (without closed loop control) and the corresponding
value of $\gamma=\gamma_{0}$ can be found using Lemma 4.

######  Lemma 4.

Given the system dynamics,in the absence of attack, the open loop (i.e., $L=0$
or $GK=0$) OCC $\|H\|_{2}$ gain $\gamma_{0}$ is given by

$\gamma_{0}=\sqrt{\frac{\mathrm{tr}\big{(}C\mathbf{P}_{x}C^{T}\big{)}+\mathrm{tr}(R_{2})}{\mathrm{tr}(R_{2})+\mathrm{tr}(R_{1})}},$

where the steady state covariance of state $\mathbf{P}_{x}$ is the solution of

$F\mathbf{P}_{x}F^{T}-\mathbf{P}_{x}+R_{1}=0.$

Proof. For both $L=0$ or $GK=0$, the state dynamics (23) are in open loop and
thus the evolution of the system is the same. When $L=0$, the state estimate
$\hat{x}_{k}$ converges to zero because the system is open loop stable and the
open loop state dynamic becomes

$x_{k+1}=Fx_{k}+\nu_{k}.$ (51)

Similarly, when $GK=0$, equation (23) becomes (51) directly. With this state
equation, the steady state covariance of the state, $\mathbf{P}_{x}$, is given
by (4) and consequently the desired performance $\gamma_{0}$ should be the
same in both cases. $\blacksquare$

Therefore, now entering into the design process, for all choices of
$\bar{\gamma}\in[\gamma^{*},\ \infty]$, the solution for performance
$\gamma<\bar{\gamma}$ always lies within the trade-off interval,

$\gamma^{*}\leq\gamma\leq\gamma_{0}.$ (52)

###  3.3 Bounding ellipsoid LMI (designing $K$ and $L$)

The goal of this paper is to construct an optimization to design $K$ and $L$
such that the impact of an attacker on the reachable states is minimized.
However, when $K$ and $L$ are considered variables of the Lemma 2
optimization, (22), and therefore, (19) contains nonlinear terms. In the
sections that follow, we impose some structure on the solution so that we can
linearize the overall design problem. Each choice will be motivated
individually, but it is also the combined effect of the these structures taken
together that yield the final linear matrix inequality. The first three (out
of four) represent a choice of how to select an outer ellipsoidal bound. While
they do impose some structure on the solution, they are best seen as one
choice out of infinitely many equally good options, so they are not very
limiting.

######  Imposed Structure 1.

There are an infinite number of tight outer ellipsoidal bounds of the stacked
state $\xi$. Of these we select one that satisfies the following structure for
the inverse of the shape matrix,

$\mathcal{P}=\begin{bmatrix}\mathcal{P}_{1}&\mathcal{P}_{12}&\\\
\mathcal{P}_{12}^{T}&\mathcal{P}_{2}&\\\ &&\mathcal{P}_{3}\end{bmatrix},$ (53)

which assumes the independence of the ellipsoidal bound on the estimation
error $e_{k}$ from the ellipsoidal bound on the combined state $x_{k}$ and
estimate $\hat{x}_{k}$. This is inspired by a similar assumption made in [10].
This selection enables us to utilize the parallel dynamics with and without
attack (see Remarks 1 and 4) and linearize the original LMI with respect to
$K$ and $L$.

This selection also permits inverting each block separately, such that
$\mathcal{P}_{3}^{-1}=\mathcal{Q}_{e}$ and

$\begin{bmatrix}\mathcal{P}_{1}&\mathcal{P}_{12}\\\
\mathcal{P}_{12}^{T}&\mathcal{P}_{2}\end{bmatrix}^{-1}=\begin{bmatrix}\mathcal{Q}_{x}&\mathcal{Q}_{x\hat{x}}\\\
\mathcal{Q}_{x\hat{x}}^{T}&\mathcal{Q}_{\hat{x}}\end{bmatrix}.$ (54)

Consider the linearizing change of coordinates used in [4, 10],

$T_{2}=\begin{bmatrix}T_{3}&&\\\ &T_{3}&\\\ &&I_{n}\end{bmatrix},\qquad
T_{3}=\begin{bmatrix}\mathcal{Q}_{x}&I&0\\\ \mathcal{Q}_{x\hat{x}}^{T}&0&0\\\
0&0&I\end{bmatrix}.$ (55)

Although (22) is not entirely linearized with this transformation, due to the
presence of term $\Sigma$ which depends on $L$, we will introduce an iterative
approach later to avoid this nonlinearity. The LMI $\mathcal{H}$, (22),
becomes

$\displaystyle\mathcal{H}_{L}=T_{2}^{T}\mathcal{H}T_{2}=\begin{bmatrix}a\mathcal{P}_{L}&A_{L}^{T}&0\\\
A_{L}&\mathcal{P}_{L}&B_{L}\\\ 0&B_{L}^{T}&\frac{1-a}{2-a}W\end{bmatrix},$
(56)

where

$\displaystyle\mathcal{P}_{L}=T_{3}^{T}\mathcal{P}T_{3}=\begin{bmatrix}\mathcal{Q}_{x}&I&0\\\
I&\mathcal{P}_{1}&0\\\ 0&0&\mathcal{P}_{3}\end{bmatrix},$ (57) $\displaystyle
B_{L}=T_{3}^{T}\mathcal{P}B=\begin{bmatrix}I&0\\\ \mathcal{P}_{1}&Y_{1}\\\
\mathcal{P}_{3}&-\mathcal{P}_{3}L\end{bmatrix},$ $\displaystyle
A_{L}=T_{3}^{T}\mathcal{P}AT_{3}=\begin{bmatrix}F\mathcal{Q}_{x}+GX_{1}&F&0\\\
Z_{1}&\mathcal{P}_{1}F+Y_{1}C&-Y_{1}C\\\ 0&0&\mathcal{P}_{3}F\end{bmatrix},$
$\displaystyle Y_{1}=\mathcal{P}_{12}L,\quad
X_{1}=K\mathcal{Q}_{x\hat{x}}^{T},$ $\displaystyle
Z_{1}=\mathcal{P}_{1}F\mathcal{Q}_{x}+\mathcal{P}_{12}LC\mathcal{Q}_{x}+\mathcal{P}_{1}GK\mathcal{Q}_{x\hat{x}}^{T}+\mathcal{P}_{12}F\mathcal{Q}_{x\hat{x}}^{T}$
$\displaystyle\qquad+\mathcal{P}_{12}GK\mathcal{Q}_{x\hat{x}}^{T}-\mathcal{P}_{12}LC\mathcal{Q}_{x\hat{x}}^{T}.$

One of the useful features of this transformation is that
$\mathcal{Q}_{x}=E_{x}^{T}\mathcal{Q}E_{x}$, the quantify used in the
objective function of Lemma 2, appears as a variable of the LMI. This section
provides the linearization necessary to separate the gains $K$ and $L$ as
variables in Lemma 2 (and could then be used as the starting point if a
different performance criteria was used, as opposed to the $\|H\|_{2}$
constraint considered in this paper).

Figure 1: The reachable set (gray) is approximated by ellipsoids with shape
matrix $\sigma\mathbf{P}$. Without relaxation effects, the role of $\sigma$ is
to identify the tight approximation (green), where $\sigma=\sigma^{*}$. When
$\sigma<\sigma^{*}$ the optimization will be infeasible or doesn’t converge,
because the red ellipsoid cannot contain the reachable set. When
$\sigma>\sigma^{*}$ the black ellipsoid loosely contains the reachable set.
Due to our techniques to linearize the optimization, the reachable set
approximation will have some extra conservatism and we do not expect the
optimal ellipsoid to be exactly tangent.

###  3.4 Combining Performance and Security

In this work, we design the controller and estimator gains to minimize the
impact of attacks on the system state, which is measured by an outer
ellipsoidal bound on the reachable states when the system is driven by the
attack and system noise. As Remark 2 states, there are an infinite number of
potential outer bounding - and tight - ellipsoids. In order to combine the LMI
constraints from the reachable set and $\|H\|_{2}$ calculations, we make a
specific choice about the outer ellipsoidal bound we select.

######  Imposed Structure 2.

We select the shape matrix of the ellipsoidal bound of the states $x_{k}$ and
estimate $\hat{x}_{k}$ under attack $E_{x\hat{x}}^{T}\mathcal{Q}E_{x\hat{x}}$
\- see (54) - to have the same orientation as the covariance of the states and
estimate without attack ($\zeta_{k}$),

$\sigma\mathbf{P}=E_{x\hat{x}}^{T}\mathcal{Q}E_{x\hat{x}},$ (58)

where $\sigma$ is a scaling factor that becomes a new variable of the method
and is a function of gains $(L,\ K)$. Since $\mathbf{Q}=\mathbf{P}^{-1}$, this
sets up a common set of variables to link the $\|H\|_{2}$ (left) and
ellipsoidal bound (right) constraints,

$\begin{aligned} \sigma\begin{bmatrix}\mathbf{P}_{x}&\mathbf{P}_{x\hat{x}}\\\
\mathbf{P}_{x\hat{x}}^{T}&\mathbf{P}_{\hat{x}}\end{bmatrix}&=\begin{bmatrix}\mathcal{Q}_{x}&\mathcal{Q}_{x\hat{x}}\\\
\mathcal{Q}_{x\hat{x}}^{T}&\mathcal{Q}_{\hat{x}}\end{bmatrix},\\\
\underbrace{\frac{1}{\sigma}\begin{bmatrix}\mathbf{Q}_{1}&\mathbf{Q}_{12}\\\
\mathbf{Q}_{12}^{T}&\mathbf{Q}_{2}\end{bmatrix}}_{\|H\|_{2}}&=\underbrace{\begin{bmatrix}\mathcal{P}_{1}&\mathcal{P}_{12}\\\
\mathcal{P}_{12}^{T}&\mathcal{P}_{2}\end{bmatrix}.}_{\mathcal{E}(\mathcal{Q})}\end{aligned}\
$ (59)

The structure above allows us to replace variables in the ellipsoidal bound
optimization $\mathcal{Q}$ and $\mathcal{P}$ with quantities from the
performance criteria, $\mathbf{P}$ and $\mathbf{Q}$, respectively.

Based on (59) we can link the variables of the bounding ellipsoid LMI with the
$\|H\|_{2}$ constraint,

$\displaystyle X$ $\displaystyle=\sigma Y_{1}=\mathbf{Q}_{12}L,\quad
Y=\frac{X_{1}}{\sigma}=K\mathbf{P}_{x\hat{x}}^{T},$ (60) $\displaystyle Z$
$\displaystyle=Z_{1}=\mathbf{Q}_{1}F\mathbf{P}_{x}+XC\mathbf{P}_{x}+\mathbf{Q}_{1}GY+\mathbf{Q}_{12}F\mathbf{P}_{x\hat{x}}^{T}$
$\displaystyle\qquad\qquad+\mathbf{Q}_{12}GY-XC\mathbf{P}_{x\hat{x}}^{T}.$

Now we can rewrite $A_{L}$, $B_{L}$ , $\mathcal{P}_{L}$ based on
$\mathbf{P}_{x}$, $\mathbf{Q}_{1}$, $\mathcal{P}_{3}$, $X$, $Y$, $Z$,

$\displaystyle A_{L}=\begin{bmatrix}\sigma(F\mathbf{P}_{x}+GY)&F&0\\\
Z&\frac{1}{\sigma}(\mathbf{Q}_{1}F+XC)&-\frac{1}{\sigma}XC\\\
0&0&\mathcal{P}_{3}F\end{bmatrix},$ (61) $\displaystyle
B_{L}=\begin{bmatrix}I&0\\\
\frac{1}{\sigma}\mathbf{Q}_{1}&\frac{1}{\sigma}X\\\
\mathcal{P}_{3}&-\mathcal{P}_{3}L\end{bmatrix},\quad\mathcal{P}_{L}=\begin{bmatrix}\sigma\mathbf{P}_{x}&I&0\\\
I&\frac{1}{\sigma}\mathbf{Q}_{1}&0\\\ 0&0&\mathcal{P}_{3}\end{bmatrix}.$

Thus the choice in (58) has facilitated integrating these optimizations.

$L\leftarrow$ such that $\mathcal{C}_{L}\geq 0$ is feasible

$\sigma\leftarrow\infty$

while  _true_ do

if  _( 62) is infeasible _ then $\mathbf{break}$

$\mathbf{P}_{x},\mathbf{Q}_{1},X,Y,Z\leftarrow$ (62)

$\mathcal{K},\mathcal{L}\leftarrow$ real solutions of Ricatti eqn (40)

$(\tilde{K},\tilde{L})\leftarrow$ select pair
$(K,L)\in(\mathcal{K},\mathcal{L})$ by smallest

Lemma 2 objective value

if  _$\|L-\tilde{L}\|\leq\epsilon$_ then $\sigma\leftarrow\sigma-\varepsilon$

else $K,L\leftarrow\tilde{K},\tilde{L}$

end while

Algorithm 1 $K,L=$ Theorem 2$(F,G,C,R_{1},R_{2})$

######  Theorem 2.

Consider a LTI system (1) with maximum allowable output covariance constrained
$\|H\|_{2}$ gain $\bar{\gamma}$ (27), chi-squared detector threshold $\alpha$
(11) and zero-alarm stealthy attacker (13). Algorithm 1 returns
(approximately) optimal controller $K^{*}$ and observer $L^{*}$ gains to
minimize the reachable set of states possible by the attacker, while
maintaining an OCC $\|H\|_{2}$ gain no bigger than $\bar{\gamma}$. Algorithm 1
uses the Ricatti equation (40) to update $L$ based on the solution of the
combined convex optimization problem which has a solution if for some
$a\in[0,1)$,

$\displaystyle\left\\{\begin{aligned}
&\min_{{a_{1},a_{2},\mathbf{P}_{x},\mathbf{Q}_{1},}\atop{X,Y,Z,\mathcal{P}_{3}}}\mathbf{tr}(\sigma\mathbf{P}_{x})\\\
&\qquad\text{s.t.}\quad 0\leq a_{1},a_{2}<1,\quad a_{1}+a_{2}\geq a,\\\
&\qquad\qquad\mathcal{H}_{L}\geq 0,\\\ &\qquad\qquad\mathcal{C}_{h}\leq 0,\
\mathcal{C}_{L}\geq 0.\end{aligned}\right.$ (62)

Proof: In the past sections we have linearized the LMIs associated with the
ellipsoidal outer bound on the reachable set ($\mathcal{H}_{L}$) and with the
$\|H\|_{2}$ constraint ($\mathcal{C}_{h}$ and $\mathcal{C}_{L}$) and finally
made a structural connection between these two optimizations (59) to use a
common set of decision variables. Because the controller gain $K$ appears
within the decision variables of the optimization, we implicitly optimize $K$.

There are two remaining challenges to be addressed by this algorithm. First,
$L$ appears, as $K$ does, implicitly in the decision variables, but also
explicitly in the nonlinear term $\mathcal{P}_{3}L$ in matrix $B_{L}$ and in
the dependence of covariance $\Sigma$ in matrix $W$ on $L$ (8), both of which
are in $\mathcal{H}_{L}$. Second, the magnification factor $\sigma$ multiplies
most of the decision variables of the optimization. Thus neither $L$ nor
$\sigma$ can be taken as variable in the optimization. We solve this by
applying an iterative algorithm over both $L$ and $\sigma$ that leverages the
structure presented in Figure 1. If we select a large enough value for the
magnification parameter $\sigma$, any choice of $L$ easily satisfies
$\mathcal{H}_{L}$ by creating a very large ellipsoidal outer bound. We satisfy
the other constraints of the optimization (the $\|H\|_{2}$ constraints) by
selecting the initial value for $L$ as the $\|H\|_{2}$ optimal $L$ (from
Theorem 1), hence satisfying $\mathcal{C}_{h}$ and $\mathcal{C}_{L}$. For
these fixed values of $\sigma$, $L$, and $\Sigma$, the optimization (62) is
solved. The solution then provides a new value of $L$ \- solved from the
Ricatti equation in (40) - that minimizes the ellipsoidal bound while
satisfying all constraints. Using the same value for $\sigma$ but the updated
$L$ (and hence updated $\Sigma$), optimization (62) is again solved, yielding
another value for $L$. This iteration is repeated until $L$ has sufficiently
converged. Once convergence is achieved, $\sigma$ is reduced and the process
is repeated with the existing $L$ as the initial value for the $L$ iteration.
At some point, the magnification factor $\sigma$ will be too small for any
choice of $L$ to permit the ellipsoidal bound to contain the reachable set,
hence the optimization (62) will become infeasible. This is the stopping
condition for the algorithm.

Note that the convergence criteria for $L$ and the decrement amount for
$\sigma$ are selected by the user. In principle, these should be small, but
making them larger will allow the algorithm to require fewer iterations. The
number of iterations for $L$ to converge tends to be quite small (typically
2-3). In practice, the decrement of $\sigma$ can also be accomplished through
a bisection algorithm which looks for the smallest value of $\sigma$ that
makes the optimization (62) feasible.

Recall that, in general, solving the Ricatti equation (40) yields up to
$2n\choose n$ solutions. We assert that there will always exist at least one
real solution due to matrices $\mathbf{P}_{L}$, $\mathbf{P}$,
$\mathcal{C}_{L}$, and $\mathcal{C}$ being positive definite. To chose between
the real solutions, we choose the one that yields the smallest ellipsoidal
outer bound, as determined by Lemma 2 (which provides a tighter approximation
of the reachable set, since it does not require the additional linearization
steps).

Finally, the stability of the closed loop system is implicitly guaranteed if
$\mathcal{H}_{L}>0$ and $\mathcal{C}_{L}>0$ (see Appendix B). $\blacksquare$

######  Remark 5.

In practice, the quality of the solution is related to the initial guess for
$L$. For the cases where $\bar{\gamma}$ is close to $\gamma^{*}$ we use $L$
from Theorem 1 as the initial guess. Otherwise when $\bar{\gamma}$ is far from
$\gamma^{*}$ we solve the problem in multiple steps increasing $\bar{\gamma}$
gradually from $\gamma^{*}$ at each step solving Theorem 2 and using its
solution as the initial values for the next step where $\bar{\gamma}$ is
increased.

##  4 Resolving the Nonlinearity

The approach offered in Theorem 2 relies on an iterative scheme to avoid the
nonlinearities surrounding the observer gain matrix $L$. In each iteration,
the solution of the convex optimization in (62) is used to find $L$ using the
Ricatti equation in (40) and the corresponding residual covariance $\Sigma$
using the Lyapunov equation in (9). To eliminate the iterative approach, both
of these nonlinearities must be linearized and absorbed into the convex
optimization. A final nonlinearity is the existance of the $-\mathcal{P}_{3}L$
term within $\mathcal{H}_{L}$ \- specifically $B_{L}$ in (61). We accomplish
this linearization through the careful selection of additional structure,
which we will show only marginally reduces the quality of the solutions found.

######  Imposed Structure 3.

Here $\mathcal{Q}_{e}=\mathcal{P}_{3}^{-1}$ is the shape matrix of the
estimation error bounding ellipsoid. Again there are an infinite number of
choices for its orientation and, following the same approach taken in Imposed
Structure 2, here we choose the orientation that matches that of the
estimation error covariance in the absence of attack
$\mathcal{P}_{3}=\frac{1}{\sigma}\mathbf{P}_{e}^{-1}$.

######  Imposed Structure 4.

The final linearization assumption is to search for solutions on the manifold
that satisfies $\mathbf{P}_{\hat{x}}=\mathbf{P}_{x\hat{x}}$. This is the most
restrictive of the four structures that have been imposed to linearize the
problem, however, we will see it does not greatly impact the conservatism of
the optimization. One of the direct effects of this structure is that the
generalized algebraic Ricatti equation (40), one of the key nonlinearities, is
directly linearized. However, it is this final structure that ultimately
removes all three nonlinearities present in Theorem 2.

From Imposed Structure 4 and the identity $\mathbf{P}\mathbf{Q}=I$ \- see (45)
- we can show that $\mathbf{P}_{e}^{-1}=\mathbf{Q}_{1}=-\mathbf{Q}_{12}$ (see
Appendix C). The relationship between $\mathbf{Q}_{1}$ and $\mathbf{Q}_{12}$
is primarily what linearizes the Ricatti equation in (40). The relationship
between $\mathbf{P}_{e}$ and $\mathbf{Q}_{1}$, combined with Imposed structure
3, enables us to replace the $-\mathcal{P}_{3}L$ nonlinearity with
$\frac{1}{\sigma}X$. Therefore we can re-write matrices
$A_{L},B_{L},\mathcal{P}_{L}$ as,

$\displaystyle A_{L}=\begin{bmatrix}\sigma(F\mathbf{P}_{x}+GY)&F&0\\\
Z&\frac{1}{\sigma}(\mathbf{Q}_{1}F+XC)&-\frac{1}{\sigma}XC\\\
0&0&\frac{1}{\sigma}\mathbf{Q}_{1}F\end{bmatrix},$ (63) $\displaystyle
B_{L}=\begin{bmatrix}[r]I\ \ &0\ \ \\\
\frac{1}{\sigma}\mathbf{Q}_{1}&\frac{1}{\sigma}X\\\
\frac{1}{\sigma}\mathbf{Q}_{1}&\frac{1}{\sigma}X\end{bmatrix},\quad\mathcal{P}_{L}=\begin{bmatrix}\sigma\mathbf{P}_{x}&I&0\\\
I&\frac{1}{\sigma}\mathbf{Q}_{1}&0\\\
0&0&\frac{1}{\sigma}\mathbf{Q}_{1}\end{bmatrix}.$

Finally, the relationship between $\mathbf{P}_{e}$ and $\mathbf{Q}_{1}$ also
help to integrate the Lyapunov equation in (9) into the convex optimization.
Lemma 5 relates the matrix $\Pi=\Sigma^{-1}/\alpha$ with the existing decision
variables. This requires us, like the system and sensor noises, to truncate
the estimation error distribution at some confidence level, measured by
$\bar{e}$.

######  Lemma 5.

Consider the truncated Gaussian system noise and measurement noise with
$\nu_{k}^{T}R_{1}^{-1}\nu_{k}\leq\bar{\nu}$ and
$\eta_{k}^{T}R_{2}^{-1}\eta_{k}\leq\bar{\eta}$ and their corresponding
Gaussian estimation error and residual with
$e_{k}^{T}\mathbf{P}_{e}^{-1}e_{k}\leq\bar{e}$ and $r_{k}^{T}\Pi r_{k}\leq 1$.
The positive definite matrix $\Pi$, can be expressed as
$\Pi=\frac{\Sigma^{-1}}{\alpha}$ where $\Sigma$ comes from (8) and (9), if,

$\displaystyle\mathcal{X}_{L}=\begin{bmatrix}\mathbf{Q}_{1}-(\bar{e}+\bar{\eta})C^{T}\Pi
C&-(\bar{e}+\bar{\eta})C^{T}\Pi\\\ -(\bar{e}+\bar{\eta})\Pi
C&R_{2}^{-1}-(\bar{e}+\bar{\eta})\Pi\end{bmatrix}\geq 0.$ (64)
$\displaystyle\mathcal{S}_{L}=\leavevmode\resizebox{352.31624pt}{}{$\begin{bmatrix}\mathbf{Q}_{1}&\mathbf{Q}_{1}F+XC&\mathbf{Q}_{1}R_{1}&XR_{2}\\\
(\mathbf{Q}_{1}F+XC)^{T}&\mathbf{Q}_{1}&0&0\\\
R_{1}\mathbf{Q}_{1}&0&R_{1}&0\\\ R_{2}X^{T}&0&0&R_{2}\end{bmatrix}$}\geq 0.$
(65)

where $\bar{e}$ is computed from the lower incomplete Gamma function such that
$\text{Pr}[e_{k}^{T}\mathbf{P}_{e}^{-1}e_{k}\leq\bar{e}]=p_{e}$, and $p_{e}$
is the desired trucation probability.

Proof: The inequality $\mathcal{H}_{L}$ already provides a lower bound on
$\Pi$. We now add an additional upper constraint to sandwich and fully
constrain $\Pi$. We exploit the fact that $r_{k}^{T}\Pi r_{k}\leq 1$, and
$r_{k}=Ce_{k}+\eta_{k}$, therefore we conclude,

$\Phi_{1}=\begin{bmatrix}e_{k}^{T}&\eta_{k}^{T}\end{bmatrix}\begin{bmatrix}C^{T}\Pi
C&C^{T}\Pi\\\ \Pi C&\Pi\end{bmatrix}\begin{bmatrix}e_{k}\\\
\eta_{k}\end{bmatrix}\leq 1.$ (66)

In addition, for the estimation error in the absence of attack we have,
$e_{k}^{T}\mathbf{P}_{e}^{-1}e_{k}\leq\bar{e}$, where $\bar{e}$ is chosen to
contain most of the Gaussian estimation error distribution (this quadratic
form follows a chi-squared distribution, therefore, $\bar{e}$ can be chosen by
the lower incomplete Gamma function). For the measurement noise we have,
$\eta_{k}^{T}R_{2}^{-1}\eta_{k}\leq\bar{\eta}$, therefore, together we have

$\Phi_{2}=\begin{bmatrix}e_{k}^{T}&\eta_{k}^{T}\end{bmatrix}\begin{bmatrix}\mathbf{P}_{e}^{-1}&0\\\
0&R_{2}^{-1}\end{bmatrix}\begin{bmatrix}e_{k}\\\
\eta_{k}\end{bmatrix}\leq\bar{e}+\bar{\eta}.$ (67)

Hence, we can provide a upper bound for $\Pi$ since $\Phi_{1}\leq 1$, then
$(\bar{e}+\bar{\eta})\Phi_{1}\leq\Phi_{2}$. Thus,

$\begin{bmatrix}\mathbf{P}_{e}^{-1}-(\bar{e}+\bar{\eta})C^{T}\Pi
C&-(\bar{e}+\bar{\eta})C^{T}\Pi\\\ -(\bar{e}+\bar{\eta})\Pi
C&R_{2}^{-1}-(\bar{e}+\bar{\eta})\Pi\end{bmatrix}\geq 0$ (68)

which can be written as (64). We now have a upper bound on $\Pi$ in terms of
$\mathbf{Q}_{1}$, however, $\mathbf{Q}_{1}$ is not itself constrained. To do
this, we relax the Lyapunov equation (9) equality to inequality,

$\mathbf{P}_{e}-(F-LC)\mathbf{P}_{e}(F-LC)^{T}-R_{1}-LR_{2}L^{T}\geq 0.$ (69)

We apply the Schur Complement to receive

$\mathcal{S}=\begin{bmatrix}\mathbf{P}_{e}-LR_{2}L^{T}-R_{1}&(F-LC)\mathbf{P}_{e}\\\
\mathbf{P}_{e}(F-LC)^{T}&\mathbf{P}_{e}\end{bmatrix}\geq 0.$ (70)

Using the transformation
$T_{4}=\text{diag}\begin{bmatrix}\mathbf{P}_{e}^{-1}&\mathbf{P}_{e}^{-1}\end{bmatrix}$
gives

$\mathcal{S}_{L}=T_{4}\mathcal{S}T_{4}^{T}=\begin{bmatrix}\mathbf{P}_{eL}-R_{eL}&A_{eL}\\\
A_{eL}^{T}&\mathbf{P}_{eL}\end{bmatrix}\geq 0,$ (71)

where,

$\displaystyle\mathbf{P}_{eL}=\mathbf{P}_{e}^{-1}=\mathbf{Q}_{1},$ (72)
$\displaystyle
R_{eL}=\mathbf{P}_{e}^{-1}(R_{1}+LR_{2}L^{T})\mathbf{P}_{e}^{-1}=\mathbf{Q}_{1}R_{1}\mathbf{Q}_{1}+XR_{2}X^{T},$
$\displaystyle A_{eL}=\mathbf{P}_{e}^{-1}(F-LC)=\mathbf{Q}_{1}F+XC.$

Applying the Schur Complement again yields the expression in (65).
$\blacksquare$

With these combined linearization steps, we can provide Theorem 3.

######  Theorem 3.

Consider a LTI system (1) with desired output covariance constrained
$\|H\|_{2}$ gain $\bar{\gamma}$ (27), chi-squared detector threshold $\alpha$
(11) and zero-alarm stealthy attacker (13). If there exist $a,a_{2}\in[0,1)$
then the solution of the following convex optimization provides the optimal
observer $L$ and controller $K$ gain matrices that minimize the set of states
reachable by an attacker while maintaining an OCC $\|H\|_{2}$ gain no bigger
than $\bar{\gamma}$,

$\displaystyle\left\\{\begin{aligned}
&\min_{{a_{1},\mathbf{P}_{x},\mathbf{Q}_{1},\Pi}\atop{X,Y,Z}}\mathbf{tr}(\sigma\mathbf{P}_{x})\\\
&\qquad\text{s.t.}\quad 0\leq a_{1}<1,\quad a_{1}+a_{2}\geq a,\\\
&\qquad\qquad\mathcal{H}_{L}\geq 0,\ \mathcal{C}_{h}\leq 0,\
\mathcal{C}_{L}\geq 0,\\\ &\qquad\qquad\mathcal{S}_{L}\geq 0,\
\mathcal{X}_{L}\geq 0.\end{aligned}\right.$ (73)

This optimization is solved for a specific value of $\sigma$, which is
optimized through a bisection algorithm.

##  5 Case study

We consider a LTI system (with matrices given below) for this study with the
chi-squared detector tuned to a false alarm rate $\mathcal{A}=0.05$ (5%),
system noise truncated with $\bar{p}=95\%$, and a worst acceptable OCC
$\|H\|_{2}$ gain of $\bar{\gamma}=8.75$. We use CVX to solve the convex
optimizations [2].

$\displaystyle F=\begin{bmatrix}[r]1.0444&-0.1409\\\
0.3001&0.6327\end{bmatrix},\ G=\begin{bmatrix}2&3\\\ 1&1\end{bmatrix},\
C=\begin{bmatrix}2&2\\\ 1&2\end{bmatrix},$ $\displaystyle
R_{1}=\begin{bmatrix}[r]0.0183&-0.0218\\\ -0.0218&0.0261\end{bmatrix},\
R_{2}=\begin{bmatrix}0.0018&0.0031\\\ 0.0031&0.0096\end{bmatrix}$

We start by computing the optimal OCC $\|H\|_{2}$ gain $\gamma^{*}$ using
Theorem 1,

$\gamma^{*}=\sqrt{\frac{\mathrm{tr}\big{(}C\ \mathbf{P}_{x}\
C^{T}\big{)}+\mathrm{tr}(R_{2})}{\mathrm{tr}(R_{2})+\mathrm{tr}(R_{1})}}=1.5673.$
(74)

The solution for gain matrices $L$ and $K$ based on the generalized algebraic
Ricatti equation (40), returns ${4\choose 2}=6$ different answers of which two
are real valued. Of the real solutions,

$L=\begin{bmatrix}[r]1.0085&-0.9780\\\ -0.0139&0.2664\end{bmatrix},\quad
K=\begin{bmatrix}[r]0.1273&-2.0544\\\ -0.4303&1.4190\end{bmatrix},$ (75)

provides the smallest reachable set according to Lemma 2. To examine the
accuracy of Theorem 1, we compare the state covariance that is the solution of
the Lyapunov equation (29) with the designed gains in (75), with the decision
variable $\mathbf{P}_{x}$ of the optimization. In this case, the (entry-wise)
mean absolute error (MAE) between these matrices is less than $0.00001$,
justifying that the relaxed inequality in Theorem 1 recovers an exact Lyapunov
equation solution.

The open loop OCC $\|H\|_{2}$ gain is found using Lemma 4 to be,

$\gamma_{0}=\sqrt{\frac{\mathrm{tr}\big{(}C^{T}C\ \mathbf{P}_{x}\
\big{)}+\mathrm{tr}(R_{2})}{\mathrm{tr}(R_{2})+\mathrm{tr}(R_{1})}}=10.1899.$
(76)

Next we use the iterative approach from Theorem 2 with convergence criteria
$\epsilon=0.03$ for the gain matrix $L$. We use bisection to find the smallest
feasible $\sigma$ (stopping tolerance of 0.01), starting with
$\sigma_{min}=0.1$ and $\sigma_{max}=10^{6}$. In addition we increase
$\bar{\gamma}$ gradually from $\gamma^{*}=1.57$ until $\bar{\gamma}=8.75$ by
steps of $0.1$ (see Remark 5). The optimal gains are,

$L=\begin{bmatrix}0.0956&-0.1248\\\ -0.1010&0.1321\end{bmatrix},\quad
K=\begin{bmatrix}0.1440&-2.0390\\\ -0.4441&1.4063\end{bmatrix},$ (77)

with corresponding magnification factor $\sigma=400.91$. Given these gains the
actual performance is computed as $\gamma=6.91$, which is computed by
substituting the optimal gains into and then solving the Lyapunov equation
(29) and then computing $\gamma$ from (36). This discrepancy between the input
$\bar{\gamma}$ and the outcome $\gamma$ is the price of the linearization,
caused by the competing LMI constraints (see Appendix E).

Finally we optimize the gains $L$ and $K$ using the fully convexified method
in Theorem 3, which results in,

$L=\begin{bmatrix}0.1274&-0.1737\\\ -0.2019&0.2872\end{bmatrix},\quad
K=\begin{bmatrix}0.1902&-1.9945\\\ -0.4757&1.3759\end{bmatrix}.$ (78)

The actual performance is $\gamma=6.81$ (as above). The bisection is
configured as above and returns $\sigma=2471.42$ as the smallest feasible
value of the magnification factor. To examine the accuracy of these
linearization method we compare the state covariance of the system which is
the solution of the Lyapunov equation (as discussed above) with the designed
gains in (78), with the decision variable $\mathbf{P}_{x}$ of the
optimization. The mean absolute error between these two matrices is $0.04$. As
before, this discrepancy is the result of the linearization (see Appendix D)

Figure 2 compares the ellipsoidal bounds provided by the solutions above
including Theorem 3 (red) and Theorem 2 (blue) for the worst allowable
performance $\bar{\gamma}=8.75$, as well as Theorem 1 (black) for optimal
performance $\gamma=\gamma^{*}$. In this figure the ellipsoidal bounds are
from Lemma 2 with optimal gains provided from Theorem 3, Theorem 2, and
Theorem 1. This figure demonstrates that despite the similarities between
Robust Contol and Resilient Control design, the $\|H\|_{2}$ optimal solution
is not the optimal security solution. It also shows that the fully convexified
solution evidences only marginally degraded solution quality, despite the
linearization steps.

Figure 2: A comparison between the ellipsoidal bounds of the optimized
reachable sets from Theorems 2 and 3 for desired performance
$\bar{\gamma}=8.75$ along with that of the OCC $\|H\|_{2}$ optimal solution.

In Figure 3 we compute the optimal gains for a wide range of worst allowable
OCC $\|H\|_{2}$ gains using Theorem 2 and Theorem 3. Given the optimal gains,
$L$ and $K$ for each, we plot the objective function of Lemma 2 with respect
to the actual performance, $\gamma$. We use the objective function of Lemma 2
as a proxy for security where, smaller values imply better security. This plot
demonstrates that there is, indeed, a trade off between security and
performance. The nonlinearity and steep slope near $\gamma^{*}$ indicates that
dramatic improvements in security can be gained by marginal concessions in
performance. The relationship between these two properties, captured by this
plot, should be used in resilient control system design.

Figure 3: Performance-Security curve for the proposed iterative semi-definite
optimization approach (dashed green) compared with fully convexified version
(green). The extra steps taken to linearize the problem make the fully convex
optimization feasible only for a portion of the entire trade off interval.

Figure 3 shows that the solution of Theorem 3 does not cover the entire trade
off interval and the problem is infeasible on the interval
$\gamma\in[\gamma^{*},\ \gamma_{c}]$, where $\gamma_{c}=5.66$ for this
scenario. The manifold $\mathbf{P}_{\hat{x}}=\mathbf{P}_{x\hat{x}}$ imposed by
Imposed Structure 2 is not able to provide a feasible solution for entire
trade off interval, satisfying both the stability constraint in both attacked
system LMI ($\mathcal{H}_{L}$) and attack free LMI ($\mathcal{C}_{L}$). There
exists an infimum for the optimal threshold on this manifold, which is
computed in Appendix D.

##  6 Conclusion

This paper presents a set of tools to design the feedback controller and
observer gains for observer-based feedback control with the aim to minimize
the effect sensor falsification attacks. As past work has observed, there is a
necessary trade-off between ensuring performance and minimizing the effect of
the attacker. Here the attacker impact is quantified as the set of states
reachable through the action of the attacker and the nominal closed-loop
performance is specified by an output covariance constrained $\|H\|_{2}$ gain.
We frame this problem as an LMI wrapped in an iterative algorithm as well as a
fully convexified optimization and a large part of the effort here is to
linearize the constraints involved. Along the way, we also contribute a convex
optimization to design optimal OCC $\|H\|_{2}$ gains.

## References

  * [1] S. Boyd, L. El Ghaoui, E. Feron, and V. Balakrishnan.  Linear Matrix Inequalities in System and Control Theory, volume 15 of Studies in Applied Mathematics.  SIAM, Philadelphia, PA, 1994.
  * [2] Inc. CVX Research.  CVX: Matlab software for disciplined convex programming, version 2.0.  http://cvxr.com/cvx, August 2012.
  * [3] Maurıcio C De Oliveira, José C Geromel, and Jacques Bernussou.  Extended h 2 and h norm characterizations and controller parametrizations for discrete-time systems.  International Journal of Control, 75(9):666–679, 2002.
  * [4] José C Geromel, Jacques Bernussou, and Maurício C De Oliveira.  H/sub 2/-norm optimization with constrained dynamic output feedback controllers: decentralized and reliable control.  IEEE Transactions on Automatic Control, 44(7):1449–1454, 1999.
  * [5] Navid Hashemi, Carlos Murguia, and Justin Ruths.  A comparison of stealthy sensor attacks on control systems.  In 2018 Annual American Control Conference (ACC), pages 973–979. IEEE, 2018.
  * [6] Navid Hashemi and Justin Ruths.  Gain design via lmis to minimize the impact of stealthy attacks.  In 2020 Annual American Control Conference (ACC). IEEE, 2020.
  * [7] Jezdimir Milosevic, David Umsonst, Henrik Sandberg, and Karl Henrik Johansson.  Quantifying the impact of cyber-attack strategies for control systems equipped with an anomaly detector.  In 2018 European Control Conference (ECC), pages 331–337. IEEE, 2018.
  * [8] Y. Mo and B. Sinopoli.  On the performance degradation of cyber-physical systems under stealthy integrity attacks.  IEEE Transactions on Automatic Control, 61:2618–2624, 2016.
  * [9] Carlos Murguia and Justin Ruths.  On model-based detectors for linear time-invariant stochastic systems under sensor attacks.  IET Control Theory & Applications, 13(8):1051–1061, 2019.
  * [10] Carlos Murguia, Iman Shames, Justin Ruths, and Dragan Nešić.  Security metrics and synthesis of secure control systems.  Automatica, 115:108757, 2020.
  * [11] Carlos Murguia, Nathan van de Wouw, and Justin Ruths.  Reachable sets of hidden cps sensor attacks: Analysis and synthesis tools.  In proceedings of the IFAC World Congress, 2016.
  * [12] Mario A Rotea.  The generalized h2 control problem.  Automatica, 29(2):373–385, 1993.
  * [13] David Umsonst and Henrik Sandberg.  Anomaly detector metrics for sensor data attacks in control systems.  In 2018 Annual American Control Conference (ACC), pages 153–158. IEEE, 2018.
  * [14] Fangwei Xu, Kwan Ho Lee, and Biao Huang.  Monitoring control performance via structured closed-loop response subject to output variance/covariance upper bound.  Journal of Process Control, 16(9):971–984, 2006.
  * [15] Guoming Zhu, MA Rotea, and R Skelton.  A convergent algorithm for the output covariance constraint control problem.  SIAM Journal on Control and Optimization, 35(1):341–361, 1997.

##  Appendix A Validity of Relaxing the Lyapunov Equation

Here we show more precisely the objective function
$\mathbf{tr}(C\mathbf{P}_{x}C^{T})$ forces the relaxed inequality constraint
(42) to its boundary (equality), hence the optimization finds a solution that
satisfies the original Lyapunov equation (29). To see this, consider the
function,

$f(\mathbf{P})=\mathbf{P}-\hat{A}\mathbf{P}\hat{A}^{T}-\hat{R}.$ (79)

and assume there exists matrices $\mathbf{P}\geq 0$ a solution to the
optimization (38) and $\bar{\mathbf{P}}\geq 0$ a solution of the Lyapunov
equation (29) such that,

$\displaystyle f(\mathbf{P})$
$\displaystyle=\mathbf{P}-\hat{A}\mathbf{P}\hat{A}^{T}-\hat{R}\geq 0,$ (80)
$\displaystyle f(\bar{\mathbf{P}})$
$\displaystyle=\bar{\mathbf{P}}-\hat{A}\bar{\mathbf{P}}\hat{A}^{T}-\hat{R}=0.$
(81)

The difference $f(\mathbf{P})-f(\bar{\mathbf{P}})\geq 0$ is positive semi-
definite, so it also satisfies a Lyapunov equation with stable $\hat{A}$ which
means,

$(\mathbf{P}-\bar{\mathbf{P}})-\hat{A}(\mathbf{P}-\bar{\mathbf{P}})\hat{A}^{T}-\bar{R}=0,$
(82)

for some $\bar{R}\geq 0$. As the solution of a Lyapunov equation the
difference $\mathbf{P}-\bar{\mathbf{P}}\geq 0$ is positive semi-definite, and
so the first block corresponding to the states
$\mathbf{P}_{x}-\bar{\mathbf{P}}_{x}\geq 0$ is positive semi-definite. Thus,

$\displaystyle C\mathbf{P}_{x}C^{T}$ $\displaystyle\geq
C\bar{\mathbf{P}}_{x}C^{T},$ (83)
$\displaystyle\mathbf{tr}(C\mathbf{P}_{x}C^{T})$
$\displaystyle\geq\mathbf{tr}(C\bar{\mathbf{P}}_{x}C^{T}).$ (84)

In the absence of no constraints other than $\mathcal{C}_{L}$, this inequality
demonstrates that if the objective function
$\mathbf{tr}(C\mathbf{P}_{x}C^{T})$ is minimized, that the solution satisfies
the Lyapunov equation (29), i.e., that the optimization forces the inequality
in (42) to equality.

##  Appendix B Proof of Stability

We show that the LMI constraints $\mathcal{H}_{L}>0$ and $\mathcal{C}_{L}>0$
imply that the attacked system and nominal system are stable.

For $\mathcal{H}_{L}$, $\mathcal{H}_{L}>0$ implies $\mathcal{H}>0$ which
implies $\mathcal{P}>0$ and $a\mathcal{P}-A^{T}\mathcal{P}A>0$. Since
$a\in[0,1)$,
$(1-a)\mathcal{P}+a\mathcal{P}-A^{T}\mathcal{P}A=\mathcal{P}-A^{T}\mathcal{P}A>0$.
Thus, $A^{T}$, and hence, $A$ is stable.

Similarly, for $\mathcal{C}_{L}$, $\mathcal{C}_{L}>0$ implies $\mathcal{C}>0$
which implies $\mathbf{P}>0$ and
$\mathbf{P}-\hat{A}\mathbf{P}\hat{A}^{T}-\hat{R}>0$. Therefore,
$\mathbf{P}-\hat{A}\mathbf{P}\hat{A}^{T}>0$, which makes $\hat{A}$ stable as
well.

##  Appendix C Manifold of $\mathbf{P}_{\hat{x}}=\mathbf{P}_{x\hat{x}}$

We know that $\mathbf{P}\mathbf{Q}=I$ which means,

$\displaystyle\mathbf{P}_{x}\mathbf{Q}_{1}+\mathbf{P}_{x\hat{x}}\mathbf{Q}_{12}^{T}=I$
(85)
$\displaystyle\mathbf{P}_{x}\mathbf{Q}_{12}+\mathbf{P}_{x\hat{x}}\mathbf{Q}_{2}=0$
(86)
$\displaystyle\mathbf{P}_{x\hat{x}}^{T}\mathbf{Q}_{1}+\mathbf{P}_{\hat{x}}\mathbf{Q}_{12}^{T}=0$
(87)
$\displaystyle\mathbf{P}_{x\hat{x}}^{T}\mathbf{Q}_{12}+\mathbf{P}_{\hat{x}}\mathbf{Q}_{2}=I.$
(88)

if we assume the points located on manifold
$(\mathbf{P}_{\hat{x}}=\mathbf{P}_{x\hat{x}})$, Based on (87) we can conclude,
$\mathbf{Q}_{1}=-\mathbf{Q}_{12}$ and based on (85) we have

$\mathbf{P}_{\hat{x}}=\mathbf{P}_{x}-\mathbf{Q}_{1}^{-1}.$ (89)

Based on $e_{k}=x_{k}-\hat{x}_{k}$,

$\mathbf{P}_{e}=\mathbf{P}_{x}+\mathbf{P}_{\hat{x}}-\mathbf{P}_{x\hat{x}}-\mathbf{P}_{x\hat{x}}^{T},$
(90)

where on the mentioned manifold,

$\mathbf{P}_{e}=\mathbf{P}_{x}-\mathbf{P}_{\hat{x}}.$ (91)

Comparing (89) and (91) we can conclude, $\mathbf{P}_{e}^{-1}=\mathbf{Q}_{1}$.

##  Appendix D Infimum value of desired performance on the manifold
$\textbf{P}_{\hat{x}}=\textbf{P}_{x\hat{x}}$

To compute this infimum value, it suffices to introduce
$\tau=\bar{\gamma}^{2}$, and obtain an optimal solution of $\tau$ from the
convex optimization,

$\displaystyle\left\\{\begin{aligned}
&\min_{{a_{1},\mathbf{P}_{x},\mathbf{Q}_{1},\Pi}\atop{X,Y,Z,\tau}}\tau\\\
&\qquad\text{s.t.}\quad 0\leq a_{1}<1,\quad a_{1}+a_{2}\geq a,\\\
&\qquad\qquad\mathcal{H}_{L}\geq 0,\ \mathcal{C}_{h}\leq 0,\
\mathcal{C}_{L}\geq 0,\\\ &\qquad\qquad\mathcal{S}_{L}\geq 0,\
\mathcal{X}_{L}\geq 0,\end{aligned}\right.$ (92)

where $a,a_{2}\in[0,\ 1]$ are introduced as a grid search variables. Here we
have replaced the objective function that minimizes the reachable set with an
objective that minimizes $\bar{\gamma}$ without penalty for increasing the
reachable set size (either through the direct decision variables or through
$\sigma$). The scaling factor $\sigma$ for this optimization is fixed as a
large value (to make the optimization feasible, but there is no need to solve
for the smallest value of $\sigma$ here). For this case, the optimal
$\bar{\gamma}=\bar{\gamma}_{c}=7.75$ and then based on the optimal gains $L$
and $K$ the actual performance is $\gamma_{c}=5.66$. The jump at $\gamma_{c}$
in Figure 3 is because the system becomes marginally stability.

##  Appendix E Evaluation of optimal performance from Theorem 2 and 3

The cost we pay for linearization in this work is that the optimizations in
Theorems 2 and 3, unlike Theorem 1, cannot guarantee that constraints
$\mathcal{C}_{h}$ and $\mathcal{C}_{L}$ meet their boundary perfectly due to
the addition of other competing inequalities. Therefore, the solution of the
decision variable $\mathbf{P}_{x}$ is slightly different from what would be
computed from the Lyapunov equation (29) using the solution gains $K$ and $L$.
Consequently the performance $\gamma$ is different from $\bar{\gamma}$. The
$\mathcal{C}_{h}$ constraint, however, ensures that the actual performance,
$\gamma$, will be smaller than the specified performance, $\bar{\gamma}$.

Figure 4 and Figure 5 depict the difference between $\gamma$ and
$\bar{\gamma}$ in the iterative and fully convexified optimizations,
respectively. In addition, these figures can be used as a tool to determine
what desired performance should be imposed on the linearized problem in order
to receive expected performance in original nonlinear problem.

Figure 4: Optimal performance versus specified performance in the solutions of
Theorem 2. In the vicinity of $\gamma=\gamma_{0}$, we see that $\bar{\gamma}$
diverges. Figure 5: Optimal performance versus specified performance in the
solutions of Theorem 3. Here $\bar{\gamma}_{c}=7.75$ and $\gamma_{c}=5.66$. In
the vicinity of $\gamma=\gamma_{0}$, we see that $\bar{\gamma}$ diverges.