# Robust Data-Enabled Predictive Control: Tractable Formulations and
Performance Guarantees

Linbin Huang†, Jianzhe Zhen†, John Lygeros, and Florian Dörfler †: The first
two authors contributed equally to this work.The authors are with the
Department of Information Technology and Electrical Engineering at ETH Zürich,
Switzerland. (Emails: {linhuang, jizhen, jlygeros, dorfler}@ethz.ch)This
research was supported by the SNSF under NCCR Automation, the ERC under the
project OCAL (grant number 787845), and ETH Zurich Funds.

###### Abstract

We introduce a general framework for robust data-enabled predictive control
(DeePC) for linear time-invariant (LTI) systems. The proposed framework
enables us to obtain model-free optimal control for LTI systems based on noisy
input/output data. More specifically, robust DeePC solves a min-max
optimization problem to compute the optimal control sequence that is resilient
to all possible realizations of the uncertainties in the input/output data
within a prescribed uncertainty set. We present computationally tractable
reformulations of the min-max problem with various uncertainty sets.
Furthermore, we show that even though an accurate prediction of the future
behavior is unattainable in practice due to inaccessibility of the perfect
input/output data, the obtained robust optimal control sequence provides
performance guarantees for the actually realized input/output cost. We further
show that the robust DeePC generalizes and robustifies the regularized DeePC
(with quadratic regularization or 1-norm regularization) proposed in the
literature. Finally, we demonstrate the performance of the proposed robust
DeePC algorithm on high-fidelity, nonlinear, and noisy simulations of a grid-
connected power converter system.

###### Index Terms:

Data-driven control, predictive control, regularization, robust control,
robust optimization.

##  I Introduction

Data-driven control seeking an optimal control strategy from data is
attracting increasing interest from both academia and industry. Compared to
conventional model-based control, data-driven control has the advantage that
it can be applied in scenarios where data is readily available, but the system
and uncertainty models are too complex to obtain or maintain, e.g., large-
scale power systems or energy efficient buildings.

There are mainly two paradigms of data-driven control: 1) indirect data-driven
control that first identifies a model and then conducts control design based
on the identified model, and 2) direct data-driven control that circumvents
the step of system identification and obtains control policies directly from
data. Indirect data-driven control has a long history in control, where many
methods have been developed for identifying the system model, e.g., prediction
error methods (PEMs), maximum likelihood methods, and subspace methods [1];
the subsequent optimal control can be conducted using, for instance, model
predictive control (MPC), linear–quadratic–Gaussian control, and system-level
synthesis [2, 3, 4, 5, 6]. Direct data-driven control gained increasing
attention and became popular thanks to iterative feedback tuning [7], virtual
reference feedback tuning [8], reinforcement learning [9], etc. The connection
between indirect data-driven control and direct data-driven control was
investigated in [10], which provides new insights from the perspective of
regularization and convex relaxations in optimization. A central promise is
that direct data-driven control may have higher flexibility and better
performance than indirect data-driven control thanks to the data-centric
representation that avoids using a specific model from identification.
Moreover, it is generally difficult to map uncertainty specifications from
system identification over to robust control in indirect data-driven control,
while, as we will show in this paper, this may become easier in direct data-
driven control.

In recent years, a result originally formulated by Willems and co-authors in
the context of behavioral system theory has received renewed attention in
direct data-driven control. This result, known as the Fundamental Lemma [11],
shows that the subspace of input/output trajectories of a linear time-
invariant (LTI) system can be obtained from the column span of a data Hankel
matrix, thereby avoiding a parametric system representation. This result was
extended in [12] and [13] to consider mosaic Hankel matrices, (Chinese) Page
matrices, and trajectory matrices as data-driven predictors. Recently,
multiple direct data-driven control methods have been proposed based on the
Fundamental Lemma, e.g., [14, 15, 16, 17, 18, 19, 20, 21].

Here we concentrate on the data-enabled predictive control (DeePC) proposed in
[14]. In the spirit of the Fundamental Lemma, the DeePC algorithm relies only
on input/output data to learn the behavior of the unknown system and perform
safe and optimal control to drive the system along a desired trajectory using
real-time feedback. The DeePC algorithm has been successfully applied in many
scenarios, including power systems [22, 23], motor drives [24], and
quadcopters [25].

When perfect (noiseless and uncorrupted) input/output data is accessible,
DeePC can accurately predict the future behavior of the system thanks to the
Fundamental Lemma. In this case, DeePC has equivalent closed-loop behavior to
conventional MPC with a model and perfect state estimation [14]. However, in
practice, perfect data is in general not accessible to the controller due to
measurement noise and process noise, which leads to inaccurate estimations and
predictions and may degrade the quality of the obtained optimal control
sequence. In fact, a key question for data-driven control is: how does the
system perform when applying control policies computed from noisy data?
Usually, uncertainty quantification is needed in system identification and
robustification is needed in the design of the control policy, e.g., through
regularization [16, 18], computation of reachable regions [26], or robust
system level synthesis [27]. In practice, the uncertainty descriptions of
identification and robust control are often incompatible. For DeePC, it has
been frequently observed that regularization is very important to ensure good
performance under noisy measurements. Theoretical support for this observation
was provided in [18] and [28], where it was shown that regularization provides
distributional robustness against stochastic disturbances. The authors of [16]
showed that a quadratic regularization is also essential for stability.
Moreover, regularization leads to a convex relaxation of the indirect (first
identify, then control) approach and accounts for an implicit identification
step [10]. In our recent work [29], we showed that including a quadratic
regularization in DeePC enables the reformulation as a min-max problem, which
provides robustness against uncertainties in the output data. Regularization
was also linked to a min-max formulation in [27] for system level synthesis.
Nonetheless, it still remains unclear whether the realized input/output cost
can be guaranteed by applying the control sequence computed from noisy data,
and whether different structural assumptions on the uncertainties (e.g., with
Hankel structures) can be taken into account to reduce the conservativeness
induced from coarse robustification.

To this end, we present a robust DeePC framework that involves solving a min-
max optimization problem to robustify the control sequence against
uncertainties in the input/output data. We justify our min-max formulation by
showing that it enables a performance guarantee for the realized input/output
cost of the system. We consider different uncertainty sets as tight estimates
of input/output uncertainties when different types of data matrices (e.g.,
Hankel or Page matrices) are used as predictors. For instance, to make the
considered set tight, one may incorporate Hankel structures on uncertainties
when Hankel matrices are employed, and may consider column-wise uncertainties
when Page matrices or trajectory matrices are used. For such uncertainty
structures we explicitly show how the min-max problems can be reduced to
tractable minimization problems. In particular, we discuss the connections
between the robust DeePC and (previously proposed and novel) regularized DeePC
algorithms, and illustrate how different uncertainty sets lead to different
regularization terms, extending our result on quadratic regularization in
[29].

The rest of this paper is organized as follows: in Section II we give a brief
review on the DeePC algorithm. Section III presents the robust DeePC
framework. Section IV derives tractable formulations for robust DeePC under
different uncertainty sets. In Section V we discuss the connections between
robustness and regularization. Section VI tests the robust DeePC with high-
fidelity simulations on a power converter system. We conclude the paper in
Section VII.

## Notation

Let $\mathbb{N}$ denote the set of positive integers, $\mathbb{R}$ denote the
set of real numbers, $\mathbb{R}_{\geq 0}$ the set of nonnegative real
numbers, $\mathbb{R}^{n}$ the $n$-dimensional Euclidean space,
$\mathbb{R}^{m\times n}$ the set of real $m$-by-$n$ matrices, and
$\mathbb{S}^{n}$ $(\mathbb{S}^{n}_{>0})$ the set of real $n$-by-$n$ symmetric
(positive definite) matrices. We denote the index set with cardinality
$n\in\mathbb{N}$ as $[n]$, that is, $[n]=\\{1,...,n\\}$. We use $\|x\|$
($\|A\|$) to denote the (induced) 2-norm of the vector $x$ (the matrix $A$);
we use $\|x\|_{1}$ to denote the 1-norm of the vector $x$; for a vector $x$,
we use ${\left\|x\right\|_{A}^{2}}$ to denote $x^{\top}Ax$; for a matrix $A$,
we use $\|A\|_{Q}$ to denote $\|Q^{\frac{1}{2}}A\|$; we use $\|A\|_{F}$ to
denote the Frobenius norm of the matrix $A$. We use $\boldsymbol{1}_{n}$ to
denote a vector of ones of length $n$, $\boldsymbol{1}_{m\times n}$ to denote
a $m$-by-$n$ matrix of ones, and $I_{n}$ to denote an $n$-by-$n$ identity
matrix (abbreviated as $I$ when the dimensions can be inferred from the
context). We use $A^{+}$ to denote the right inverse of $A$, and
$A^{\bot}=I-A^{+}A$ to denote the orthogonal projector onto the kernel of $A$.
We use ${\rm col}(Z_{0},Z_{1},...,Z_{\ell})$ to denote the matrix
$[Z_{0}^{\top}\;Z_{1}^{\top}\;\cdots\;Z_{\ell}^{\top}]^{\top}$. We use
$\otimes$ to denote the Kronecker product.

##  II Data-Enabled Predictive Control

###  II-A Notation and Preliminaries on the Fundamental Lemma

Consider a discrete-time linear time-invariant (LTI) system

$\left\\{\begin{array}[]{l}{x_{t+1}}=A{x_{t}}+B{u_{t}}\\\
{y_{t}}=C{x_{t}}+D{u_{t}},\end{array}\right.\,$ (1)

where $A\in\mathbb{R}^{n\times n}$, $B\in\mathbb{R}^{n\times m}$,
$C\in\mathbb{R}^{p\times n}$, $D\in\mathbb{R}^{p\times m}$,
$x_{t}\in\mathbb{R}^{n}$ is the state of the system at $t\in\mathbb{Z}_{\geq
0}$, $u_{t}\in\mathbb{R}^{m}$ is the input vector, and
$y_{t}\in\mathbb{R}^{p}$ is the output vector. Recall the respective extended
observability matrix and convolution (impulse-response) matrices

$\mathscr{O}_{\ell}(A,C):={\rm col}(C,CA,...,CA^{\ell-1})\,,\quad{\rm
and}\vspace{-0.5mm}$ $\mathcal{T}_{N}=\begin{bmatrix}D&0&0&\cdots&0\\\
CB&D&0&\cdots&0\\\ CAB&CB&D&\cdots&0\\\ \vdots&\vdots&\vdots&\ddots&\vdots\\\
CA^{N-2}B&CA^{N-3}B&CA^{N-4}B&\cdots&D\\\ \end{bmatrix}\,.$ (2)

The lag of the system (1) is defined by the smallest integer
$\ell\in\mathbb{Z}_{\geq 0}$ such that the observability matrix
$\mathscr{O}_{\ell}(A,C)$ has rank $n$, i.e., the state can be reconstructed
from $\ell$ measurements. In a data-driven setting, $\ell$ and $n$ are
generally unknown, but upper bounds can usually be inferred from knowledge of
the system. Consider $L,T\in\mathbb{Z}_{\geq 0}$ with $T\geq L>\ell$ and
length-$T$ input and output trajectories of (1): $u={\rm
col}(u_{0},u_{1},\dots u_{T-1})\in\mathbb{R}^{mT}$ and $y={\rm
col}(y_{0},y_{1},\dots y_{T-1})\in\mathbb{R}^{pT}$. For the inputs $u$, define
the Hankel matrix of depth $L$ as

$\mathscr{H}_{L}(u):=\left[{\begin{array}[]{*{20}{c}}{{u_{0}}}&{{u_{1}}}&\cdots&{{u_{T-L}}}\\\
{{u_{1}}}&{{u_{2}}}&\cdots&{{u_{T-L+1}}}\\\ \vdots&\vdots&\ddots&\vdots\\\
{{u_{L-1}}}&{{u_{L}}}&\cdots&{{u_{T-1}}}\end{array}}\right]\,.$ (3)

Accordingly, for the outputs define the Hankel matrix $\mathscr{H}_{L}(y)$.
Consider the stacked matrix
$\mathscr{H}_{L}(u,y)=\left[\begin{smallmatrix}\mathscr{H}_{L}(u)\\\
\mathscr{H}_{L}(y)\end{smallmatrix}\right]$. By [12, Corollary 19], the
restricted behavior (of length $L$) equals the image of $\mathscr{H}_{L}(u,y)$
if and only if rank$\left(\mathscr{H}_{L}(u,y)\right)=mL+n$. Note that this
result extends and includes the original Fundamental Lemma [11, Theorem 1]
which requires controllability and persistency of excitation of order $L+n$
(i.e., $\mathscr{H}_{L+n}(u)$ must have full row rank) as sufficient
conditions.

These behavioral results can be leveraged for data-driven prediction and
estimation as follows. Consider $T_{\rm ini},N,T\in\mathbb{Z}_{\geq 0}$, as
well as an input/output time series ${\rm
col}(u^{\rm{d}},y^{\rm{d}})\in\mathbb{R}^{(m+p)T}$ so that
rank$\left(\mathscr{H}_{T_{\rm ini}+N}(u^{\rm{d}},y^{\rm{d}})\right)=m(T_{\rm
ini}+N)+n$. Here the superscript “d” denotes data collected offline, and the
rank condition is met by choosing $u^{\rm{d}}$ to be persistently exciting of
sufficiently high order. The Hankel matrix $\mathscr{H}_{T_{\rm
ini}+N}(u^{\rm{d}},y^{\rm{d}})$ can be partitioned into

$\left[{\begin{array}[]{*{20}{c}}{{\bar{U}_{\rm P}}}\\\ {{\bar{U}_{\rm
F}}}\end{array}}\right]:=\mathscr{H_{c}}_{T_{\rm
ini}+N}(u^{\rm{d}})\quad\text{and}\quad\left[{\begin{array}[]{*{20}{c}}{{\bar{Y}_{\rm
P}}}\\\ {{\bar{Y}_{\rm F}}}\end{array}}\right]:=\mathscr{H_{c}}_{T_{\rm
ini}+N}(y^{\rm{d}})\,,$

where $\bar{U}_{\rm P}\in\mathbb{R}^{mT_{\rm ini}\times H_{c}}$, $\bar{U}_{\rm
F}\in\mathbb{R}^{mN\times H_{c}}$, $\bar{Y}_{\rm P}\in\mathbb{R}^{pT_{\rm
ini}\times H_{c}}$, $\bar{Y}_{\rm F}\in\mathbb{R}^{pN\times H_{c}}$, and
$H_{c}=T-T_{\rm ini}-N+1$. In the sequel, the data in the partition with
subscript P (for “past”) will be used to implicitly estimate the initial
condition of the system, whereas the data with subscript F will be used to
predict the “future” trajectories. In this case, $T_{\rm ini}$ is the length
of an initial trajectory measured in the immediate past during on-line
operation, and $N$ is the length of a predicted trajectory starting from the
initial trajectory. Recall that the image of $\mathscr{H}_{T_{\rm
ini}+N}(u^{\rm{d}},y^{\rm{d}})$ spans all length-$(T_{\rm ini}+N)$
trajectories, that is, ${\rm{col}}(u_{\rm ini},u,y_{\rm
ini},y)\in\mathbb{R}^{(m+p)(T_{\rm ini}+N)}$ is a trajectory of (1) if and
only if there exists $g\in\mathbb{R}^{H_{c}}$ so that

$\left[{\begin{array}[]{*{20}{c}}{{\bar{U}_{\rm P}}}\\\ {{\bar{Y}_{\rm P}}}\\\
{{\bar{U}_{\rm F}}}\\\ {{\bar{Y}_{\rm
F}}}\end{array}}\right]g=\left[{\begin{array}[]{*{10}{c}}{{\bar{u}_{\rm
ini}}}\\\ {{\bar{y}_{\rm ini}}}\\\ u\\\ y\end{array}}\right]\,.$ (4)

The initial trajectory ${\rm{col}}(\bar{u}_{\rm ini},\bar{y}_{\rm
ini})\in\mathbb{R}^{(m+p)T_{\rm ini}}$ can be thought of as setting the
initial condition for the future (to be predicted) trajectory
${\rm{col}}(u,y)\in\mathbb{R}^{(m+p)N}$. In particular, if $T_{\rm
ini}\geq\ell$, for every given future input trajectory $u$, the future output
trajectory $y$ is uniquely determined through (4) [30].

In addition to Hankel matrices, one can also use more input/output data to
construct (Chinese) Page matrices or trajectory matrices as data-driven
predictors [12], where the Page matrix (of depth $L$) for a signal
$u\in\mathbb{R}^{mLT^{\prime}}$ is

$\mathscr{P}_{L}(u):=\left[{\begin{array}[]{*{20}{c}}{{u_{0}}}&{{u_{L}}}&\cdots&{{u_{L(T^{\prime}-1)}}}\\\
{{u_{1}}}&{{u_{L+1}}}&\cdots&{{u_{L(T^{\prime}-1)+1}}}\\\
\vdots&\vdots&\ddots&\vdots\\\
{{u_{L-1}}}&{{u_{2L-1}}}&\cdots&{{u_{LT^{\prime}-1}}}\end{array}}\right]\,,$
(5)

and a trajectory matrix (of depth $L$) can be constructed from $T^{\prime}$
independent trajectories $u^{i}\in\mathbb{R}^{mL}$ ($i\in[T^{\prime}]$) as

$\mathscr{T}_{L}(u):=\left[{\begin{array}[]{*{20}{c}}{{u^{1}}}&{{u^{2}}}&\cdots&{{u^{T^{\prime}}}}\end{array}}\right]\,.$
(6)

Notice that the elements in a Hankel matrix are correlated (structured), while
the elements in a Page matrix or a trajectory matrix are independent
(unstructured).

###  II-B Review of the DeePC algorithm

The DeePC algorithm proposed in [14] directly uses input/output data collected
from the unknown system to predict the future behaviour, and performs optimal
and safe control without identifying a parametric system representation. More
specifically, DeePC solves the following optimization problem to obtain the
optimal future control inputs

$\begin{array}[]{cl}\displaystyle\mathop{{\rm{min}}}\limits_{g\atop(u,y)\in\mathcal{C}}\left\\{{\left\|u\right\|_{R}^{2}}+{\left\|{y-r}\right\|_{Q}^{2}}\,\Bigg{|}\,\,{\begin{bmatrix}{{\bar{U}_{\rm
P}}}\\\ {{\bar{Y}_{\rm P}}}\\\ {{\bar{U}_{\rm F}}}\\\ {{\bar{Y}_{\rm
F}}}\end{bmatrix}}g={\begin{bmatrix}{{\bar{u}_{\rm ini}}}\\\ {{\bar{y}_{\rm
ini}}}\\\ u\\\ y\end{bmatrix}}\right\\},\end{array}$ (7)

where the set of input/output constraints is defined as
$\mathcal{C}=\\{(u,y)\in\mathbb{R}^{(m+p)N}\ |\ W{\rm col}(u,y)\leq w\\}$ for
$W\in\mathbb{R}^{n_{w}\times(m+p)N}$ and $w\in\mathbb{R}^{n_{w}}$. The
positive definite matrix $R\in\mathbb{S}^{mN\times mN}_{>0}$ and positive
semi-definite matrix $Q\in\mathbb{S}^{pN\times pN}_{\geq 0}$ are the cost
matrices. The vector $r\in\mathbb{R}^{pN}$ is a prescribed reference
trajectory for the future outputs. DeePC involves solving the convex quadratic
problem (7) in a receding horizon fashion, that is, after calculating the
optimal control sequence $u^{\star}$, we apply
$(u_{t},...,u_{t+k-1})=(u_{0}^{\star},...,u_{k-1}^{\star})$ to the system for
$k\leq N$ time steps, then, reinitialize the problem (7) by updating ${\rm
col}(\bar{u}_{\rm ini},\bar{y}_{\rm ini})$ to the most recent input and output
measurements, and setting $t$ to $t+k$, to calculate the new optimal control
for the next $k\leq N$ time steps. As in MPC, the control horizon $k$ is a
design parameter.

The standard DeePC algorithm in (7) assumes perfect (noiseless and
uncorrupted) input/output data generated from the unknown system (1). However,
in practice, the perfect data is not accessible to the controller due to, for
example, measurement noise, process noise, or noise that enters through the
input channels. Instead, the corresponding measured output data and the
corresponding recorded input data are used in the control algorithm.

Perfect Data and Noisy Data. Throughout the paper, we use $\bar{U}_{\rm P}$,
$\bar{U}_{\rm F}$, $\bar{u}_{\rm ini}$, $\bar{Y}_{\rm P}$, $\bar{Y}_{\rm F}$,
and $\bar{y}_{\rm ini}$ to denote the perfect data generated from the system
(1), which accurately captures the system dynamics according to (4); we use
$\hat{U}_{\rm P}$, $\hat{U}_{\rm F}$, and $\hat{u}_{\rm ini}$ to denote the
corresponding recorded input data, and use $\hat{Y}_{\rm P}$, $\hat{Y}_{\rm
F}$, and $\hat{y}_{\rm ini}$ to denote the corresponding measured output data.

##  III Robust DeePC

In this section, we propose a novel framework for robust DeePC to incorporate
different geometry and tighter estimates of uncertainty sets in the
input/output data.

###  III-A Robust DeePC

The perfect input/output data generated from the unknown system (1) precisely
captures the system dynamics and accurately predicts the future behavior of
the system thanks to Fundamental Lemma [11]. However, in practice, the perfect
data is not accessible to the controller due to, for instance, process noise,
or measurement noise. This results in inaccurate prediction and may degrade
performance once the control sequence is applied to the system. To deal with
this problem, we propose a robust optimization framework that incorporates
uncertainties in the data matrices and the initial trajectory of (7) by
considering the following robust counterpart

$\begin{array}[]{l}\displaystyle\mathop{{\rm{min}}}\limits_{g,\sigma_{u},\sigma_{y}\atop(u,y)\in\mathcal{C}}{\left\|u\right\|_{R}^{2}}+{\left\|{y-r}\right\|_{Q}^{2}}+h(\sigma_{u},\sigma_{y})\\\
\;\;{\rm s.t.}\;\;{\begin{bmatrix}{{U_{\rm P}(\xi)}}\\\ {{Y_{\rm P}(\xi)}}\\\
{{U_{\rm F}(\xi)}}\\\ {{Y_{\rm
F}(\xi)}}\end{bmatrix}}g={\begin{bmatrix}{{u_{\rm ini}(\xi)}}\\\ {{y_{\rm
ini}(\xi)}}\\\ u\\\ y\end{bmatrix}}+{\begin{bmatrix}\sigma_{u}\\\
\sigma_{y}\\\ 0\\\ 0\end{bmatrix}},\ \forall\xi\in\mathcal{D}\end{array}$ (8)

where the prescribed uncertainty set $\mathcal{D}\subset\mathbb{R}^{n_{\xi}}$
is a compact and convex set with nonempty relative interior, and
$(\sigma_{u},\sigma_{y})$ captures the violation of the first two sets of
equality constraints. The estimation regularizer $h(\sigma_{u},\sigma_{y})$
penalizes violation of the first two sets of equality constraints due to
uncertainties, which can be seen as a penalty on the estimation uncertainty.
We assume that
$h:\mathbb{R}^{mN}\times\mathbb{R}^{pN}\rightarrow\mathbb{R}_{\geq 0}$ is a
convex function with $h(\sigma_{u},\sigma_{y})=0$ if and only if
$\sigma_{u}=0$ and $\sigma_{y}=0$.

We assume that the elements of the Hankel matrices and the initial trajectory
are affected by uncertainties in an affine way. More specifically, the
matrices ${{U_{\rm P}(\xi)}},{{Y_{\rm P}(\xi)}},{{U_{\rm F}(\xi)}},{{Y_{\rm
F}(\xi)}}$, and vectors ${{u_{\rm ini}(\xi)}}$ and ${{y_{\rm ini}(\xi)}}$ are
affine functions of the uncertain parameter $\xi\in\mathcal{D}$, with
${{U_{\rm P}(0)}}=\hat{U}_{\rm P},{{Y_{\rm P}(0)}}=\hat{Y}_{\rm P},{{U_{\rm
F}(0)}}=\hat{U}_{\rm F},{{Y_{\rm F}(0)}}=\hat{Y}_{\rm F}$, ${{u_{\rm
ini}(0)}}=\hat{u}_{\rm ini}$ and ${{y_{\rm ini}(0)}}=\hat{y}_{\rm ini}$. For
instance,

$Y_{\rm P}(\xi)=\hat{Y}_{\rm P}+\sum_{j=1}^{n_{\xi}}Y_{\rm
P}^{(j)}\xi_{j}\quad\text{and}\quad y_{\rm ini}(\xi)=\hat{y}_{\rm
ini}+\sum_{j=1}^{n_{\xi}}y^{(j)}_{\rm ini}\xi_{j},$

and $Y_{\rm P}^{(j)}\in\mathbb{R}^{pT_{\rm ini}\times H_{c}}$, $y^{(j)}_{\rm
ini}\in\mathbb{R}^{pT_{\rm ini}}$ for every $j\in[n_{\xi}]$.

In problem (8), one seeks a vector $(u,y,g,\sigma_{u},\sigma_{y})$ that
minimizes the objective function value and is feasible for all possible
uncertainties $\xi$ residing in $\mathcal{D}$.

Problem (8) is overly conservative, or even infeasible because there may not
exist a $(u,y,g,\sigma_{u},\sigma_{y})$ that satisfies the equality
constraints for all possible realization of $\xi$ in $\mathcal{D}$. To reduce
the conservativeness of (8), we propose the following two-stage robust DeePC
decision problem

$\begin{array}[]{l}\displaystyle\mathop{{\rm{min}}}\limits_{g}\max_{\xi\in\mathcal{D}}\min_{(u,y)\in\mathcal{C}\atop\sigma_{u},\sigma_{y}}{\left\|u\right\|_{R}^{2}}+{\left\|{y-r}\right\|_{Q}^{2}}+h(\sigma_{u},\sigma_{y})\\\
\;{\rm s.t.}\;\;{\begin{bmatrix}{{U_{\rm P}(\xi)}}\\\ {{Y_{\rm P}(\xi)}}\\\
{{U_{\rm F}(\xi)}}\\\ {{Y_{\rm
F}(\xi)}}\end{bmatrix}}g={\begin{bmatrix}{{u_{\rm ini}(\xi)}}\\\ {{y_{\rm
ini}(\xi)}}\\\ u\\\ y\end{bmatrix}}+{\begin{bmatrix}\sigma_{u}\\\
\sigma_{y}\\\ 0\\\ 0\end{bmatrix}}.\end{array}$ (9)

Here, we assume $g\in\mathbb{R}^{H_{c}}$ to be a first stage variable that is
decided before the realization of the disturbances $\xi$, and
$(u,y,\sigma_{u},\sigma_{y})$ to be the second stage variable which is
determined after the value of $\xi$ is revealed. Therefore, the second stage
variable $(u,y,\sigma_{u},\sigma_{y})$ is a general function of $\xi$ that
minimizes the objective function after $\xi$ is realized. Problem (9) should
be read as: the designer chooses $g$, then nature strikes and chooses an
adversarial $\xi$. In the end, the control sequence $(u,y)$ and estimation
errors $(\sigma_{u},\sigma_{y})$ realize. We will later show that this two-
stage problem provides robust control input sequence to the system and
guarantees a performance certificate for the realized input/output cost.

The coefficient of the second stage variable $(u,y,\sigma_{u},\sigma_{y})$ is
constant in (9), which corresponds to the stochastic programming format known
as fixed recourse [31]. In this case, since $(u,y,\sigma_{u},\sigma_{y})$ are
fully characterized by the uncertain linear equality constraints in (9),
linear decision rules are optimal for $(u,y,\sigma_{u},\sigma_{y})$. Without
loss of generality, one can eliminate all the second stage variables
$(u,y,\sigma_{u},\sigma_{y})$ in (9), which is equivalent to imposing linear
decision rules on $(u,y,\sigma_{u},\sigma_{y})$ [32, Lemma 2], and problem (9)
can be equivalently reformulated as the following min-max robust optimization
problem

$\begin{array}[]{cl}\displaystyle\mathop{{\rm{min}}}\limits_{g\in\mathcal{G}}\max_{\xi\in\mathcal{D}}&\left\|U_{\rm
F}(\xi)g\right\|_{R}^{2}+\left\|{Y_{\rm F}(\xi)g-r}\right\|_{Q}^{2}\\\
&+h(U_{\rm P}(\xi)g-u_{\rm ini}(\xi),Y_{\rm P}(\xi)g-y_{\rm
ini}(\xi)),\end{array}$ (10)

where $\mathcal{G}=\\{g\ |\ \forall\xi\in\mathcal{D}:\,(U_{\rm F}(\xi)g,Y_{\rm
F}(\xi)g)\in\mathcal{C}\\}$.

TABLE I: Tractable robust counterparts for $\forall\xi\in\mathcal{D}:(a+B^{\top}\xi)^{\top}x\leq b$. LO: linear optimization problems. CQO: convex quadratic optimization problems. Note that $\nu$ is an additional optimization variable in the resulting tractable robust counterparts. This table is adapted from [33]. Type (Tractability) | $\mathcal{D}$ | Robust Counterpart
---|---|---
Box (LO) | $\|\xi\|_{\infty}\leq\rho$ |  $a^{\top}x+\rho\|Bx\|_{1}\leq b$
Ellipsoidal (CQO) | $\|\xi\|\leq\rho$ |  $a^{\top}x+\rho\|Bx\|\leq b$
Budget (LO) | $D\xi\leq d$ | $\begin{cases}a^{\top}x+d^{\top}\nu\leq b\\\ D^{\top}\nu=Bx\\\ \nu\geq 0\end{cases}$
Polyhedral (LO) | $\begin{cases}\|\xi\|_{\infty}\leq\rho\\\ \|\xi\|_{1}\leq\tau\end{cases}$ | $a^{\top}x+\rho\|\nu\|_{1}+\tau\|Bx-\nu\|_{\infty}\leq b$

Note that since $\mathcal{C}$ is a polytope, and the uncertainty set
$\mathcal{D}$ is a compact and convex set with nonempty relative interior,
then each semi-infinite constraint in $\mathcal{G}$ admits a tractable robust
counterpart, that is, it can be reformulated into a finite set of convex
constraints by using standard techniques from robust optimization [34]; see
Table I for examples. Note that if some of the linear equality constraints in
(8) are certain, that is, not affected by uncertainties, then one can include
such equalities in the set $\mathcal{G}$, in a similar way as in [29], instead
of including them in the objective function of (10).

Figure 1: Interaction between the system and the robust DeePC algorithm. The
data matrices $\hat{U}_{\rm P}$, $\hat{Y}_{\rm P}$, $\hat{U}_{\rm F}$, and
$\hat{Y}_{\rm F}$ are obtained offline.

The minimizer $g^{\star}$ of (10) can be used to compute the optimal control
sequence $u^{\star}=\hat{U}_{\rm F}g^{\star}$. We illustrate the relationship
between the unknown system and the robust DeePC algorithm in Fig. 1. Note that
the realized control actions $u_{\rm sys}$ of the system may be different from
$u^{\star}$ due to input disturbances. The realized output trajectory $y_{\rm
sys}$ of the system (in response to $u_{\rm sys}$) is in general different
from the predicted trajectory $y^{\star}=\hat{Y}_{\rm F}g^{\star}$ in the
algorithm because the perfect data is not accessible to the controller. Note
that $u^{\star}$ can also be applied in a receding horizon fashion as done for
DeePC.

In the remainder of this paper, we focus on (10) with a square quadratic
estimation regularizer, that is,

$h(\sigma_{u},\sigma_{y})=\lambda_{u}\|\sigma_{u}\|^{2}+\lambda_{y}\|\sigma_{y}\|^{2},$
(11)

where $\lambda_{u}$ and $\lambda_{y}$ are positive scalars.

###  III-B Performance Guarantees Induced by Robust DeePC

As a first result, we link the realized cost of the system, defined by $c_{\rm
realized}=\left\|u_{\rm sys}\right\|_{R}^{2}+\left\|y_{\rm
sys}-r\right\|_{Q}^{2}$, to the optimization cost of (10), and show how the
realized cost can be certified by applying the robust DeePC. Here “realized
cost” refers to the cost accrued by the real system if the optimizer of (10)
is applied in open loop for the whole horizon $N$ (see Fig. 1). We assume that
the uncertainty set $\mathcal{D}$ contains the realization of the uncertain
input/output data deviations.

######  Assumption 1.

There exists a $\bar{\xi}\in\mathcal{D}$ such that ${{U_{\rm
P}(\bar{\xi})}}=\bar{U}_{\rm P}$, ${{Y_{\rm P}(\bar{\xi})}}=\bar{Y}_{\rm P}$,
${{U_{\rm F}(\bar{\xi})}}=\bar{U}_{\rm F}$, ${{Y_{\rm
F}(\bar{\xi})}}=\bar{Y}_{\rm F}$, ${{u_{\rm ini}(\bar{\xi})}}=\bar{u}_{\rm
ini}$ and ${{y_{\rm ini}(\bar{\xi})}}=\bar{y}_{\rm ini}$.

Assumption 1 ensures that the perfect input/output data is captured by the
considered uncertainty set $\mathcal{D}$. However, even with Assumption 1 in
place, it is still not obvious that the realized cost can be guaranteed by
solving (10) and then applying $u^{\star}=\hat{U}_{\rm F}g^{\star}$ to the
system. This is because the violations of the initial conditions are penalized
in (10), and thus we generally have $\bar{U}_{\rm P}g\neq\bar{u}_{\rm ini}$
and $\bar{Y}_{\rm P}g\neq\bar{y}_{\rm ini}$, hence, the decision variable $g$
cannot be used to accurately predict the realized future trajectory. In the
following theorem, we relate the optimal cost of the proposed robust DeePC
algorithm (10) to the realized cost on the unknown system.

######  Theorem III.1.

If Assumption 1 holds and $\|u_{\rm sys}-\bar{U}_{\rm
F}g^{\star}\|\leq\eta_{p}$ where $g^{\star}$ is the minimizer of (10), then
there exists sufficiently large
$(\lambda_{u},\lambda_{y})\in\mathbb{R}^{2}_{>0}$ such that

$2\sqrt{c_{\rm
opt}}+\eta_{p}\left(\sqrt{2}\left\|I_{mN}\right\|_{R}+\left\|{\mathcal{T}_{N}}\right\|_{Q}\right)\geq\sqrt{c_{\rm
realized}}\,,$

where $c_{\rm opt}$ denotes the minimum of (10), and $\mathcal{T}_{N}$ is
defined in (2).

###### Proof.

If Assumption 1 holds, it follows that

$\begin{split}c_{\rm opt}\geq&\left\|\bar{U}_{\rm
F}g^{\star}\right\|_{R}^{2}+\left\|\bar{Y}_{\rm
F}g^{\star}-r\right\|_{Q}^{2}+\lambda_{u}\left\|\bar{U}_{\rm
P}g^{\star}-\bar{u}_{\rm ini}\right\|^{2}\\\ &+\lambda_{y}\left\|\bar{Y}_{\rm
P}g^{\star}-\bar{y}_{\rm ini}\right\|^{2}\,.\end{split}$ (12)

Since the realized $(u_{\rm sys},y_{\rm sys})$ is an input/output trajectory
departing from $(\bar{u}_{\rm ini},\bar{y}_{\rm ini})$, according to the
Fundamental Lemma [11], there exists a $\bar{g}$ that satisfies

$\left[{\begin{array}[]{*{20}{c}}{\bar{U}_{\rm P}}\\\ {\bar{Y}_{\rm P}}\\\
{\bar{U}_{\rm F}}\\\ {\bar{Y}_{\rm
F}}\end{array}}\right]\bar{g}=\left[{\begin{array}[]{*{20}{c}}{\bar{u}_{\rm
ini}}\\\ {\bar{y}_{\rm ini}}\\\ u_{\rm sys}\\\ y_{\rm
sys}\end{array}}\right].$ (13)

By defining $\Delta_{g}=\bar{g}-g^{\star}$, we have

$\left[{\begin{array}[]{*{20}{c}}{\bar{U}_{\rm P}}\\\ {\bar{Y}_{\rm P}}\\\
{\bar{U}_{\rm F}}\\\ {\bar{Y}_{\rm
F}}\end{array}}\right]\Delta_{g}=\left[{\begin{array}[]{*{20}{c}}{\bar{u}_{\rm
ini}-\bar{U}_{\rm P}g^{\star}}\\\ {\bar{y}_{\rm ini}-\bar{Y}_{\rm
P}g^{\star}}\\\ u_{\rm sys}-\bar{U}_{\rm F}g^{\star}\\\ y_{\rm
sys}-\bar{Y}_{\rm
F}g^{\star}\end{array}}\right]=:\left[{\begin{array}[]{*{20}{c}}\epsilon_{u_{\rm
ini}}\\\ \epsilon_{y_{\rm ini}}\\\ \epsilon_{u_{\rm sys}}\\\ y_{\rm
sys}-\bar{Y}_{\rm F}g^{\star}\end{array}}\right]\,,$ (14)

which implies that

$\begin{split}\Delta_{g}&=\begin{bmatrix}\bar{U}_{\rm P}\\\ \bar{Y}_{\rm P}\\\
\bar{U}_{\rm F}\end{bmatrix}^{+}\begin{bmatrix}\epsilon_{u_{\rm ini}}\\\
\epsilon_{y_{\rm ini}}\\\ \epsilon_{u_{\rm
sys}}\end{bmatrix}+\begin{bmatrix}\bar{U}_{\rm P}\\\ \bar{Y}_{\rm P}\\\
\bar{U}_{\rm F}\end{bmatrix}^{\bot}x\,,\\\ \text{and}\quad\bar{Y}_{\rm
F}\Delta_{g}&=\bar{Y}_{\rm F}\begin{bmatrix}\bar{U}_{\rm P}\\\ \bar{Y}_{\rm
P}\\\ \bar{U}_{\rm F}\end{bmatrix}^{+}\begin{bmatrix}\epsilon_{u_{\rm ini}}\\\
\epsilon_{y_{\rm ini}}\\\ \epsilon_{u_{\rm sys}}\end{bmatrix}=:K\epsilon_{\rm
ini}+{\mathcal{T}_{N}}\epsilon_{u_{\rm sys}}\,,\end{split}$

where $x\in\mathbb{R}^{H_{c}}$ and $\epsilon_{\rm ini}={\rm
col}(\epsilon_{u_{\rm ini}},\epsilon_{y_{\rm ini}})$ , while the second
equality follows from $\bar{Y}_{\rm F}[{\rm col}(\bar{U}_{\rm P},\bar{Y}_{\rm
P},\bar{U}_{\rm F})]^{\bot}=0$ because the output trajectory is uniquely
determined given an initial trajectory and a future input trajectory [30].
Note that $\begin{bmatrix}K&\mathcal{T}_{N}\end{bmatrix}$ is the $N$-step
auto-regressive matrix with extra input (ARX) of the system (1). Since
$\mathcal{T}_{N}$ maps the $N$-step future inputs to the $N$-step future
outputs, it can be obtained as (2).

By definding $\Lambda:={\rm diag}(\lambda_{u}I_{mT_{\rm
ini}},\lambda_{y}I_{pT_{\rm ini}})$, it follows from (12) and (14) that

$\begin{split}c_{\rm opt}\geq&\left\|u_{\rm sys}-\epsilon_{u_{\rm
sys}}\right\|_{R}^{2}\\\ &+\left\|y_{\rm sys}-K\epsilon_{\rm
ini}-{\mathcal{T}_{N}}\epsilon_{u_{\rm
sys}}-r\right\|_{Q}^{2}+\left\|\epsilon_{\rm ini}\right\|^{2}_{\Lambda}\\\
\geq&\left\|u_{\rm sys}-\epsilon_{u_{\rm sys}}\right\|_{R}^{2}\\\
&+\frac{1}{2}\left(\left\|y_{\rm sys}-r-K\epsilon_{\rm
ini}-{\mathcal{T}_{N}}\epsilon_{u_{\rm sys}}\right\|_{Q}+\left\|\epsilon_{\rm
ini}\right\|_{\Lambda}\right)^{2}\\\ \geq&\left\|u_{\rm sys}-\epsilon_{u_{\rm
sys}}\right\|_{R}^{2}\\\ &+\frac{1}{2}\left(\left\|y_{\rm
sys}-r-{\mathcal{T}_{N}}\epsilon_{u_{\rm
sys}}\right\|_{Q}-\left\|K\epsilon_{\rm ini}\right\|_{Q}+\left\|\epsilon_{\rm
ini}\right\|_{\Lambda}\right)^{2}\\\ \geq&\left\|u_{\rm sys}-\epsilon_{u_{\rm
sys}}\right\|_{R}^{2}+\frac{1}{2}\left\|y_{\rm
sys}-r-{\mathcal{T}_{N}}\epsilon_{u_{\rm sys}}\right\|_{Q}^{2}\,.\end{split}$
(15)

Here, the second inequality follows from the simple fact that
$a^{2}+b^{2}\geq\frac{1}{2}(a+b)^{2}$, the third inequality holds thanks to
the reverse triangle inequality. The last inequality is satisfied if we take
$\lambda_{u}$ and $\lambda_{y}$ large enough (by assumption) to ensure

$\Lambda\succeq K^{\top}QK\quad\implies\quad\left\|\epsilon_{\rm
ini}\right\|_{\Lambda}\geq\left\|K\epsilon_{\rm ini}\right\|_{Q}.$

From (15), we have

$c_{\rm opt}\geq\frac{1}{4}\left(\sqrt{2}\left\|u_{\rm sys}-\epsilon_{u_{\rm
sys}}\right\|_{R}+\left\|y_{\rm sys}-r-{\mathcal{T}_{N}}\epsilon_{u_{\rm
sys}}\right\|_{Q}\right)^{2}\,,$

considering that $a^{2}+b^{2}\geq\frac{1}{2}(a+b)^{2}$. Further, we have

$\begin{split}&2\sqrt{c_{\rm opt}}\geq\sqrt{2}\left\|u_{\rm
sys}-\epsilon_{u_{\rm sys}}\right\|_{R}+\left\|y_{\rm
sys}-r-{\mathcal{T}_{N}}\epsilon_{u_{\rm sys}}\right\|_{Q}\\\
&\geq\sqrt{2}\left\|u_{\rm sys}\right\|_{R}-\sqrt{2}\left\|\epsilon_{u_{\rm
sys}}\right\|_{R}+\left\|y_{\rm
sys}-r\right\|_{Q}-\left\|{\mathcal{T}_{N}}\epsilon_{u_{\rm
sys}}\right\|_{Q}\\\ &\geq\sqrt{2}\left\|u_{\rm sys}\right\|_{R}+\left\|y_{\rm
sys}-r\right\|_{Q}-\eta_{p}\left(\sqrt{2}\left\|I_{mN}\right\|_{R}+\left\|{\mathcal{T}_{N}}\right\|_{Q}\right)\end{split}$

where the second inequality follows from the reverse triangle inequality, and
the third due to the Cauchy–Schwarz inequality. Notice that
$\sqrt{2}\left\|u_{\rm sys}\right\|_{R}+\left\|y_{\rm
sys}-r\right\|_{Q}\geq\sqrt{c_{\rm realized}}$ due to the simple fact that
$a+b\geq\sqrt{a^{2}+b^{2}}$ when $a,b\geq 0$. This completes the proof. ∎

TABLE II: Comparisons of different uncertainty sets for robust DeePC. |  Uncertainty type | Tractability |  Setting different bounds for different columns |  Removing the effects of scaling matrices |  Hankel structure
---|---|---|---|---|---
Section IV-A | Unstructured | Conic Quadratic | ✗ | ✗ | ✗
Section IV-B | Column-wise | Conic Quadratic | ✓ | ✗ | ✗
Section IV-C | Interval | Convex Quadratic | ✓ | ✓ | ✗
Section IV-D | Structured | Semi-definite Program | ✗ | ✓ | ✓

Theorem III.1 shows that by solving the min-max problem in (10), one can
obtain robust control input sequences that guarantee the realized input/output
cost of the system, even with the initial conditions violated (i.e., we
generally have $\bar{U}_{\rm P}g\neq\bar{u}_{\rm ini}$ and $\bar{Y}_{\rm
P}g\neq\bar{y}_{\rm ini}$) in the min-max problem. The choice of
$(\lambda_{u},\lambda_{y})$ depends only on $Q$ and $K$, where $K$ is a
constant matrix that is implicitly determined by the state-space matrices $A$,
$B$, $C$, and $D$. Since Theorem III.1 only requires
$(\lambda_{u},\lambda_{y})$ to be sufficiently large, one can use a guess on
the upper bound of the largest eigenvalue of $K^{\top}QK$ and choose
$(\lambda_{u},\lambda_{y})$ accordingly. The value of $\eta_{p}$ increases
with a higher level of input disturbance. In the absence of input disturbance,
that is, $u_{\rm sys}=\bar{U}_{\rm F}g^{\star}$, we construct a new upper
bound on the realized input/output cost in the following result.

######  Corollary III.2.

If Assumption 1 holds and the unknown system (1) does not have input
disturbance, then there exist sufficiently large
$(\lambda_{u},\lambda_{y})\in\mathbb{R}^{2}_{>0}$ such that

$2c_{\rm opt}\geq c_{\rm realized}\,.$

###### Proof.

Since $\epsilon_{u_{\rm sys}}=u_{\rm sys}-\bar{U}_{\rm F}g^{\star}=0$ by
assumption, the claim directly follows from (15). ∎

Note that with $\eta_{p}=0$, the inequality in Corollary III.2 is at least as
tight as the one from Theorem III.1.

In the absence of uncertainties in input/output data, one can consider
$\mathcal{D}=\\{0\\}$, and thus (10) reduces to

$\begin{array}[]{cl}\mathop{\min}\limits_{g}&\|\bar{U}_{\rm
F}g\|_{R}^{2}+\|\bar{Y}_{\rm F}g-r\|_{Q}^{2}\\\ &+\lambda_{u}\|\bar{U}_{\rm
P}g-\bar{u}_{\rm ini}\|^{2}+\lambda_{y}\|\bar{Y}_{\rm P}g-\bar{y}_{\rm
ini}\|^{2}\\\ {\rm s.t.}&(\bar{U}_{\rm F}g,\bar{Y}_{\rm
F}g)\in\mathcal{C}\,,\end{array}$

which recovers (7) by further considering $\bar{U}_{\rm P}g=\bar{u}_{\rm ini}$
and $\bar{Y}_{\rm P}g=\bar{y}_{\rm ini}$ as hard constraints. In this case and
in the absence of disturbances, we have $c_{\rm opt}=c_{\rm realized}$.

The minimum of (10) will in general be increasing when a larger uncertainty
set is needed to cover possible realized uncertainties, which, according to
Theorem III.1, indicates that the (worst-case) realized input/output cost may
also increase, as $c_{\rm opt}$ is always greater than $c_{\rm realized}$. In
other words, the (worst-case) performance of robust DeePC degrades with
increasing of noise level, and the best performance is achieved when perfect
data is accessible. Moreover, one should also make the uncertainty set tight
to reduce the conservativeness.

##  IV Tractable Reformulations for Robust DeePC

In what follows, we derive tractable reformulations for (10) when various
uncertainty sets are considered, and we will compare the conservativeness and
the resulting robust counterparts with the different uncertainty sets. Note
that (10) can be compactly and equivalently represented as

$\displaystyle\mathop{{\rm{min}}}\limits_{g\in\mathcal{G}}\max_{\xi\in\mathcal{D}}\
\ \left\|A(\xi)g-b(\xi)\right\|^{2}\,,$ (16)

where the elements of $A(\xi)$ and $b(\xi)$ are affine functions of $\xi$,
namely,

$\displaystyle A(\xi)=A^{(0)}+\sum_{j=1}^{n_{\xi}}A^{(j)}\xi_{j},\ {\rm
and}\,\,\,b(\xi)=b_{0}+\sum_{j=1}^{n_{\xi}}b_{j}\xi_{j}\,,$

where

$\begin{split}A^{(j)}&={\rm col}\left(\lambda_{u}^{\frac{1}{2}}U_{\rm
P}^{(j)},\lambda_{y}^{\frac{1}{2}}Y_{\rm P}^{(j)},R^{\frac{1}{2}}U_{\rm
F}^{(j)},Q^{\frac{1}{2}}Y_{\rm F}^{(j)}\right)\in\mathbb{R}^{H_{r}\times
H_{c}},\\\ b^{(j)}&={\rm col}\left(\lambda_{u}^{\frac{1}{2}}u_{\rm
ini}^{(j)},\lambda_{y}^{\frac{1}{2}}y_{\rm
ini}^{(j)},0,0\right)\in\mathbb{R}^{H_{r}},\ \forall j\in[n_{\xi}]\,,\\\
A^{(0)}&={\rm col}\left(\lambda_{u}^{\frac{1}{2}}\hat{U}_{\rm
P},\lambda_{y}^{\frac{1}{2}}\hat{Y}_{\rm P},R^{\frac{1}{2}}\hat{U}_{\rm
F},Q^{\frac{1}{2}}\hat{Y}_{\rm F}\right)\in\mathbb{R}^{H_{r}\times H_{c}},\\\
b^{(0)}&={\rm col}\left(\lambda_{u}^{\frac{1}{2}}\hat{u}_{\rm
ini},\lambda_{y}^{\frac{1}{2}}\hat{y}_{\rm
ini},0,Q^{\frac{1}{2}}r\right)\in\mathbb{R}^{H_{r}},\end{split}$

and $H_{r}=(m+p)(T_{\rm ini}+N)$.

Below we discuss various uncertainty quantifications and tractable
reformulations sorted from coarse to fine, as listed in Table II. We will
discuss the tractability and capabilities of the robust DeePC formulation (16)
when these uncertainty sets are considered.

###  IV-A Unstructured Uncertainties

Consider the following variant of (16) with unstructured uncertainties in the
following form

$\displaystyle\mathop{{\rm{min}}}\limits_{g\in\mathcal{G}}\max_{(\Delta_{A},\Delta_{b})\in\mathcal{D}_{\rm
uns}}\ \ \left\|(A^{(0)}+\Delta_{A})g-(b^{(0)}+\Delta_{b})\right\|,$ (17)

where the unstructured uncertainty set is defined through

$\mathcal{D}_{\rm uns}=\\{(\Delta_{A},\Delta_{b})\mid\|[\Delta_{A}\
\Delta_{b}]\|_{F}\leq\rho_{u}\\},$

with $\rho_{u}\in\mathbb{R}_{>0}$. Here we take the square-root of the
quadratic objective function of (16), which will not affect the minimizer of
the problem.

Note that (17) is a special (albeit unstructured) case of (16) with $\xi={\rm
col}\left((\Delta_{A})_{1},(\Delta_{A})_{2},\dots,(\Delta_{A})_{H_{c}},\Delta_{b}\right)\in\mathbb{R}^{H_{r}(H_{c}+1)}$
where $(\Delta_{A})_{i}$ is the $i$-th column of $\Delta_{A}$, and it can be
easily verified that the elements of $A(\xi)=A^{(0)}+\Delta_{A}$ and
$b(\xi)=b^{(0)}+\Delta_{b}$ are affine functions of $\xi$. Moreover,
Assumption 1 can be easily satisfied by choosing a sufficiently large
$\rho_{u}$, which leads to a performance certificate according to Theorem
III.1. The following result shows that (17) can be reformulated into a
tractable conic quadratic problem.

######  Proposition IV.1.

A vector $g^{\star}$ is a minimizer of (17) if and only if $g^{\star}$ also
minimizes the following conic quadratic problem

$\begin{array}[]{l}\mathop{{\rm{min}}}\limits_{g\in\mathcal{G}}\;\;\|A^{(0)}g-b^{(0)}\|+\rho_{u}\sqrt{\|g\|^{2}+1}.\end{array}$
(18)

Moreover, the minima of (17) and (18) coincide.

###### Proof.

The proof is adapted from [35, Theorem 3.1]. Consider a fixed $g$ in (17). It
follows from triangle inequality that

$\displaystyle\max_{(\Delta_{A},\Delta_{b})\in\mathcal{D}_{\rm uns}}\ \
\left\|(A^{(0)}+\Delta_{A})g-(b^{(0)}+\Delta_{b})\right\|$ $\displaystyle\leq$
$\displaystyle\left\|A^{(0)}g-b^{(0)}\right\|+\max_{(\Delta_{A},\Delta_{b})\in\mathcal{D}_{\rm
uns}}\left\|\Delta_{A}g-\Delta_{b}\right\|$ $\displaystyle=$
$\displaystyle\left\|A^{(0)}g-b^{(0)}\right\|+\rho_{u}\sqrt{\|g\|^{2}+1}.$

Choose $\hat{\Delta}=[\hat{\Delta}_{A}\ \hat{\Delta}_{b}]\in\mathcal{D}_{\rm
uns}$ such that

$[\hat{\Delta}_{A}\
\hat{\Delta}_{b}]=\frac{\rho_{u}\omega}{\sqrt{\|g\|^{2}+1}}[g^{\top}\
1],\quad\text{where}$
$\omega=\begin{cases}\frac{A^{(0)}g-b^{(0)}}{{\|A^{(0)}g-b^{(0)}\|}}&\text{if
}A^{(0)}g\neq b^{(0)}\\\ \text{any unit-norm
vector}&\text{otherwise}.\end{cases}$

Since $\hat{\Delta}$ is a rank-one matrix, we have
$\|\hat{\Delta}\|_{F}=\|\hat{\Delta}\|=\rho_{u}$, and

$\displaystyle\left\|(A^{(0)}+\hat{\Delta}_{A})g-(b^{(0)}+\hat{\Delta}_{b})\right\|$
$\displaystyle=$
$\displaystyle\left\|A^{(0)}g-b^{(0)}\right\|+\left\|\hat{\Delta}_{A}g-\hat{\Delta}_{b}\right\|$
$\displaystyle=$
$\displaystyle\left\|A^{(0)}g-b^{(0)}\right\|+\rho_{u}\sqrt{\|g\|^{2}+1}$
$\displaystyle\leq$
$\displaystyle\max_{(\Delta_{A},\Delta_{b})\in\mathcal{D}_{\rm uns}}\ \
\left\|(A^{(0)}+\Delta_{A})g-(b^{(0)}+\Delta_{b})\right\|,$

where the first equality follows from triangle inequality and the fact that
$(A^{(0)}g-b^{(0)})$ is a nonnegative scalar of
$(\hat{\Delta}_{A}g-\hat{\Delta}_{b})$. Since the lower bound coincides with
the upper bound on the inner maximization problem (17) for any fixed $g$, the
claim follows. ∎

Similarly, if the vector $b$ in (17) is not affected by uncertainty, that is,
$\Delta_{b}=0$, then the corresponding tractable reformulation becomes [36,
Theorem 2]

$\begin{array}[]{l}\mathop{{\rm{min}}}\limits_{g\in\mathcal{G}}\;\;\|A^{(0)}g-b^{(0)}\|+\rho_{u}\|g\|.\end{array}$

The min-max formulation (17) is equivalent to the regularization (18), which
reminiscent of other DeePC regularizations studied in [10, 16, 18, 29]. We
will explore these connections in Section V.

We remark that the considered uncertainties in (17) can be overly conservative
for several reasons: 1) the constructed worst-case realization of the
uncertainties satisfies rank$([\Delta_{A}\;\Delta_{b}])=1$ [35], which cannot
be achieved in general when Hankel matrices are used as predictors; 2) Hankel
structures cannot be imposed on $\Delta_{A}$; 3) the uncertainties are
indirectly added to the input/output data through the scaling matrices, e.g.,
$Q^{\frac{1}{2}}$ and $\lambda_{y}^{\frac{1}{2}}I$; 4) the columns in
$[\Delta_{A}\;\Delta_{b}]$ are uncorrelated when Page matrices or trajectory
matrices are used, but $\mathcal{D}_{\rm uns}$ considers the Frobenius norm of
$[\Delta_{A}\;\Delta_{b}]$ to be bounded which may impose correlations among
the columns. In the following sections we alleviate these sources of
conservatism. First we consider column-wise uncertainties for
$[\Delta_{A}\;\Delta_{b}]$.

###  IV-B Generalized Column-wise Uncertainties

Consider the following variant of (16) with generalized column-wise
uncertainties, that is,

$\displaystyle\mathop{{\rm{min}}}\limits_{g\in\mathcal{G}}\max_{(\Delta_{A},\Delta_{b})\in\mathcal{D}_{\rm
gcol}}\ \ \left\|(A^{(0)}+\Delta_{A})g-(b^{(0)}+\Delta_{b})\right\|.$ (19)

Here the generalized column-wise uncertainty set [37, Section 3] is defined
column-by-column as

$\begin{split}\mathcal{D}_{\rm
gcol}=\\{&(\Delta_{A},\Delta_{b})\mid\exists\rho_{A}\in\mathbb{R}^{H_{c}}:\|(\Delta_{A})_{i}\|\leq(\rho_{A})_{i},\forall
i\in[H_{c}],\\\ &f_{j}(\rho_{A})\leq 0,\ \forall
j\in[J],\Delta_{b}\leq\rho_{b}\\},\end{split}$

where the function
$f_{j}:\mathbb{R}^{H_{c}}\rightarrow\mathbb{[}-\infty,+\infty]$ is a proper,
closed and convex function for each $j\in[J]$, and $(\Delta_{A})_{i}$ is the
$i$-th column of the matrix $\Delta_{A}$.

######  Proposition IV.2.

If $\mathcal{D}_{\rm gcol}$ is nonempty and admits a Slater point, then a
vector $g^{\star}$ is a minimizer of (19) if and only if there exists a
$(\lambda^{\star},\\{y_{j}^{\star}\\}_{j=1}^{J})$ such that $g^{\star}$ also
minimizes the following convex optimization problem

$\begin{array}[]{cll}\displaystyle\min_{g\in\mathcal{G},\lambda\atop\\{y_{j}\\}_{j=1}^{J}}&\displaystyle\|A^{(0)}g-b^{(0)}\|+\sum_{j\in[J]}\psi_{j}(y_{j},\lambda_{j})+\rho_{b}\\\
{\rm s.t.}&\lambda_{j}\geq 0\quad\forall j\in[J]\\\
&\displaystyle\sum_{j\in[J]}y_{j}\leq-|g|,\end{array}$ (20)

where $\psi:\mathbb{R}^{n_{x}}\times\mathbb{R}_{\geq
0}\rightarrow[-\infty,+\infty]$ denotes the perspective function111We adopt
the definition of the perspective function from [38], where
$\psi(x,t)=tf^{*}(x/t)$ if $t>0$, and $\psi(x,t)=\delta^{*}_{{\rm dom}(f)}(x)$
if $t=0$, and $\delta^{*}_{\mathcal{S}}$ denotes the support function over the
set $\mathcal{S}$. The perspective function $\psi$ is proper, convex and
closed if and only if $f^{*}$ is proper, convex and closed. of the conjugate
function $f^{*}:\mathbb{R}^{n_{x}}\rightarrow[-\infty,+\infty]$ of $f$.
Moreover, the minima of (19) and (20) coincide.

###### Proof.

The proof is adapted from [37, Theorem 1]. Fix $g$ in (19). We first prove
that

$\displaystyle\max_{(\Delta_{A},\Delta_{b})\in\mathcal{D}_{\rm gcol}}\ \
\left\|(A^{(0)}+\Delta_{A})g-(b^{(0)}+\Delta_{b})\right\|$ $\displaystyle=$
$\displaystyle\|A^{(0)}g-b^{(0)}\|+\displaystyle\max_{\rho_{A}\in\mathcal{V}}\
\rho_{A}^{\top}|g|+\rho_{b},$

where $\mathcal{V}=\\{\rho_{A}\in\mathbb{R}^{H_{c}}_{>0}\mid\
f_{j}(\rho_{A})\leq 0,\ \forall j\in[J]\\}$. Note that

$\displaystyle\max_{(\Delta_{A},\Delta_{b})\in\mathcal{D}_{\rm gcol}}\ \
\|(A^{(0)}+\Delta_{A})g-(b^{(0)}+\Delta_{b})\|$ $\displaystyle=$
$\displaystyle\max_{(\Delta_{A},\Delta_{b})\in\mathcal{D}_{\rm gcol}}\ \
\|A^{(0)}g-b^{(0)}+\sum_{i\in[H_{r}]}(\Delta_{A})_{i}g_{i}-\Delta_{b}\|$
$\displaystyle\leq$
$\displaystyle\|A^{(0)}g-b^{(0)}\|+\max_{(\Delta_{A},\Delta_{b})\in\mathcal{D}_{\rm
gcol}}\sum_{i\in[H_{r}]}\|(\Delta_{A})_{i}g_{i}\|+\|\Delta_{b}\|$
$\displaystyle\leq$
$\displaystyle\|A^{(0)}g-b^{(0)}\|+\max_{\rho_{A}\in\mathcal{V}}\
\rho_{A}^{\top}|g|+\rho_{b}.$

To prove the converse inequality, choose a $\hat{\Delta}=[\hat{\Delta}_{A}\
\hat{\Delta}_{b}]\in\mathcal{D}_{\rm gcol}$ such that

$(\hat{\Delta}_{A})_{i}=-(\rho_{A})_{i}{\rm sgn}(g_{i})\omega\ {\rm and}\
\hat{\Delta}_{b}=\rho_{b}\omega,\quad\text{where}$
$\omega=\begin{cases}\frac{A^{(0)}g-b^{(0)}}{{\|A^{(0)}g-b^{(0)}\|}}&\text{if}A^{(0)}g\neq
b^{(0)}\\\ \text{any unit-norm vector}&\text{otherwise}.\end{cases}$

Then, we have

$\displaystyle\max_{(\Delta_{A},\Delta_{b})\in\mathcal{D}_{\rm gcol}}\ \
\|(b^{(0)}+\Delta_{b})-(A^{(0)}+\Delta_{A})g\|$ $\displaystyle\geq$
$\displaystyle\max_{\rho_{A}\in\mathcal{V}}\|b^{(0)}+\hat{\Delta}_{b}-A^{(0)}g-\sum_{i\in[H_{r}]}(\hat{\Delta}_{A})_{i}g_{i}\|$
$\displaystyle=$
$\displaystyle\max_{\rho_{A}\in\mathcal{V}}\|b^{(0)}+\rho_{b}\omega-A^{(0)}g+\sum_{i\in[H_{r}]}g_{i}(\rho_{A})_{i}{\rm
sgn}(g_{i})\omega\|$ $\displaystyle=$
$\displaystyle\max_{\rho_{A}\in\mathcal{V}}\|b^{(0)}+\rho_{b}\omega-A^{(0)}g+\sum_{i\in[H_{r}]}\rho_{A}^{\top}|g|\omega\|$
$\displaystyle=$
$\displaystyle\|A^{(0)}g-b^{(0)}\|+\displaystyle\max_{\rho_{A}\in\mathcal{V}}\
\rho_{A}^{\top}|g|+\rho_{b}.$

Finally, one can reformulate the inner maximization into its dual and obtain
(20) thanks to the strong duality of convex programs, which applies because
$\mathcal{D}_{\rm gcol}$ is nonempty and admits a Slater point. ∎

We remark that if the column-wise uncertainty is assumed, one can
conservatively approximate (19) via (17) with a large enough $\mathcal{D}_{\rm
uns}$, e.g., if $\rho_{u}^{2}\geq\|\rho_{A}\|^{2}+\rho_{b}^{2}$ in
$\mathcal{D}_{\rm uns}$, which implies that the minimum of (19) is smaller or
equal than that of (17). Problem (19) is useful if unstructured Page matrices
or trajectory matrices are considered as predictors. Furthermore, different
bounds for different columns in the data matrices can be incorporated in (19),
which comes in handy when the data comes from a time-varying system or the
columns are from independent experiments (i.e., trajectory matrices are used).
We relegate this direction to future research.

The following result shows that for a special case of (19), where
$\rho_{A}\in\mathbb{R}^{H_{c}}$ is a constant vector in $\mathcal{D}_{\rm
gcol}$, then the objective function of (19) can be reformulated into a least-
square problem with a 1-norm regularization.

######  Corollary IV.3.

Consider $\mathcal{D}_{\rm gcol}=\mathcal{D}_{\rm col}$ in (19), where

$\mathcal{D}_{\rm
col}:=\\{(\Delta_{A},\Delta_{b})\mid\|(\Delta_{A})_{i}\|\leq(\rho_{A})_{i},\forall
i\in[H_{c}],\Delta_{b}\leq\rho_{b}\\},$ (21)

then $g^{\star}$ is a minimizer of (19) if and only if $g^{\star}$ minimizes

$\mathop{{\rm{min}}}\limits_{g\in\mathcal{G}}\;\;\|A^{(0)}g-b^{(0)}\|+\rho_{A}^{\top}|g|+\rho_{b}\,.$
(22)

Moreover, the minima of (19) and (22) coincide.

###### Proof.

The claimed result follows from Proposition IV.2. ∎

Since $(A^{(0)})_{i}$ is a weighted input/output trajectory (with scaling
matrices $\lambda_{y}^{\frac{1}{2}}I$, $Q^{\frac{1}{2}}$, etc.), it could be
conservative to consider a column-wise uncertainty vector $(\Delta_{A})_{i}$
(with 2-norm bound) on $(A^{(0)})_{i}$ as the uncertainties are indirectly
added to the input/output data through the scaling matrices. We will next
demonstrate how to use interval uncertainties to incorporate the effects of
the scaling matrices.

###  IV-C Interval Uncertainties

Consider the following special case of (16)

$\displaystyle\mathop{{\rm{min}}}\limits_{g\in\mathcal{G}}\max_{(\Delta_{A},\Delta_{b})\in\mathcal{D}_{\rm
int}}\ \ \left\|(A^{(0)}+\Delta_{A})g-(b^{(0)}+\Delta_{b})\right\|$ (23)

with the nonempty interval uncertainty set

$\mathcal{D}_{\rm
int}=\\{(\Delta_{A},\Delta_{b})\mid|(\Delta_{A})_{ij}|\leq\bar{A}_{ij},\forall(i,j),|(\Delta_{b})_{i}|\leq\bar{b}_{i},\forall
i\\},$

where $\bar{A}\in\mathbb{R}^{H_{r}\times H_{c}}_{\geq 0}$,
$\bar{b}\in\mathbb{R}^{H_{r}}_{\geq 0}$, and $(\Delta_{A})_{ij}$ is the
element in the $i$-th row and $j$-th column of $\Delta_{A}$. The element of
$\bar{A}$ and $\bar{b}$ can be set to $0$ if there is no uncertainty in the
corresponding element of the data matrices. The corresponding $\mathcal{G}$ in
(23) can be represented as a polyhedron. Note that if the interval uncertainty
is assumed, one can still (conservatively) capture it with a large enough
$\mathcal{D}_{\rm col}$, e.g., $(\rho_{A})_{i}\geq\|(\bar{A})_{i}\|\ \forall
i\in[H_{c}]$ and $\rho_{b}\geq\|\bar{b}\|$ in $\mathcal{D}_{\rm col}$. In this
case, the minimum of (23) is less or equal than that of (19), suggesting a
less conservative result (better performance in the worst-case scenario) by
applying $\mathcal{D}_{\rm int}$. The following result shows that (23) can be
reformulated as a convex quadratic problem.

######  Proposition IV.4.

A vector $g^{\star}$ is a minimizer of (23) if and only if there exists a
$(\gamma^{\star},\nu^{\star})$ such that
$(g^{\star},\gamma^{\star},\nu^{\star})$ also minimizes the following convex
quadratic problem

$\begin{array}[]{cl}\displaystyle\min_{g\in\mathcal{G},\gamma,\nu}&\;\;\|\gamma+\bar{b}+\bar{A}\nu\|^{2}\\\
{\rm s.t.}&\;\;\displaystyle-\gamma\leq A^{(0)}g-b^{(0)}\leq\gamma\\\
&\;\;-\nu\leq g\leq\nu.\end{array}$ (24)

Moreover, the minimum of (23) coincides with
$\|\gamma^{\star}+\bar{b}+\bar{A}\nu^{\star}\|$.

###### Proof.

This proof is adapted from [34, Section 6.2]. One may observe that (23) can be
equivalently reformulated as

$\begin{array}[]{cl}\displaystyle\min_{g\in\mathcal{G},\tau}&\tau^{\top}\tau\\\
{\rm
s.t.}&\displaystyle|(A^{(0)}g-b^{(0)})_{i}|+\bar{b}_{i}+\sum_{j=1}^{H_{c}}\bar{A}_{ij}|g_{j}|\leq\tau_{i},\
\forall i\in[H_{r}],\end{array}$

where $(x)_{i}$ denotes the $i$-th element of the vector $x$, and the minimum
of (23) coincides with $\|\tau^{\star}\|$. One can now replace $\tau$ by the
introduced auxiliary variables $\gamma_{i}$ and $\nu_{j}$ for
$|(A^{(0)}g-b^{(0)})_{i}|$ and $|g_{j}|$, respectively, for every
$i\in[H_{r}]$ and $j\in[H_{c}]$. ∎

The interval uncertainty set allows us to consider independent upper and lower
bounds for the possible deviation in each entry of $A^{(0)}$ and $b^{(0)}$.
Hence, one can directly consider uncertainties on the input/output data and
cancel the effects of the scaling matrices. For example, one may consider the
measurement error of the $i{\rm th}$ output to be bounded by
$\tilde{y}_{i}\in\mathbb{R}_{\geq 0}$, and the process noise of the $i{\rm
th}$ input to be bounded by $\tilde{u}_{i}\in\mathbb{R}_{\geq 0}$ (which can
also be zero if no uncertainty occurs). In this case, we have

$\bar{A}=\begin{bmatrix}\lambda_{u}^{\frac{1}{2}}{\bf 1}_{T_{\rm ini}\times
H_{c}}\otimes\tilde{u}\\\ \lambda_{y}^{\frac{1}{2}}{\bf 1}_{T_{\rm ini}\times
H_{c}}\otimes\tilde{y}\\\ R^{\frac{1}{2}}{\bf 1}_{N\times
H_{c}}\otimes\tilde{u}\\\ Q^{\frac{1}{2}}{\bf 1}_{N\times
H_{c}}\otimes\tilde{y}\end{bmatrix}\quad{\rm
and}\quad\bar{b}=\begin{bmatrix}\lambda_{u}^{\frac{1}{2}}{\bf 1}_{T_{\rm
ini}}\otimes\tilde{u}\\\ \lambda_{y}^{\frac{1}{2}}{\bf 1}_{T_{\rm
ini}}\otimes\tilde{y}\\\ 0\\\ 0\end{bmatrix}\,.$ (25)

In the above setting, the possible deviation in each entry of $A^{(0)}$ and
$b^{(0)}$ is assumed to reside within an independent interval, and thus Hankel
structures for uncertainties in $A(\xi)$ cannot be enforced. For this reason,
in contrast to the case in which Page matrices or trajectory matrices are
considered, if Hankel matrices are used in (24) as predictors, then the
obtained solution could be overly conservative, similar to the cases with
unstructured uncertainties or column-wise uncertainties. To this end, in what
follows we illustrate how the uncertainties with Hankel structures can be
taken into account in the robust DeePC algorithm with a tractable formulation.

###  IV-D Structured Uncertainties

Consider (16) with the uncertainties residing within the following uncertainty
set

${\mathcal{D}}=\mathcal{D}_{\rm struct}:=\\{\xi\in\mathbb{R}^{n_{\xi}}\ |\
\|\xi\|\leq\rho_{s}\\}.$

In this case, the corresponding $\mathcal{G}$ is a convex set that consists of
linear and conic quadratic constraints. Recall that the vector of
uncertainties $\xi$ enters the elements of $A(\xi)$ and $b(\xi)$ with affine
structures, and we will show below how these affine structures allow us to
construct Hankel structures for the uncertainties.

We start by noting that (16) can be reformulated as

$\displaystyle\mathop{{\rm{min}}}\limits_{g\in\mathcal{G}}\max_{\xi\in\mathcal{D}}\
\ \left\|D(g)\xi-c(g)\right\|^{2},$ (26)

where the elements of $D(g)$ and $c(g)$ that are affine functions of $g$,
namely,

$D(g)=D^{(0)}+\sum_{\ell=1}^{H_{c}}D^{(\ell)}g_{\ell},\ \text{and}\
c(g)=c_{0}+\sum_{\ell=1}^{H_{c}}c_{\ell}g_{\ell}.$

The equivalence between (16) and (26) directly follows from

$\displaystyle(A^{(0)})_{\ell}=$ $\displaystyle c_{\ell}\
\forall\ell\in[H_{c}],\ (D^{(0)})_{j}=b_{j}\ \forall j\in[n_{\xi}],\
b_{0}=c_{0},$ $\displaystyle\quad(A^{(j)})_{\ell}=$
$\displaystyle(D^{(\ell)})_{j}\ \forall\ell\in[H_{c}]\ \forall j\in[n_{\xi}],$

where, for instance, $(A^{(0)})_{\ell}$ is the $\ell$-th column of $A^{(0)}$.

######  Example 1 (Additive Uncertainties with Hankel Structures in (26)).

Consider uncertainties directly on the time series recorded input data and
measured output data (i.e., $\hat{u}^{\rm d}$, $\hat{u}_{\rm ini}$,
$\hat{y}^{\rm d}$, and $\hat{y}_{\rm ini}$). Let $\xi={\rm
col}(\xi_{1}\in\mathbb{R}^{mT},\xi_{2}\in\mathbb{R}^{pT},\xi_{3}\in\mathbb{R}^{mT_{\rm
ini}},\xi_{4}\in\mathbb{R}^{pT_{\rm ini}})\in\mathcal{D}_{\rm struct}$ such
that

$\begin{bmatrix}U_{\rm P}(\xi)\\\ U_{\rm
F}(\xi)\end{bmatrix}=\mathscr{H}_{T_{\rm ini}+N}(\hat{u}^{\rm
d}+\alpha_{1}\xi_{1})=\begin{bmatrix}\hat{U}_{\rm P}\\\ \hat{U}_{\rm
F}\end{bmatrix}+\alpha_{1}\mathscr{H}_{T_{\rm ini}+N}(\xi_{1})\,,$
$\begin{bmatrix}Y_{\rm P}(\xi)\\\ Y_{\rm
F}(\xi)\end{bmatrix}=\mathscr{H}_{T_{\rm ini}+N}(\hat{y}^{\rm
d}+\alpha_{2}\xi_{2})=\begin{bmatrix}\hat{Y}_{\rm P}\\\ \hat{Y}_{\rm
F}\end{bmatrix}+\alpha_{2}\mathscr{H}_{T_{\rm ini}+N}(\xi_{2})\,,$ $u_{\rm
ini}(\xi)=\hat{u}_{\rm ini}+\alpha_{3}\xi_{3}\,,y_{\rm ini}(\xi)=\hat{y}_{\rm
ini}+\alpha_{4}\xi_{4}\,,$

where $\alpha_{1}$, $\alpha_{2}$, $\alpha_{3}$, and $\alpha_{4}$ are the
scaling coefficients. Note that $\alpha_{1}=\alpha_{3}=0$ in the absence of
input disturbance.

We now show that (26) with additive uncertainties in the input/output data
matrices (that obey the prescribed Hankel structures) admit a compact
representation. Notice that

$\begin{split}\begin{bmatrix}U_{\rm P}(\xi)\\\ U_{\rm
F}(\xi)\end{bmatrix}g=&\begin{bmatrix}\hat{U}_{\rm P}\\\ \hat{U}_{\rm
F}\end{bmatrix}g+\alpha_{1}\mathscr{H}_{T_{\rm ini}+N}(\xi_{1})g\\\
=&\begin{bmatrix}\hat{U}_{\rm P}\\\ \hat{U}_{\rm
F}\end{bmatrix}g+\alpha_{1}\mathscr{M}_{T_{\rm ini}+N}(g)\otimes
I_{m}\xi_{1}\,,\end{split}$ $\begin{split}\begin{bmatrix}Y_{\rm P}(\xi)\\\
Y_{\rm F}(\xi)\end{bmatrix}g=&\begin{bmatrix}\hat{Y}_{\rm P}\\\ \hat{Y}_{\rm
F}\end{bmatrix}g+\alpha_{2}\mathscr{H}_{T_{\rm ini}+N}(\xi_{2})g\\\
=&\begin{bmatrix}\hat{Y}_{\rm P}\\\ \hat{Y}_{\rm
F}\end{bmatrix}g+\alpha_{2}\mathscr{M}_{T_{\rm ini}+N}(g)\otimes
I_{p}\xi_{2}\,,\end{split}$

where $\mathscr{M}_{L}(x)\in\mathbb{R}^{L\times(n_{x}+L-1)}$
($x\in\mathbb{R}^{n_{x}}$) is defined as

$\mathscr{M}_{L}(x)=\begin{bmatrix}x_{1}&x_{2}&\cdots&x_{L}&\cdots&x_{n_{x}}&&&\\\
&x_{1}&x_{2}&\cdots&x_{L}&\cdots&x_{n_{x}}&&\\\
&&\ddots&\ddots&&\ddots&&\ddots&\\\
&&&x_{1}&x_{2}&\cdots&x_{L}&\cdots&x_{n_{x}}\end{bmatrix}.$

One can further partition $\mathscr{M}_{T_{\rm ini}+N}(g)$ into

$\begin{bmatrix}\mathscr{M}_{T_{\rm ini}+N}^{\rm P}(g)\in\mathbb{R}^{T_{\rm
ini}\times T}\\\\[6.0pt] \mathscr{M}_{T_{\rm ini}+N}^{\rm
F}(g)\in\mathbb{R}^{N\times T}\end{bmatrix}:=\mathscr{M}_{T_{\rm
ini}+N}(g)\,.$

Since

$A(\xi)=\begin{bmatrix}\lambda_{u}^{\frac{1}{2}}U_{\rm P}(\xi)\\\
\lambda_{y}^{\frac{1}{2}}Y_{\rm P}(\xi)\\\ R^{\frac{1}{2}}U_{\rm F}(\xi)\\\
Q^{\frac{1}{2}}Y_{\rm F}(\xi)\end{bmatrix}\quad{\text{a}nd}\quad
b(\xi)=\begin{bmatrix}\lambda_{u}^{\frac{1}{2}}u_{\rm ini}(\xi)\\\
\lambda_{y}^{\frac{1}{2}}y_{\rm ini}(\xi)\\\ 0\\\
Q^{\frac{1}{2}}r\end{bmatrix}\,$

by definition, we then have

$\begin{split}D(g)=&\begin{bmatrix}D_{1}(g)&D_{2}(g)\end{bmatrix}\,,\\\
c(g)=&\begin{bmatrix}\lambda_{u}^{\frac{1}{2}}\hat{U}_{\rm P}\\\
\lambda_{y}^{\frac{1}{2}}\hat{Y}_{\rm P}\\\ R^{\frac{1}{2}}\hat{U}_{\rm F}\\\
Q^{\frac{1}{2}}\hat{Y}_{\rm
F}\end{bmatrix}g-\begin{bmatrix}\lambda_{u}^{\frac{1}{2}}\hat{u}_{\rm ini}\\\
\lambda_{y}^{\frac{1}{2}}\hat{y}_{\rm ini}\\\ 0\\\
Q^{\frac{1}{2}}r\end{bmatrix}\,,\end{split}$ (27)

where

$\displaystyle
D_{1}(g)=\begin{bmatrix}\lambda_{u}^{\frac{1}{2}}\alpha_{1}\mathscr{M}^{\rm
P}_{T_{\rm ini}+N}(g)\otimes I_{m}&0\\\
0&\lambda_{y}^{\frac{1}{2}}\alpha_{2}\mathscr{M}^{\rm P}_{T_{\rm
ini}+N}(g)\otimes I_{p}\\\ R^{\frac{1}{2}}\alpha_{1}\mathscr{M}^{\rm
F}_{T_{\rm ini}+N}(g)\otimes I_{m}&0\\\
0&Q^{\frac{1}{2}}\alpha_{2}\mathscr{M}^{\rm F}_{T_{\rm ini}+N}(g)\otimes
I_{p}\end{bmatrix}$ $\displaystyle\text{and}\quad
D_{2}(g)=\begin{bmatrix}\lambda_{u}^{\frac{1}{2}}\alpha_{3}I_{mT_{\rm
ini}}&0\\\ 0&\lambda_{y}^{\frac{1}{2}}\alpha_{4}I_{pT_{\rm ini}}\\\ 0&0\\\
0&0\end{bmatrix}\,.$

Example 1 shows that the formulation in (26) admits a compact representation
of uncertainties with Hankel structures, as derived in (27). The following
result explicitly shows how (26) with $\mathcal{D}=\mathcal{D}_{\rm struct}$
can be solved by considering a tractable semi-definite programming problem.

######  Proposition IV.5.

A vector $g^{\star}$ is a minimizer of (26) with $\mathcal{D}=\mathcal{D}_{\rm
struct}$ if and only if there exists a $(\tau^{\star},\lambda^{\star})$ such
that $g^{\star}$ also minimizes the following semi-definite programming
problem

$\begin{array}[]{cll}\displaystyle\min_{g\in\mathcal{G},\tau,\lambda}&\;\;\tau\\\
{\rm s.t.}&\;\;\begin{bmatrix}\tau-\lambda\rho_{s}^{2}&0&c(g)^{\top}\\\
0&\lambda I&D(g)^{\top}\\\ c(g)&D(g)&I\end{bmatrix}\succeq 0.\end{array}$ (28)

Moreover, the minima of (26) and (28) coincide.

###### Proof.

The proof is adapted from [35]. For any fixed $g$, we have

$\displaystyle\tau\geq\max_{\xi\in\mathcal{D}}\left\|D(g)\xi-c(g)\right\|^{2}$
$\displaystyle\iff\quad$
$\displaystyle\tau\geq\left\|D(g)\xi-c(g)\right\|^{2}\quad\forall\xi:\|\xi\|\leq\rho_{s}.$

This semi-infinite constraint is satisfied if and only if

$\exists\lambda\geq
0:\begin{bmatrix}\tau-\lambda\rho_{s}^{2}-c(g)^{\top}c(g)&-c(g)^{\top}D(g)\\\
-D(g)^{\top}c(g)&\lambda I-D(g)^{\top}D(g)\end{bmatrix}\succeq 0,$

thanks to S-lemma [39, §2.3], which applies because
$\rho_{s}\in\mathbb{R}_{>0}$. The inner maximization problem in (26) can be
reformulated as

$\begin{array}[]{cll}\displaystyle\min_{\tau,\lambda}&\;\;\tau\\\ {\rm
s.t.}&\;\;\begin{bmatrix}\tau-\lambda\rho_{s}^{2}-c(g)^{\top}c(g)&-c(g)^{\top}D(g)\\\
-D(g)^{\top}c(g)&\lambda I-D(g)^{\top}D(g)\end{bmatrix}\succeq 0.\end{array}$

It then follows from Schur complement [40, §A.5.5] that

$\displaystyle\begin{bmatrix}\tau-\lambda\rho_{s}^{2}&0\\\ 0&\lambda
I_{H_{r}}\end{bmatrix}-\begin{bmatrix}c(g)^{\top}\\\
D(g)^{\top}\end{bmatrix}I_{H_{r}}\begin{bmatrix}c(g)&D(g)\end{bmatrix}\succeq
0$

if and only if

$\begin{bmatrix}\tau-\lambda\rho_{s}^{2}&0&c(g)^{\top}\\\ 0&\lambda
I&D(g)^{\top}\\\ c(g)&D(g)&I_{H_{r}}\end{bmatrix}\succeq 0,$

and the claim follows. ∎

Note that S-lemma [39, §2.3] enabling the proof still applies if a general
quadratic uncertainty set $\mathcal{D}_{\rm
quad}=\\{\xi\in\mathbb{R}^{n_{\xi}}\ |\
\xi^{\top}V\xi+2v^{\top}\xi\leq\rho_{s}^{2}\\}$ is considered in Proposition
IV.5, where $V\in\mathbb{R}^{n_{\xi}\times n_{\xi}}$ is a symmetric matrix
(not necessarily positive semi-definite) and $v\in\mathbb{R}^{n_{\xi}}$.

We remark that the structured uncertainty set in Example 1 is tighter than the
other uncertainty sets introduced before in modeling the uncertainties on
Hankel matrices, and thus it leads to less conservative results when Hankel
matrices are used as predictors. However, it requires one to solve a semi-
definite program and generally needs more computational effort. As will be
discussed in the simulation section, the advantages of the structured
uncertainty set come at the cost of more computational effort, and we will
also show that though the unstructured uncertainty set is conservative, it can
lead to satisfactory performance. One can analogously also consider $\xi$ on
the input/output data to construct Pages matrices and derive $D(g)$ and
$c(g)$, though it may result in a very high-dimension semi-definite program,
as Page matrices require much longer trajectories.

It is worth noting that for a special case of (26), where $D(g)$ is
independent of $g$, the semi-definite programming problem in Proposition IV.5
can be reduced to a second-order cone program, as shown in the following
result. This scenario is useful if the Hankel matrix data has no noise (e.g.,
the data is generated from a prescribed true model in simulations) and only
the initial trajectory $(u_{\rm ini},y_{\rm ini})$ is corrupted by noise.

######  Corollary IV.6.

If $D(g)=D$ in (26), then $g^{\star}$ is a minimizer of (26) if and only if
there exists a $(\nu^{\star},\lambda^{\star})$ such that $g^{\star}$ also
minimizes the second-order cone program

$\begin{array}[]{cll}\displaystyle\min_{g\in\mathcal{G},\nu,\lambda}&\;\;c(g)^{\top}c(g)+{\bf
1}^{\top}\nu+\lambda\rho_{s}^{2}\\\ {\rm
s.t.}&\;\;\left\|\begin{matrix}(2c(g)^{\top}DS)_{\ell}\\\
\nu_{\ell}+\alpha_{\ell}-\lambda\end{matrix}\right\|\leq\nu_{\ell}-\alpha_{\ell}+\lambda&\quad\forall\ell\in[H_{c}],\end{array}$

where $S$ is a nonsingular matrix that simultaneously diagonalizes $D^{\top}D$
and $I_{H_{c}}\in\mathbb{R}^{H_{c}\times H_{c}}$, namely,

$S^{\top}D^{\top}DS=:{\rm diag}(\alpha)\ ,\quad S^{\top}S=I_{H_{c}},$

and $\alpha=[\alpha_{1},\cdots,\alpha_{H_{c}}]$.

###### Proof.

For any fixed $g$, it follows from S-lemma [39, §2.3] that the inner
maximization problem in (26) can be equivalently reformulated as

$\begin{array}[]{l}\displaystyle\min_{\tau,\lambda}\;\;\tau\\\ {\rm
s.t.}\;\begin{bmatrix}\tau-\lambda\rho_{s}^{2}-c(g)^{\top}c(g)&-c(g)^{\top}D\\\
-D^{\top}c(g)&\lambda I_{H_{c}}-D^{\top}D\end{bmatrix}\succeq 0.\end{array}$

Since there exists a nonsingular matrix $S$ that simultaneously diagonalizes
$D^{\top}D$ and $I_{H_{c}}\in\mathbb{R}^{H_{c}\times H_{c}}$, the linear
matrix inequality (LMI) constraint above is satisfied if and only if

$\begin{bmatrix}\tau-\lambda\rho_{s}^{2}-c(g)^{\top}c(g)&-c(g)^{\top}DS\\\
-(DS)^{\top}c(g)&\lambda I_{H_{c}}-{\rm diag}(\alpha)\end{bmatrix}\succeq 0.$
(29)

Note that such a matrix $S$ always exists. It then follows from Schur
complement [40, §A.5.5] that (29) can be reduced to a set of second-order cone
constraints. ∎

##  V Robustness Induced by Regularization

Proposition IV.1 and Corollary IV.3 show that when unstructured or column-wise
uncertainties are considered, the min-max formulation in (10) reduces to a
minimization problem with an additional regularization term on $g$ (see (18)
and (22)). The robustness induced by different geometry of uncertainties can
thus be interpreted as different types of regularization on $g$. This is
consistent with the observation that regularization is instrumental for
ensuring good performance when the system is subject to disturbances (see also
the discussion in Section I).

Notice that the 2-norm is considered in the objective functions of (18) and
(22), which differ from the existing literature where quadratic costs are
generally considered (see, e.g., [14, 16, 29]). For instance, in the spirit of
regularized DeePC, one may consider quadratic costs on the input/output
signals and a quadratic regularization on $g$ as

$\begin{array}[]{cl}\displaystyle\mathop{{\rm{min}}}\limits_{g\in\mathcal{G}}&\|\hat{U}_{\rm
F}g\|_{R}^{2}+\|{\hat{Y}_{\rm F}g-r}\|_{Q}^{2}+\lambda_{u}\|\hat{U}_{\rm
P}g-\hat{u}_{\rm ini}\|^{2}\\\ &+\lambda_{y}\|\hat{Y}_{\rm P}g-\hat{y}_{\rm
ini}\|^{2}+\lambda_{g}\|g\|^{2},\end{array}$ (30)

or a 1-norm regularization on $g$ as

$\begin{array}[]{cl}\displaystyle\mathop{{\rm{min}}}\limits_{g\in\mathcal{G}}&\|\hat{U}_{\rm
F}g\|_{R}^{2}+\|{\hat{Y}_{\rm F}g-r}\|_{Q}^{2}+\lambda_{u}\|\hat{U}_{\rm
P}g-\hat{u}_{\rm ini}\|^{2}\\\ &+\lambda_{y}\|\hat{Y}_{\rm P}g-\hat{y}_{\rm
ini}\|^{2}+\lambda_{g}\|g\|_{1},\end{array}$ (31)

where the set $\mathcal{G}$ is the same as that in (10) to further robustify
the input/output constraints. The following results show how (30) and (31) are
respectively related to the tractable formulations of robust DeePC (18) and
(22), and thus can be considered as special cases of (10).

######  Theorem V.1 (Robustness induced by quadratic regularization).

If $g^{\star}\neq 0$ is a minimizer of (30), then $g^{\star}$ minimizes (17)
(and equivalently, (18)) with

$\rho_{u}=\begin{cases}\frac{\lambda_{g}\sqrt{\left\|g^{\star}\right\|^{2}+1}}{\left\|A^{(0)}g^{\star}-b^{(0)}\right\|}&\text{if
$A^{(0)}g^{\star}\neq b^{(0)}$}\\\
\lambda_{g}\sqrt{\left\|g^{\star}\right\|^{2}+1}&\text{otherwise.}\end{cases}$
(32)

Moreover, $\rho_{u}$ in (32) is strictly monotonically increasing with
$\lambda_{g}$ chosen in (30).

###### Proof.

See Appendix A. ∎

Theorem V.1 shows that by incorporating unstructured uncertainties, robust
DeePC reduces to regularized DeePC with quadratic regularization where the
input/output constraints are also robustified. Note that Theorem V.1 extends
the results in [29] by further incorporating robust constraints. It highlights
the importance of the regularization from a robust optimization perspective,
namely, quadratic regularization of $g$ in (30) is equivalent to a robust
reformulation as in (17) with an implicit bound for the uncertainty set.
Hence, by choosing a sufficiently large $\lambda_{g}$ in (30), a sufficiently
large bound $\rho_{u}$ can be obtained in (17) such that Assumption 1 holds
and Theorem III.1 applies.

It has been observed in [14] and [10] that the 1-norm regularization promotes
sparsity, selecting the most informative noisy trajectories to predict the
future behaviour. The following result provides new interpretations of 1-norm
regularization in DeePC from a min-max optimization perspective.

######  Theorem V.2 (Robustness induced by 1-norm regularization).

If $g^{\star}\in\mathbb{R}^{H_{c}}$ is a minimizer of (31), then $g^{\star}$
minimizes (19) with $\mathcal{D}_{\rm gcol}=\\{(\Delta_{A},\Delta_{b})\ |\
\|(\Delta_{A})_{i}\|\leq\rho_{c},\forall
i\in[H_{c}],\Delta_{b}\leq\rho_{c}\\}$ (and equivalently, (22)), where

$\rho_{c}=\begin{cases}\frac{\lambda_{g}}{2\left\|A^{(0)}g^{\star}-b^{(0)}\right\|}&\text{if
$A^{(0)}g^{\star}\neq b^{(0)}$}\\\ \lambda_{g}/2&\text{otherwise.}\end{cases}$
(33)

Moreover, $\rho_{c}$ in (33) is strictly monotonically increasing with
$\lambda_{g}$ chosen in (31).

###### Proof.

See Appendix B. ∎

Theorem V.2 shows that the 1-norm regularization in (31) provides robustness
to column-wise uncertainties in the data matrices and the initial trajectory.
One can further conclude that a sufficiently large $\lambda_{g}$ in (39) leads
to a sufficiently large $\rho_{c}$ in (19) such that Assumption 1 is satisfied
and Theorem III.1 applies. The 1-norm regularization on $g$ may be more
appropriate than the quadratic regularization when Page matrices or trajectory
matrices are used as predictors, because the 1-norm regularization implies a
column-wise uncertainty set that assumes no correlation among different
columns. Moreover, (31) can be solved by quadratic programming since it can be
reformulated as (41) in the proof.

##  VI Simulation Results

Figure 2: Control scheme for a grid-connected converter.

In this section, we provide simulation results to illustrate the effectiveness
of the proposed robust DeePC algorithm.

###  VI-A Simulation Case Study

We consider a grid-connected three-phase power converter represented as in
Fig. 2, and apply robust DeePC to regulate the active and reactive power.
Conventionally, power regulation of grid-connected converters can be achieved
by PI controllers. However, the power grid is ever-changing and in general
unknown from the perspective of a converter, which significantly affects the
performance of (fixed) PI controllers and may even result in instabilities
[41]. As a remedy, we employ robust DeePC to perform model-free, robust, and
optimal power control for converters. Thanks to the data-centric
representation, input/output data of the converter can be collected to capture
the system dynamics, predict the future behaviors, and calculate optimal
control sequences, without assuming a specific model.

As shown in Fig. 2, we choose the active power $P_{E}=V_{d}I_{d}+V_{q}I_{q}$
and the reactive power $Q_{E}=V_{q}I_{d}-V_{d}I_{q}$ to be the output signals
of the converter system, and the robust DeePC algorithm provides optimal
control inputs for the current references $I_{d}^{\rm ref}$, $I_{q}^{\rm
ref}$. Note that $u_{i,t}$ is the $i{\rm th}$ element of $u_{t}$ and $y_{i,t}$
is the $i{\rm th}$ element of $y_{t}$ in Fig.2. The power regulation is
achieved by setting the reference vector in (10) to $r=I_{N}\otimes{\rm
col}(P_{0},Q_{0})$, where $P_{0}$ is the active power reference value and
$Q_{0}$ is the reactive power reference value.

Throughout our simulations, we use the nonlinear converter model; the same
trends in the results are, however, also observed if one uses a linearization.
We use the base values $f_{\rm b}=50{\rm Hz}$, $S_{\rm b}=1.5{\rm kW}$, and
$U_{\rm b}=280{\rm V}$ for per-unit calculations of the converter system. The
LCL parameters are: $L_{F}=0.05{\rm(p.u.)}$ (with resistance
$R_{F}=0.01{\rm(p.u.)}$), $L_{g}=0.05{\rm(p.u.)}$ (with resistance
$R_{g}=0.01{\rm(p.u.)}$), and capacitor $C_{F}=0.05{\rm(p.u.)}$. The local
load is $R_{\rm Load}=2{\rm(p.u.)}$. The PI parameters of the current control
are $\\{0.2,10\\}$. The sampling time for the robust DeePC algorithm is $1{\rm
ms}$. The parameters for the robust DeePC are: $T_{\rm ini}=5$, $N=25$,
$T=120$, $R=I_{2N}$, $Q=10^{5}I_{2N}$, $\lambda_{u}=\lambda_{y}=10^{5}$. The
control horizon is $k=N$. In what follows, we assume that the output data (in
the Hankel matrices and the initial trajectory) are corrupted by measurement
noise, while the input data are known exactly, and we show the realized
trajectory of the system outputs (without showing the noise). Before the
robust DeePC algorithm is activated, persistently exciting white noise signals
are injected into the system through $I_{d}^{\rm ref}$ and $I_{d}^{\rm ref}$
for $0.12{\rm s}$ to collect the input/output data. Fig. 3 plots the
input/output responses during this data-collection period, which implies that
the active power and the reactive power are perturbed, but within an
acceptable range.

Figure 3: Input and output responses during the data-collection period.

###  VI-B Comparison of Conservativeness

Figure 4: Comparison of robust DeePC with structured set and unstructured set
when the uncertainty samples are drawn from a structured set. Hankel matrices
are used as predictors.

In this subsection, we compare the performance and conservativeness of robust
DeePC when incorporating different uncertainty sets. We activate the robust
DeePC algorithm at $t=0{\rm s}$ with $P_{0}=0$ and $Q_{0}=0$ and change
$P_{0}$ from $0$ to $0.1{\rm(p.u.)}$ at $t=0.2{\rm s}$. We consider the
realized cost of applying the optimal control sequence obtained by solving
(10) at $t=0.2{\rm s}$ in open loop for the whole horizon $N$. Though the
realized cost is related to the optimization cost in Theorem III.1, it is
still unclear if a smaller (albeit tight) uncertainty set in robust DeePC
leads to a better average realized cost.

We start with a study where Hankel matrices are used as predictors.
Theoretical intuition from Section IV suggests that when Hankel matrices are
used, robust DeePC with structured uncertainties should be less conservative
than that with unstructured uncertainties. To test this intuition, we assume
that the output data $\hat{y}^{\rm d}$ and $\hat{y}_{\rm ini}$ are affected by
uncertainties residing in the structured set with $\rho_{s}=0.016$, and draw
1000 samples that are uniformly-distributed over this structured set. Fig. 4
plots the corresponding realized costs of applying robust DeePC with the
structured set (by solving (28) where $\rho_{s}=0.016$) and the unstructured
set (by solving (18) where $\rho_{u}=27.39$, the smallest value such that the
unstructured set contains the structured set), respectively. Fig. 4 shows that
robust DeePC with structured set performs better than that with unstructured
set, which confirms our intuition. Realized cost becomes even higher with the
column-wise and interval uncertainty sets when scaled to contain the
structured uncertainty set (data not shown).

Figure 5: Comparison of robust DeePC with interval set, column-wise set, and
unstructured set when the uncertainty samples are drawn from a interval set.
Page matrices are used as predictors.

The analysis in Section IV indicates that column-wise set and interval set are
more appropriate choices when Page matrices are used. To validate this, we
next consider Page matrices as predictors and compare the performance of
unstructured set, column-wise set, and interval set. We collect a longer
input/output trajectory such that the constructed Page matrices have the same
dimensions as the Hankel matrices. We first assume that the output data is
affected by uncertainties residing in an interval set with the upper and lower
bounds for the output signals being $\pm\rho_{i}$ ($\rho_{i}=0.001$). We
consider 1000 samples uniformly-distributed over this interval set and plot in
Fig. 5 the corresponding realized costs of applying robust DeePC with
respectively the same interval set (by solving (24)), the column-wise set (by
solving (22) where for all $k\in[H_{c}]$: $(\rho_{A})_{k}=\rho_{c}=2.45$), and
the unstructured set (by solving (18) where $\rho_{u}=23.49$); we again choose
the smallest $\rho_{c}$ and $\rho_{u}$ such that the considered interval set
is contained in the column-wise set and the unstructured set. Fig. 5 shows
that robust DeePC with the interval set achieves better average and worst-case
performance than the other two methods, because the interval set is tight in
this case. However, the cases with the column-wise set did not show superior
performance over the unstructured set, as the column-wise set is also
conservative under the above setting.

Figure 6: Comparison of robust DeePC with column-wise set and unstructured set
when the uncertainty samples are drawn from a column-wise set. Trajectory
matrices are used as predictors.

One advantage of incorporating column-wise uncertainties in robust DeePC is
that one can compactly consider different bounds for different columns in the
data matrices. To demonstrate this, we use trajectory matrices as predictors
and test the performance of robust DeePC when the columns in the trajectory
matrices are subject to different levels of uncertainties, which could be the
case when the columns are from different experiments. We consider 1000
uncertainty samples uniformly-distributed over the column-wise set (21) where
for all $k\in[H_{c}]$: $(\rho_{A})_{k}=0.63k$ and $\rho_{b}=0.63H_{c}$. Fig. 6
shows the performance of robust DeePC with respectively the considered column-
wise set (by solving (22)) and the unstructured set that contains the column-
wise set (by solving (18) where $\rho_{u}=324.73$). Under this setting, we
observe that robust DeePC with the column-wise set has a better average
performance than that with the unstructured set, thanks to the tighter
uncertainty representation.

The average runtime of these simulations (on an Intel Core i7 9750H CPU with
16GB RAM) is shown in Table III. As expected, it takes significantly more time
to solve the semi-definite program in (28).

TABLE III: Average runtime when applying different uncertainty sets. Uncertainty sets | Tractable formulations | Average runtime
---|---|---
Unstructured set | Eq. (18) | $0.0146{\rm s}$
Column-wise set | Eq. (22) | $0.0269{\rm s}$
Interval set | Eq. (24) with (25) | $0.1063{\rm s}$
Structured set | Eq. (28) with (27) | $4.3234{\rm s}$

###  VI-C Comparison of On-line Performance

In the previous section, we consider uncertainty samples that are uniformly-
distributed over some prescribed uncertainty set to test the conservativeness
of robust DeePC. However, in practice, measurement noise may be better
represented by white noise. Hence, for a more realistic comparison, we then
turn to band-limited Gaussian white noise [42] in the output measurements and
test the on-line performance of robust DeePC under two different noise levels:
1) noise power: $4\times 10^{-10}{\rm(p.u.)}$; 2) noise power: $1.6\times
10^{-7}{\rm(p.u.)}$. Two example sequences are given in Fig. 7 to illustrate
the range of uncertainties with different noise power.

In what follows, we consider Hankel matrices as predictors. After the robust
DeePC algorithm is activated at $t=0{\rm s}$ with $P_{0}=0$ and $Q_{0}=0$, we
change $P_{0}$ from $0$ to $0.1{\rm(p.u.)}$ at $t=0.2{\rm s}$ and back to $0$
at $t=0.4{\rm s}$. Fig. 8 (a) plots the active power responses of applying
robust DeePC with different uncertainty sets, where the noise power is
$4\times 10^{-10}{\rm(p.u.)}$. For the interval set, we consider the upper and
lower bounds of the output signals to be $\pm 0.001$ (compare Fig. 7 (a)),
then choose the smallest $\rho_{c}$ such that the interval set is contained in
the column-wise set; we choose the smallest $\rho_{s}$ such that the
structured set contains the interval uncertainty $[-0.001,0.001]$ in the
outputs; finally, we choose the smallest $\rho_{u}$ such that the unstructured
set contains the structured set. In all cases, the converter has satisfactory
tracking performance. Moreover, the performance is improved by reducing the
control horizon from $25$ to $5$, thanks to the faster feedback.

The above simulations with $k=25$ (i.e., the whole input sequence is applied,
in line with Theorem III.1) are repeated 100 times with different data sets to
construct the Hankel matrices and different random seeds to generate the
measurement noise. The histogram in Fig. 8 (b) shows the realized input/output
cost from $0$ to $0.6{\rm s}$ (i.e., $\sum_{i=0}^{600}\|u_{{\rm
sim},i}\|^{2}_{R}+\|y_{{\rm sim},i}-r_{{\rm sim},i}\|^{2}_{Q}$ where $u_{{\rm
sim},i}$, $y_{{\rm sim},i}$, and $r_{{\rm sim},i}$ are respectively the input,
output, and reference data at time $i$ obtained from the simulations). All the
cases have satisfactory performance, and the cases with the structured set and
the interval set have better average performance than the cases with the
unstructured set or the column-wise set, consistent with the analysis in
Section IV.

Figure 7: Two example sequences of white noise signals. (a) Noise power:
$4\times 10^{-10}{\rm(p.u.)}$. (b) noise power: $1.6\times
10^{-7}{\rm(p.u.)}$. Figure 8: Comparison of on-line performance with noise
power being $4\times 10^{-10}{\rm(p.u.)}$. (a) Time-domain responses. (b)
Realized costs. Figure 9: Comparison of on-line performance with noise power
being $4\times 1.6^{-7}{\rm(p.u.)}$. (a) Time-domain responses. (b) Realized
costs.

We next test the performance of robust DeePC when a higher measurement noise
level is considered. Fig. 9 plots the active power responses when the noise
power is $4\times 1.6^{-7}{\rm(p.u.)}$, with the bounds for the uncertainty
sets scaled accordingly. For instance, the upper and lower bounds of the
output signals are $\pm 0.02$ in the interval set (compare Fig. 7 (b)). All
the cases still have satisfactory tracking performances, but are worse than
those in Fig. 8 (a) as the impact of noise increases. The above simulations
with $k=25$ are repeated 100 times, and the histogram in Fig. 9 (b) compares
the performance of robust DeePC with different uncertainty sets. Surprisingly,
the cases with the unstructured set have the best average performance. This
counter-intuitive effect is due to the fact that the noise is not uniformly-
distributed over any of the uncertainty sets, and it suggests that robust
DeePC with unstructured set may have superior performance in practice.

###  VI-D Resilience in Presence of Outliers

To test the resilience of robust DeePC, we consider bad data points in the
output trajectory $\hat{y}^{\rm d}$ used to construct the Hankel matrices. The
output data is set as 0 in these bad data points, to reflect, for example, the
loss of a measurement. The histogram in Fig. 10 displays the realized
input/output cost from $0$ to $0.6{\rm s}$ in the presence of $5$ such bad
data points (appeared in the middle of the data sequence of $P_{E}$). The
robust DeePC with the structured set has the best average performance in this
case, which we attribute to the Hankel structure of the uncertainties. Robust
DeePC also leads to satisfactory performance when the interval set or the
column-wise set are considered. By comparison, the performance degrades
significantly when the unstructured set is used, suggesting that this may
result in unacceptable resilience. This problem can possibly be resolved by
setting a sufficiently large bound, but at the cost of being overly
conservative in normal situations.

Figure 10: Comparison of on-line performance with $5$ bad data points.

##  VII Conclusions

This paper proposed a robust DeePC framework to perform robust and optimal
model-free predictive control. The robust DeePC involves solving a min-max
optimization problem to robustify the optimal control sequence against
uncertainties in the input/output data that used for predictions. We showed
that by applying robust DeePC, the realized input/output cost can be bounded
if the considered uncertainty set captures the perfect input/output data. We
explicitly derived tractable formulations for robust DeePC when different
geometries of uncertainty sets are incorporated. In particular, when Hankel
matrices are used as predictors, we illustrated how uncertainties with Hankel
structures can be taken into account in a structured uncertainty set to reduce
the conservativeness. By incorporating appropriate geometries of uncertainty
sets, the robust DeePC recovers existing regularized DeePC algorithms with
robustified constraints. In particular, a quadratic regularization corresponds
to considering an unstructured uncertainty set, while a 1-norm regularization
corresponds to a column-wise uncertainty set. The robust DeePC algorithm with
different uncertainty sets was tested by high-fidelity simulations on a grid-
connected converter system, which shows satisfactory performance even with
noisy measurements and bad data.

##  Appendix A Proof of Theorem V.1

We start by compactly rewriting (30) as

$\mathop{\min}\limits_{g\in\mathcal{G}}\;\|A^{(0)}g-b^{(0)}\|^{2}+\lambda_{g}\|g\|^{2}\,.$
(34)

Note that the set $\mathcal{G}$ admits a tractable reformulation as the semi-
infinite constraints can be reformulated into a finite number of convex
quadratic inequalities using standard robust optimization techniques, as shown
in Table I. Accordingly, we compactly rewrite $\mathcal{G}$ as
$\mathcal{G}=\\{g\ |\ Gg+\mathcal{N}_{V}(g)\leq q\\}$ where
$\mathcal{N}_{V}(g)=[\|V^{(1)}g\|\;\|V^{(2)}g\|\;\cdots\;\|V^{(n_{q})}g\|]^{\top}$
for some $G\in\mathbb{R}^{n_{q}\times H_{c}}$, $v,q\in\mathbb{R}^{n_{q}}$, and
$V^{(i)}\in\mathbb{R}^{n_{q}\times n_{q}}\;\forall i\in[n_{q}]$.

Consider the Lagrangian of (34)

$\mathcal{L}_{Q}(g,\mu_{Q})=\|A^{(0)}g-b^{(0)}\|^{2}+\lambda_{g}\|g\|^{2}+\mu_{Q}^{\top}(Gg+\mathcal{N}_{V}(g)-q),$
(35)

where $\mu_{Q}$ is the vector of the dual variables. Since $\mathcal{G}$ is
nonempty, there exist a solution $(g^{\star},\mu_{Q}^{\star})$ to the
Karush–Kuhn–Tucker (KKT) conditions of (34)

$\displaystyle\begin{array}[]{rl}2{(A^{(0)})^{\top}(A^{(0)}g-b^{(0)})}+2{\lambda_{g}g}+G^{\top}\mu_{Q}&\\\
+\sum_{i\in[n_{q}]}(\mu_{Q})_{i}\frac{(V^{(i)})^{\top}V^{(i)}g}{\|V^{(i)}g\|}&=0\,,\end{array}$
(36c) $\displaystyle\mu_{Q}^{\top}(Gg+\mathcal{N}_{V}(g)-q)=0\,,$ (36d)
$\displaystyle Gg+\mathcal{N}_{V}(g)\leq q\,,$ (36e) $\displaystyle\mu_{Q}\geq
0\,.$ (36f)

Therefore, the vector $g^{\star}$ is a minimizer of (34).

Following Proposition IV.1, the min-max problem (17) can be equivalently
reformulated as (18). Consider the Lagrangian of (18)

$\mathcal{L}(g,\mu){}={}\|A^{(0)}g-b^{(0)}\|+{\rho_{u}}\sqrt{\left\|g\right\|^{2}+1}+\mu^{\top}(Gg+\mathcal{N}_{V}(g)-q),$
(37)

where $\mu$ is the vector of the dual variables. By choosing $\rho_{u}$ in
(32), it can be verified that $(g^{\star},\mu^{\star},y^{\star})$, where

$(\mu^{\star},y^{\star})=\begin{cases}\begin{array}[]{r}\left(\frac{\mu_{Q}^{\star}}{2\|A^{(0)}g^{\star}-b^{(0)}\|},\frac{(A^{(0)})^{\top}(A^{(0)}g^{\star}-b^{(0)})}{\|A^{(0)}g^{\star}-b^{(0)}\|}\right)\\\
\text{if $A^{(0)}g^{\star}\neq b^{(0)}$,}\end{array}\\\
\left(\frac{\mu_{Q}^{\star}}{2},0\right)\quad\text{otherwise,}\end{cases}$

and $g^{\star}$ as before, satisfy the KKT conditions of (37):

$\displaystyle\begin{array}[]{ll}y+\frac{\rho_{u}g}{\sqrt{\|g\|^{2}+1}}+G^{\top}\mu\\\
\hskip
2.84526pt+\sum_{i\in[n_{q}]}\mu_{i}\frac{(V^{(i)})^{\top}V^{(i)}g}{\|V^{(i)}g\|}=0\,,\end{array}$
(38c)
$\displaystyle\begin{array}[]{r}\|A^{(0)}\bar{g}-b^{(0)}\|\geq\|A^{(0)}g-b^{(0)}\|+y^{\top}(\bar{g}-g)\\\
\forall\bar{g}\in\mathbb{R}^{H_{c}}\,,\end{array}$ (38f)
$\displaystyle\mu^{\top}(Gg+\mathcal{N}_{V}(g)-q)=0\,,$ (38g) $\displaystyle
Gg+\mathcal{N}_{V}(g)\leq q\,,$ (38h) $\displaystyle\mu\geq 0\,.$ (38i)

Thus, the vector $g^{\star}$ is also a minimizer of (18) (and (17)).

Next we prove the monotonic relationship between $\rho_{u}$ and $\lambda_{g}$.
Let $g_{1}$ be the minimizer of (34) with $\lambda_{g}=\lambda_{g1}>0$ (and
the minimizer of (18) with $\rho_{u}=\rho_{1}$), and let $g_{2}$ be the
minimizer of (34) with $\lambda_{g}=\lambda_{g2}>\lambda_{g1}$ (and the
minimizer of (18) with $\rho_{u}=\rho_{2}$). If $g_{1}\neq g_{2}$, according
to the definitions of $g_{1}$ and $g_{2}$, we have

$\|A^{(0)}g_{1}-b^{(0)}\|^{2}+\lambda_{g1}{\left\|g_{1}\right\|^{2}}-\lambda_{g1}{\left\|g_{2}\right\|^{2}}<\|A^{(0)}g_{2}-b^{(0)}\|^{2}\,,$
$\|A^{(0)}g_{2}-b^{(0)}\|^{2}<\|A^{(0)}g_{1}-b^{(0)}\|^{2}+\lambda_{g2}{\left\|g_{1}\right\|^{2}}-\lambda_{g2}{\left\|g_{2}\right\|^{2}}\,,$

leading to

$\lambda_{g1}({\left\|g_{1}\right\|^{2}}-{\left\|g_{2}\right\|^{2}})<\lambda_{g2}({\left\|g_{1}\right\|^{2}}-{\left\|g_{2}\right\|^{2}})\,,$

which indicates that ${\left\|g_{1}\right\|^{2}}>{\left\|g_{2}\right\|^{2}}$
because $0<\lambda_{g1}<\lambda_{g2}$. Then, we have
${\left\|g_{1}\right\|}>{\left\|g_{2}\right\|}$. Since $g_{1}$ minimizes (18)
with $\rho=\rho_{1}$ and $g_{2}$ minimizes (18) with $\rho=\rho_{2}$, we also
have

$\begin{split}&\|A^{(0)}g_{1}-b^{(0)}\|+\rho_{1}\sqrt{\left\|g_{1}\right\|^{2}+1}-\rho_{1}\sqrt{\left\|g_{2}\right\|^{2}+1}\\\
<&\|A^{(0)}g_{2}-b^{(0)}\|\\\
<&\|A^{(0)}g_{1}-b^{(0)}\|+\rho_{2}\sqrt{\left\|g_{1}\right\|^{2}+1}-\rho_{2}\sqrt{\left\|g_{2}\right\|^{2}+1}\,.\end{split}$

It can then be deduced that $\rho_{1}<\rho_{2}$. If $g_{1}=g_{2}$, we have
$\rho_{1}<\rho_{2}$ according to (32). Hence, $\rho_{u}$ is increasing with
the increase of $\lambda_{g}$ if $g^{\star}\neq 0$. This completes the proof.

##  Appendix B Proof of Theorem V.2

Notice that (31) can be compactly rewritten as

$\mathop{\min}\limits_{g\in\mathcal{G}}\;\|A^{(0)}g-b^{(0)}\|^{2}+\lambda_{g}\|g\|_{1}\,.$
(39)

When $(\rho_{A})_{i}=\rho_{b}=\rho_{c}\in\mathbb{R}_{>0}\ \forall i$,
Corollary IV.3 implies that $g^{\star}$ is a minimizer of (19) if and only if
$g^{\star}$ is a minimizer of the following conic quadratic optimization
problem

$\mathop{\min}\limits_{g\in\mathcal{G}}\;\|A^{(0)}g-b^{(0)}\|+\rho_{c}\|g\|_{1}\,.$
(40)

In what follows, we show that If $g^{\star}$ is a minimizer of (39), then
$g^{\star}$ also minimizes (40) by choosing $\rho_{c}$ according to (33). To
this end, first observe that both (39) and (40) have the same feasible region,
and they can be equivalently reformulated into the following compact forms

$\displaystyle\min_{\left(g,\nu\right)\in\mathcal{G}^{\prime}}\;\;\|A^{(0)}g-b^{(0)}\|^{2}+\lambda_{g}{\bf
1}_{H_{c}}^{\top}\nu,$ $\displaystyle{\rm and}$ (41)
$\displaystyle\min_{\left(g,\nu\right)\in\mathcal{G}^{\prime}}\;\;\|A^{(0)}g-b^{(0)}\|+\rho_{c}{\bf
1}_{H_{c}}^{\top}\nu,$ $\displaystyle{\rm respectively,}$ (42)

where $\mathcal{G}^{\prime}=\\{(g,\nu)\ |\ g\in\mathcal{G},\ -\nu\leq
g\leq\nu\\}$ is a nonempty polyhedron. Note that the set $\mathcal{G}$ can be
compactly rewritten as $\mathcal{G}=\\{g\ |\ G_{1}g+V_{1}|g|\leq q_{1}\\}$ for
some $G_{1},V_{1}\in\mathbb{R}^{n_{q1}\times H_{c}}$, and
$q_{1}\in\mathbb{R}^{n_{q1}}$. Then, one can replace $|g|$ by $\nu$ such that
the set $\mathcal{G}^{\prime}$ can be compactly rewritten as
$\mathcal{G}^{\prime}=\\{(g,\nu)\ |\ G_{\nu}{\rm col}(g,\nu)\leq q_{\nu}\\}$
for some $G_{\nu}\in\mathbb{R}^{n_{q_{\nu}}\times H_{c}}$ and
$q_{\nu}\in\mathbb{R}^{n_{q_{\nu}}}$.

Then, analogous to the proof of Theorem V.1, one can derive the KKT conditions
of (41) and (42), and then check that by choosing $\rho_{c}$ according to
(33), if $g^{\star}$ is a minimizer of (41), then it is also a minimizer of
(42).

It remains to prove the monotonic relationship between $\rho_{c}$ and
$\lambda_{g}$. Note that (39) and (40) have unique solutions thanks to [43,
Lemmas 3 and 4], which apply because $A^{(0)}$ contains noisy data so that
each element of $A^{(0)}$ can be seen as a random variable drawn from a
continuous distribution. Based on this fact, one can analogously prove the
monotonic relationship by referring to the proof of Theorem V.1.

## References

  * [1] L. Ljung, “System identification,” _Wiley Encyclopedia of Electrical and Electronics Engineering_ , pp. 1–19, 1999.
  * [2] M. Gevers, “Identification for control: From the early achievements to the revival of experiment design,” _European journal of control_ , vol. 11, no. 4-5, pp. 335–352, 2005.
  * [3] W. Favoreel, B. De Moor, and M. Gevers, “SPC: Subspace predictive control,” _IFAC Proceedings Volumes_ , vol. 32, no. 2, pp. 4004–4009, 1999.
  * [4] M. Morari and J. H. Lee, “Model predictive control: past, present and future,” _Computers & Chemical Engineering_, vol. 23, no. 4-5, pp. 667–682, 1999.
  * [5] G. Pillonetto, F. Dinuzzo, T. Chen, G. De Nicolao, and L. Ljung, “Kernel methods in system identification, machine learning and function estimation: A survey,” _Automatica_ , vol. 50, no. 3, pp. 657–682, 2014.
  * [6] R. Boczar, N. Matni, and B. Recht, “Finite-data performance guarantees for the output-feedback control of an unknown system,” in _2018 IEEE Conference on Decision and Control (CDC)_. IEEE, 2018, pp. 2994–2999.
  * [7] H. Hjalmarsson, M. Gevers, S. Gunnarsson, and O. Lequin, “Iterative feedback tuning: theory and applications,” _IEEE control systems magazine_ , vol. 18, no. 4, pp. 26–41, 1998.
  * [8] M. C. Campi, A. Lecchini, and S. M. Savaresi, “Virtual reference feedback tuning: a direct method for the design of feedback controllers,” _Automatica_ , vol. 38, no. 8, pp. 1337–1346, 2002.
  * [9] B. Recht, “A tour of reinforcement learning: The view from continuous control,” _Annual Review of Control, Robotics, and Autonomous Systems_ , vol. 2, pp. 253–279, 2019.
  * [10] F. Dörfler, J. Coulson, and I. Markovsky, “Bridging direct & indirect data-driven control formulations via regularizations and relaxations,” _arXiv preprint arXiv:2101.01273_ , 2021.
  * [11] J. C. Willems, P. Rapisarda, I. Markovsky, and B. L. De Moor, “A note on persistency of excitation,” _Systems & Control Letters_, vol. 54, no. 4, pp. 325–329, 2005.
  * [12] I. Markovsky and F. Dörfler, “Identifiability in the behavioral setting,” 2020, available online.
  * [13] H. J. van Waarde, C. De Persis, M. K. Camlibel, and P. Tesi, “Willems’ fundamental lemma for state-space systems and its extension to multiple datasets,” _IEEE Control Systems Letters_ , vol. 4, no. 3, pp. 602–607, 2020\.
  * [14] J. Coulson, J. Lygeros, and F. Dörfler, “Data-enabled predictive control: In the shallows of the DeePC,” in _2019 18th European Control Conference (ECC)_. IEEE, 2019, pp. 307–312.
  * [15] C. De Persis and P. Tesi, “Formulas for data-driven control: Stabilization, optimality, and robustness,” _IEEE Trans. Autom. Control_ , vol. 65, no. 3, pp. 909–924, 2019.
  * [16] J. Berberich, J. Köhler, M. A. Muller, and F. Allgower, “Data-driven model predictive control with stability and robustness guarantees,” _IEEE Trans. Autom. Control_ , vol. 66, no. 4, pp. 1702–1717, 2021.
  * [17] H. J. Van Waarde, J. Eising, H. L. Trentelman, and M. K. Camlibel, “Data informativity: a new perspective on data-driven analysis and control,” _IEEE Transactions on Automatic Control_ , vol. 65, no. 11, pp. 4753–4768, 2020.
  * [18] J. Coulson, J. Lygeros, and F. Dörfler, “Distributionally robust chance constrained data-enabled predictive control,” _arXiv preprint arXiv:2006.01702_ , 2020.
  * [19] N. Wieler, J. Berberich, A. Koch, and F. Allgöwer, “Data-driven controller design via finite-horizon dissipativity,” _arXiv preprint arXiv:2101.06156_ , 2021.
  * [20] Y. Lian and C. N. Jones, “Nonlinear data-enabled prediction and control,” _arXiv preprint arXiv:2101.03187_ , 2021.
  * [21] A. B. Alexandru, A. Tsiamis, and G. J. Pappas, “Data-driven control on encrypted data,” _arXiv preprint arXiv:2008.12671_ , 2020.
  * [22] L. Huang, J. Coulson, J. Lygeros, and F. Dörfler, “Decentralized data-enabled predictive control for power system oscillation damping,” _arXiv preprint arXiv:1911.12151_ , 2019.
  * [23] ——, “Data-enabled predictive control for grid-connected power converters,” in _2019 IEEE 58th Conference on Decision and Control (CDC)_. IEEE, 2019.
  * [24] P. G. Carlet, A. Favato, S. Bolognani, and F. Dörfler, “Data-driven predictive current control for synchronous motor drives,” in _2020 IEEE Energy Conversion Congress and Exposition (ECCE)_. IEEE, 2020, pp. 5148–5154.
  * [25] E. Elokda, J. Coulson, P. Beuchat, J. Lygeros, and F. Dörfler, “Data-enabled predictive control for quadcopters,” _ETH Zurich, Research Collection_ , 2019.
  * [26] A. Alanwar, Y. Stürz, and K. H. Johansson, “Robust data-driven predictive control using reachability analysis,” _arXiv preprint arXiv:2103.14110_ , 2021.
  * [27] A. Xue and N. Matni, “Data-driven system level synthesis,” _arXiv preprint arXiv:2011.10674_ , 2020.
  * [28] J. Coulson, J. Lygeros, and F. Dörfler, “Regularized and distributionally robust data-enabled predictive control,” in _2019 IEEE 58th Conference on Decision and Control (CDC)_. IEEE, 2019\.
  * [29] L. Huang, J. Zhen, J. Lygeros, and F. Dörfler, “Quadratic regularization of data-enabled predictive control: Theory and application to power converter experiments,” _arXiv preprint arXiv:2012.04434_ , 2020.
  * [30] I. Markovsky and P. Rapisarda, “Data-driven simulation and control,” _International Journal of Control_ , vol. 81, no. 12, pp. 1946–1959, 2008\.
  * [31] J. Birge and F. Louveaux, _Introduction to Stochastic Programming_. Springer Science & Business Media, 2011.
  * [32] J. Zhen and D. den Hertog, “Computing the maximum volume inscribed ellipsoid of a polytopic projection,” _INFORMS Journal on Computing_ , vol. 30, no. 1, pp. 31–42, 2017.
  * [33] A. Ben-Tal, D. den Hertog, and J.-P. Vial, “Deriving robust counterparts of nonlinear uncertain inequalities,” _Mathematical Programming_ , vol. 149, no. 1, pp. 265–299, 2015.
  * [34] A. Ben-Tal, L. El Ghaoui, and A. Nemirovski, _Robust Optimization_. Princeton Series in Applied Mathematics. Princeton University Press, Princeton, NJ, 2009.
  * [35] L. El Ghaoui and H. Lebret, “Robust solutions to least-squares problems with uncertain data,” _SIAM Journal on matrix analysis and applications_ , vol. 18, no. 4, pp. 1035–1064, 1997.
  * [36] D. Bertsimas and M. S. Copenhaver, “Characterization of the equivalence of robustification and regularization in linear and matrix regression,” _European Journal of Operational Research_ , vol. 270, no. 3, pp. 931–942, 2018.
  * [37] H. Xu, C. Caramanis, and S. Mannor, “Robust regression and lasso,” _IEEE Transactions on Information Theory_ , vol. 56, no. 7, pp. 3561–3574, 2010.
  * [38] R. Rockafellar, _Convex Analysis._ Princeton University Press, 1970.
  * [39] V. A. Yakubovich, “S-procedure in nonlinear control theory,” _Vestnik Leningradskogo Universiteta_ , vol. 1, pp. 62–77, 1971, (in Russian).
  * [40] S. Boyd and L. Vandenberghe, _Convex Optimization_. Cambridge University Press, 2004.
  * [41] L. Huang, H. Xin, Z. Li, P. Ju, H. Yuan, Z. Lan, and Z. Wang, “Grid-synchronization stability analysis and loop shaping for pll-based power converters with different reactive power control,” _IEEE Trans. Smart Grid_ , vol. 11, no. 1, pp. 501–516, 2019.
  * [42] Mathworks, _Band-Limited White Noise_ , 2021. [Online]. Available: www.mathworks.com/help/simulink/slref/bandlimitedwhitenoise.html
  * [43] R. J. Tibshirani _et al._ , “The lasso problem and uniqueness,” _Electronic Journal of statistics_ , vol. 7, pp. 1456–1490, 2013.