# Designing Experiments for Data-Driven Control of Nonlinear Systems

Claudio De Persis  Pietro Tesi  ENTEG, University of Groningen, Nijenborgh 4,
9747AG Groningen, The Netherlands (e-mail: c.de.persis@rug.nl).  DINFO,
University of Florence, 50139 Florence, Italy (e-mail<EMAIL_ADDRESS>
###### Abstract

In a recent paper we have shown that data collected from linear systems
excited by persistently exciting inputs during low-complexity experiments, can
be used to design state- and output-feedback controllers, including optimal
Linear Quadratic Regulators (LQR), by solving linear matrix inequalities (LMI)
and semidefinite programs. We have also shown how to stabilize in the first
approximation unknown nonlinear systems using data. In contrast to the case of
linear systems, however, in the case of nonlinear systems the conditions for
learning a controller directly from data may not be fulfilled even when the
data are collected in experiments performed using persistently exciting
inputs. In this paper we show how to design experiments that lead to the
fulfilment of these conditions.

###### keywords:

Nonlinear systems; Nonlinear control; Control system design; Data-driven
control; Convex programming.

##  1 Introduction

Recent advances in learning have prompted a renewed interest in the use of
data for control of complex dynamical systems in at least two different ways.
On one hand, the use of new learning-based identification techniques for
system identification followed by off-the-shelf top-notch robust control
design (Chiuso, 2016; Recht, 2019). On the other hand, the design of control
policies directly form data, skipping altogether any attempt of identifying
the system’s model (Campi et al., 2002; Campi and Savaresi, 2006; Formentin et
al., 2013; Novara et al., 2016; Tanaskovic et al., 2017). Related
contributions are the iterative feedback tuning (Hjalmarsson and Gevers,
1998), correlation-based tuning (Karimi et al., 2007), and the design of
controllers for all the systems compatible with the measured data via
polynomial optimization (Dai and Sznaier, 2018). Optimal (Gonçalves da Silva
et al., 2019; Baggio et al., 2019) and nonlinear (Wabersich and Zeilinger,
2018) control problems have also been investigated.

Inspired by the the so-called Fundamental Lemma of Willems et al. (2005), a
novel approach to simulate and control unknown dynamical systems without any
identification, but merely using finite-length input-output data, was
introduced in Markovsky and Rapisarda (2008). The spotlight on the ideas of
Markovsky et al. (2005); Markovsky and Rapisarda (2008) and their important
role in data-driven control design was turned on again by the recent paper by
Coulson et al. (2019a). See Coulson et al. (2019b) and Huang et al. (2019) for
recent follow-ups.

Also relying on the results of Willems et al. (2005), recently De Persis and
Tesi (2019b, a) have shown that data collected from linear systems excited by
persistently exciting inputs during low-complexity finite-horizon off-line
experiments, can be used to design state- and output-feedback controllers,
including optimal Linear Quadratic Regulators (LQR), by solving linear matrix
inequalities (LMI) and semidefinite programs. The authors have also shown how
these methods are robust to the use of data corrupted by bounded-but-unknown
deterministic noise, meaning that they return a stabilizing controller for
noisy data that satisfy a quantified signal-to-noise ratio.

In a continuation of De Persis and Tesi (2019b, a), the necessity of the
persistency-of-excitation conditions has been throughly discussed in van
Waarde et al. (2020), the design of controllers robust to process disturbances
has been further explored in Berberich et al. (2019) and data-based guarantees
for polyhedral set-invariance properties (safe controllers) have been studied
in Bisoffi et al. (2019).

Another result of De Persis and Tesi (2019a) was to extend the stabilization
result to the case of nonlinear systems (see Theorem 1 in Section 2 below).
Namely, if the system is nonlinear and if the data collected during a finite-
length off-line experiment satisfy suitable conditions (Assumptions 1 and 2
below) then we can design a controller from these data that stabilizes the
system in the first approximation. In contrast with what happens with linear
systems, however, where it can be shown that these conditions are always
satisfied if persistently exciting inputs are applied in the experiment, in
the case of nonlinear systems the fulfilment of Assumptions 1 and 2 is a more
challenging task.

The purpose of this note is threefold. First we show that experiments that are
performed on nonlinear systems may fail to generate data that satisfy
Assumptions 1 and 2, even when persistently exciting inputs are applied.
Second, we show how to design experiments that lead to the fulfilment of these
assumptions. Third, we provide some heuristic considerations on how to verify
that the stabilizing controller that we design is actually stabilizing when no
a priori knowledge about the nonlinear system is available.

In Section 2 we recall the results of De Persis and Tesi (2019a). The main
results are in Section 3. Conclusions are drawn in Section 4.

##  2 Preliminaries

We start by recalling the result given in De Persis and Tesi (2019a).

Consider a smooth nonlinear system

$x(k+1)=f(x(k),u(k))$ (1)

and let $(\overline{x},\overline{u})$ be a _known_ equilibrium pair, that is
such that $\overline{x}=f(\overline{x},\overline{u})$. Throughout his note,
for ease of notation we will assume without loss of generality that the
equilibrium point of interest is $(\overline{x},\overline{u})=(0,0)$.

Let us rewrite the nonlinear system as

$\displaystyle x(k+1)=Ax(k)+Bu(k)+d(k)$ (2)

where

$\displaystyle A:=\left.\frac{\partial f}{\partial
x}\right|_{(x,u)=(0,0)},\quad B:=\left.\frac{\partial f}{\partial
u}\right|_{(x,u)=(0,0)}\,.$ (3)

The quantity $d$ accounts for higher-order terms and it has the property that
it goes to zero faster than $x$ and $u$, namely we have

$\displaystyle d=R(x,u)\begin{bmatrix}x\\\ u\end{bmatrix}$

with $R(x,u)$ a matrix of smooth functions with the property that

$\displaystyle\lim_{{\tiny\begin{bmatrix}x\\\
u\end{bmatrix}\to{0}}}R(x,u)={0}$ (4)

It is known that if the pair $(A,B)$ defining the linearized system is
stabilizable then the controller $K$ rendering $A+BK$ stable exponentially
stabilizes the equilibrium of the nonlinear system. Let

$\displaystyle X_{0,T}$ $\displaystyle:=$
$\displaystyle\left[\begin{array}[]{cccc}x(0)&x(1)&\cdots&x(T-1)\end{array}\right]$
$\displaystyle X_{1,T}$ $\displaystyle:=$
$\displaystyle\left[\begin{array}[]{cccc}x(1)&x(2)&\cdots&x(T)\end{array}\right]$
$\displaystyle U_{0,T}$ $\displaystyle:=$
$\displaystyle\left[\begin{array}[]{cccc}u(0)&u(1)&\cdots&u(T-1)\end{array}\right]$
$\displaystyle D_{0,T}$ $\displaystyle:=$
$\displaystyle\left[\begin{array}[]{cccc}d(0)&d(1)&\cdots&d(T-1)\end{array}\right]$

be data resulting from an experiment carried out on the nonlinear system (1).
Note that the matrices $X_{0,T}$, $X_{1,T}$ and $U_{0,T}$ are available from
data.

Consider the following assumptions.

###### Assumption 1

The matrices

$\displaystyle\begin{bmatrix}U_{0,T}\\\ \hline\cr
X_{0,T}\end{bmatrix},\;X_{1,T}$ (9)

have full row rank.

###### Assumption 2

It holds that

$\displaystyle D_{0,T}D_{0,T}^{\top}\preceq\gamma X_{1,T}X_{1,T}^{\top}$ (10)

for some $\gamma>0$.

The following result holds.

######  Theorem 1 (De Persis and Tesi (2019a))

Consider a nonlinear system as in (1), along with an equilibrium pair
$(\overline{x},\overline{u})$. Suppose that Assumptions 1 and 2 hold. Then,
any solution $(Q,\alpha)$ to

$\displaystyle\begin{array}[]{l}\max_{(Q,\alpha)}\alpha\\\\[11.38092pt]
\textrm{subject to}\\\\[8.5359pt]
\begin{array}[]{rl}\begin{bmatrix}X_{0,T}\,Q-\alpha
X_{1,T}X_{1,T}^{\top}&X_{1,T}Q\\\
Q^{\top}X_{1,T}^{\top}&X_{0,T}\,Q\end{bmatrix}\succeq&0\\\\[11.38109pt]
\begin{bmatrix}I_{T}&Q\\\
Q^{\top}&X_{0,T}\,Q\end{bmatrix}\succeq&0\end{array}\end{array}$ (16)

such that $\gamma<\alpha^{2}/(4+2\alpha)$, $\alpha>0$, returns a state-
feedback controller $K=U_{0,T}Q(X_{0,T}Q)^{-1}$ that locally stabilizes the
equilibrium $(\overline{x},\overline{u})$. $\blacksquare$

The reason for considering the maximization of $\alpha$ is to render condition
$\gamma<\alpha^{2}/(4+2\alpha)$ easier to fulfil. Notice that (16) is a
semidefinite program.

Theorem 1 rests on the fact that Assumptions 1 and 2 hold and that (16) admits
a solution with $\gamma<\alpha^{2}/(4+2\alpha)$. For _linear_ controllable
systems, this is always the case provided that the experiments are carried out
with persistently exciting inputs of a sufficiently high order (see (De Persis
and Tesi, 2019a, Definition 1) for a definition of persistency of excitation).
In fact, under these conditions: (i) Assumptions 1 and 2 are satisfied; (ii)
problem (16) is feasible, and any solution is such that
$\gamma<\alpha^{2}/(4+2\alpha)$ (since $D_{0,T}=0$ Assumption 2 holds with an
arbitrary $\gamma$). We summarize this fact.

######  Theorem 2 (De Persis and Tesi (2019a))

Let (1) be a linear controllable system with state-space dimension $n$, and
consider an experiment carried out with a persistently exciting input of order
$n+1$. Then, Assumptions 1 and 2 hold. Further, problem (16) is feasible and
any solution is such that $\gamma<\alpha^{2}/(4+2\alpha)$. Thus, any solution
returns a stabilizing state-feedback gain $K=U_{0,T}Q(X_{0,T}Q)^{-1}$.
$\blacksquare$

In the sequel, it will be useful to relate the evolution of the nonlinear
system with the one of the corresponding linearized system, which we express
as

$\displaystyle x^{l}(k+1)=Ax^{l}(k)+Bu^{l}(k)$ (17)

We also let

$\displaystyle X^{l}_{0,T}$ $\displaystyle:=$
$\displaystyle\left[\begin{array}[]{cccc}x^{l}(0)&x^{l}(1)&\cdots&x^{l}(T-1)\end{array}\right]$
$\displaystyle X^{l}_{1,T}$ $\displaystyle:=$
$\displaystyle\left[\begin{array}[]{cccc}x^{l}(1)&x^{l}(2)&\cdots&x^{l}(T)\end{array}\right]$
$\displaystyle U^{l}_{0,T}$ $\displaystyle:=$
$\displaystyle\left[\begin{array}[]{cccc}u^{l}(0)&u^{l}(1)&\cdots&u^{l}(T-1)\end{array}\right]$

be the data resulting from an (hypothetical) experiment made on the linearized
system (17). Note that the matrices $X^{l}_{0,T}$, $X^{l}_{1,T}$ and
$U^{l}_{0,T}$ are not available from data.

##  3 Main results

Since around an equilibrium point nonlinear systems behave as linear systems,
it is tempting to conclude that Theorem 2 can be extended to the nonlinear
case as long as the experiments are carried out sufficiently close to the
equilibrium point. As we will see, this is true _only if_ the experiments are
carried out in a certain manner. More precisely, we will see that there indeed
exist experiments that ensure the same properties as in Theorem 2, but not all
experiments (even with persistently exciting inputs and arbitrarily close to
the equilibrium) guarantee these properties. In the remainder of this section:

  1. 1.

we show that there exist experiments (carried out with persistently exciting
inputs and arbitrarily close to the equilibrium) for which Assumptions 1 and 2
fail to hold;

  2. 2.

we show that there exist experiments (carried out with persistently exciting
inputs and sufficiently close to the equilibrium) ensuring the same properties
as in Theorem 2, and we characterize a class of experiments ensuring these
properties;

  3. 3.

we provide some heuristic considerations on how to verify that a solution to
(16) returns a stabilizing controller. (In general, one cannot determine how
close to the equilibrium the experiments should be run unless we have some
prior knowledge on $D_{0,T}$, that is on the type of nonlinearity; thus, in
general one cannot assess the fulfilment of the condition
$\gamma<\alpha^{2}/(4+2\alpha)$ since $\gamma$ is unknown.)

###  3.1 Experiment design issues for nonlinear systems

We start with showing that there exist experiments (carried out with
persistently exciting inputs and arbitrarily close to the equilibrium) for
which Assumptions 1 and 2 fail to hold. We show this fact through an example.

Example. Consider the nonlinear system

$x(k+1)=x(k)^{2}+u(k)$

with equilibrium $(\overline{x},\overline{u})=(0,0)$. The system can be
written as

$x(k+1)=u(k)+d(k),\quad d(k)=x(k)^{2}$

with $A=0$ and $B=1$, and the corresponding linearized system is given by

$x^{l}(k+1)=u^{l}(k)$

Consider now the experiment given by

$\displaystyle x(0)=x_{0}$ $\displaystyle
U_{0,T}=\left[\begin{array}[]{ccc}u(0)&u(1)&u(2)\end{array}\right]$

where

$x_{0}=u(0)=\theta,\,\,u(1)=\theta+\theta^{2},\,\,u(2)=\theta+\theta^{2}+(\theta+\theta^{2})^{2}$

with $\theta$ real. The input is persistently exciting of order $n+1=2$. Now,
If we hypothetically apply this experiment to the linearized system, that is
we let $x^{l}(0)=x_{0}$ and
${\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}U^{l}_{0,T}=U_{0,T}}$,
we obtain

$\left[\begin{array}[]{c}U^{l}_{0,T}\\\
X^{l}_{0,T}\end{array}\right]=\left[\begin{array}[]{ccc}\theta&\theta+\theta^{2}&\theta+\theta^{2}+(\theta+\theta^{2})^{2}\\\
\theta&\theta&\theta+\theta^{2}\end{array}\right]$

which is full row rank for all $\theta\neq 0$, and the same holds for
$X^{l}_{1,T}$. However,

$\begin{array}[]{rl}\left[\begin{array}[]{c}U_{0,T}\\\
X_{0,T}\end{array}\right]=&\left[\begin{array}[]{ccc}u(0)&u(1)&u(2)\\\
x(0)&u(0)+x(0)^{2}&u(1)+x(1)^{2}\end{array}\right]\\\\[11.38109pt]
=&\left[\begin{array}[]{ccc}\theta&\theta+\theta^{2}&\theta+\theta^{2}+(\theta+\theta^{2})^{2}\\\
\theta&\theta+\theta^{2}&\theta+\theta^{2}+(\theta+\theta^{2})^{2}\end{array}\right]\end{array}$

has rank $1$ for all $\theta\neq 0$. We conclude that Assumption 1 does not
hold. Note that this holds true even though the input sequence if persistently
exciting and $\theta$ is arbitrarily small. $\blacksquare$

This simple example shows that for nonlinear systems there might exist
experiments for which Assumption 1 is not satisfied, even if these experiments
originate from persistently exciting inputs and are carried out arbitrarily
close to the equilibrium.

An intuitive explanation for this fact can be obtained by relating the
matrices associated with the linear and nonlinear system. Consider an
experiment on the nonlinear system with initial condition $x_{0}$ and input
sequence $\\{u(0),u(1),\ldots,u(T-1)\\}$ resulting in state matrices $X_{0,T}$
and $X_{1,T}$. It is simple to see that

$x(k)=x^{l}(k)+\sum_{i=0}^{k-1}A^{k-i-1}d(i)$

This implies that the matrices in Assumption 1 can be written as

$\displaystyle\begin{bmatrix}U_{0,T}\\\ \hline\cr
X_{0,T}\end{bmatrix}=\begin{bmatrix}U^{l}_{0,T}\\\ \hline\cr
X^{l}_{0,T}\end{bmatrix}+\begin{bmatrix}0\\\
\hline\cr\Xi\end{bmatrix},\,\,X_{1,T}=X^{l}_{1,T}+\Psi$ (22)

where

$\Xi:=\begin{bmatrix}0&d(0)&Ad(0)+d(1)&A^{2}d(0)+Ad(1)+d(2)&\cdots&\end{bmatrix}$
(23)

and

$\Psi:=\begin{bmatrix}d(0)&Ad(0)+d(1)&A^{2}d(0)+Ad(1)+d(2)&\cdots&\end{bmatrix}$
(24)

In connection with the previous example, this translates to

$\begin{array}[]{rl}\left[\begin{array}[]{c}U^{l}_{0,T}\\\
X^{l}_{0,T}\end{array}\right]=&\left[\begin{array}[]{ccc}\theta&\theta+\theta^{2}&\theta+\theta^{2}+(\theta+\theta^{2})^{2}\\\
\theta&\theta&\theta+\theta^{2}\end{array}\right]\\\\[11.38109pt]
\begin{bmatrix}0\\\ \Xi\end{bmatrix}=&\begin{bmatrix}0&0&0\\\
0&\theta^{2}&(\theta+\theta^{2})^{2}\end{bmatrix}\end{array}$

This means that there exist trajectories (i.e. experiments) for which the
perturbation that causes

$\left[\begin{array}[]{c}U^{l}_{0,T}\\\ X^{l}_{0,T}\end{array}\right]$ (25)

to lose rank is of the same order ($\theta^{2}$ in this example) as the
nonlinear terms, thus of the same order as $\Xi$. In the sequel, we will show
that one can nonetheless design experiments so as to satisfy Assumptions 1 and
2.

###  3.2 Scaling the experiments ensures Assumptions 1 and 2

The basic idea is to show that there also exist experiments such that $\Xi$
vanishes faster than the perturbation causing the matrix (25) to lose rank.

Consider an experiment with initial condition $x_{0}$ and input sequence

$u_{0,T}=\\{u(0),u(1),\ldots,u(T-1)\\}$

persistently exciting of order $n+1$. For this experiment, relation (22)
holds. Consider now a _scaled_ version of this experiment with initial
condition $\varepsilon x_{0}$ and input sequence $\varepsilon u_{0,T}$ with
$\varepsilon>0$ real. Denote by $\tilde{U}_{0,T}$ the corresponding matrix
which satisfies

$\displaystyle\tilde{U}_{0,T}$ $\displaystyle=$
$\displaystyle\begin{bmatrix}\tilde{u}(0)&\tilde{u}(1)&\cdots&\tilde{u}(T-1)\end{bmatrix}$
$\displaystyle=$ $\displaystyle\varepsilon U_{0,T}=\varepsilon U^{l}_{0,T}$

Also, denote by $\tilde{X}_{0,T}$ and $\tilde{X}_{1,T}$ the state matrices of
the nonlinear system resulting from this new experiment,

$\begin{array}[]{rl}\tilde{X}_{0,T}=&\begin{bmatrix}\tilde{x}(0)&\tilde{x}(1)&\cdots&\tilde{x}(T-1)\end{bmatrix}\\\\[14.22636pt]
\tilde{X}_{1,T}=&\begin{bmatrix}\tilde{x}(1)&\tilde{x}(2)&&\cdots&\tilde{x}(T)\end{bmatrix}\end{array}$

where $\tilde{x}$ satisfies

$\displaystyle\tilde{x}(k+1)$ $\displaystyle=$ $\displaystyle
f(\tilde{x}(k),\tilde{u}(k))$ $\displaystyle=$ $\displaystyle
A\tilde{x}(k)+B\tilde{x}(k)+\tilde{d}(k)$

with $\tilde{x}(0)=\varepsilon x_{0}$, and where, as before,

$\displaystyle\tilde{d}=R(\tilde{x},\tilde{u})\begin{bmatrix}\tilde{x}\\\
\tilde{u}\end{bmatrix}$

with $R(\tilde{x},\tilde{u})$ a matrix of smooth functions with the property
that

$\displaystyle\lim_{{\tiny\begin{bmatrix}\tilde{x}\\\
\tilde{u}\end{bmatrix}\to{0}}}R(\tilde{x},\tilde{u})={0}$ (26)

Now, for this new experiment it holds that

$\displaystyle\begin{bmatrix}\tilde{U}_{0,T}\\\
\hline\cr\tilde{X}_{0,T}\end{bmatrix}=\varepsilon\begin{bmatrix}U^{l}_{0,T}\\\
\hline\cr X^{l}_{0,T}\end{bmatrix}+\begin{bmatrix}0\\\
\hline\cr\tilde{\Xi}\end{bmatrix},\,\,\tilde{X}_{1,T}=\varepsilon
X^{l}_{1,T}+\tilde{\Psi}$ (27)

where the matrices $\tilde{\Xi}$ and $\tilde{\Psi}$ are defined as in (23) and
(24) with $d$ replaced by $\tilde{d}$. We will only show the first of (27)
since the reasoning for the second is analogous.

In order to prove the first of (27) it is sufficient to show that
$\tilde{X}_{0,T}=\varepsilon X^{l}_{0,T}+\tilde{\Xi}$ since the relation
$\tilde{U}_{0,T}=\varepsilon U^{l}_{0,T}$ holds by construction. The result
can be proved by induction on the entries $\tilde{x}$ of the matrix
$\tilde{X}_{0,T}$. The claim holds for $k=0$. Suppose that the claim holds up
to a certain $k\geq 0$. We have

$\displaystyle\tilde{x}(k+1)$ $\displaystyle=$ $\displaystyle
A\tilde{x}(k)+B\tilde{u}(k)+\tilde{d}(k)$ (28) $\displaystyle=$ $\displaystyle
A\left(\varepsilon
x^{l}(k)+\sum_{i=0}^{k-1}A^{k-i-1}\tilde{d}(i)\right)+B\tilde{u}(k)+\tilde{d}(k)$
$\displaystyle=$
$\displaystyle\varepsilon\left(Ax^{l}(k)+Bu^{l}(k)\right)+\sum_{i=0}^{k-1}A^{k-i}\tilde{d}(i)+\tilde{d}(k)$
$\displaystyle=$ $\displaystyle\varepsilon
x^{l}(k+1)+\sum_{i=0}^{k}A^{k-i}\tilde{d}(i)$

which gives the claim.

The important fact resulting from these relations is that the matrices
$X^{l}_{0,T}$ and $X^{l}_{1,T}$ are _fixed_ , and the perturbation which cause
them to lose rank depends on $\varepsilon$. To this end, we recall the
following result.

###### Theorem 3

(Dahleh et al., 2003, Theorem 5.1) Suppose $M\in\mathbb{C}^{m\times n}$ is
full column rank. Then,

$\min_{\Delta\in\mathbb{C}^{m\times
n}}\left\\{\|\Delta\|\,\,|\,\,M+\Delta\textrm{ has
rank}<n\right\\}=\sigma_{n}(M)$

where $\|\cdot\|$ denotes the Euclidean norm, and where $\sigma_{n}$ is the
smallest singular value of the matrix $M$. $\blacksquare$

The reasoning for full row rank matrices is analogous (one can alternately
consider (27) with the transpose).

Recalling that the singular values of a matrix $M$ are the square roots of the
eigenvalues of the matrix $MM^{\top}$, in connection with (27) it follows that
the smallest perturbation that causes $\varepsilon M$ with

$M:=\begin{bmatrix}U^{l}_{0,1,T}\\\ \hline\cr X^{l}_{0,T}\end{bmatrix}$

to lose rank is of order $\varepsilon$. In fact, since $M$ is full row rank
then the smallest eigenvalue of the matrix $MM^{\top}$ is strictly positive.
Denoting this eigenvalue by $\underline{\lambda}$, it follows that the
smallest singular value of $\varepsilon M$ is
$\varepsilon\sqrt{\underline{\lambda}}$. Hence, to prove that the perturbed
matrix in (27) is full row rank, it is sufficient to show that the matrix
$\tilde{\Xi}$ goes to zero faster than $\varepsilon$ since this implies the
existence of a value $\overline{\varepsilon}$ such that

$\|\tilde{\Xi}\|<\varepsilon\sqrt{\underline{\lambda}}\quad\forall\varepsilon\in(0,\overline{\varepsilon})$
(29)

This fact can be proven by looking at the various elements $\tilde{d}$ which
form the entries of the matrix $\tilde{\Xi}$, since

$\|\tilde{\Xi}\|^{2}\leq\displaystyle\sum_{k=1}^{T-1}{\rm
trace}\left[\left(\sum_{i=0}^{k-1}A^{k-1-i}\tilde{d}(i)\right)^{\top}\left(\sum_{i=0}^{k-1}A^{k-1-i}\tilde{d}(i)\right)\right]$

with the right-hand side being the Frobenius norm $\|\tilde{\Xi}\|_{F}^{2}$.
In fact, if one can prove that $\tilde{d}$ converges to zero faster than
$\varepsilon$, then the bound (29) on the maximum allowable perturbation
holds.

Claim. For all $k\geq 0$, $\tilde{d}(k)$ converges to zero faster than
$\varepsilon$.

Proof of the Claim. Recall that, for all $k$,

$\tilde{d}(k)=R(\tilde{x}(k),\tilde{u}(k))\begin{bmatrix}\tilde{x}(k)\\\
\tilde{u}(k)\end{bmatrix}$ (30)

with $R(\tilde{x}(k),\tilde{u}(k))$ a matrix of smooth functions with the
property that

$\displaystyle\lim_{{\tiny\begin{bmatrix}\tilde{x}(k)\\\
\tilde{u}(k)\end{bmatrix}\to{0}}}R(\tilde{x}(k),\tilde{u}(k))={0}$ (31)

For $k=0$, we have

$\begin{bmatrix}\tilde{x}(0)\\\
\tilde{u}(0)\end{bmatrix}=\varepsilon\begin{bmatrix}x(0)\\\
u(0)\end{bmatrix},$

hence, by (31), $\lim_{\varepsilon\to{0}}R(\tilde{x}(0),\tilde{u}(0))={0}$,
i.e. $R(\tilde{x}(0),\tilde{u}(0))$ converges to zero as fast as
$\varepsilon$. The two just established facts and (30) written for $k=0$ allow
us to conclude that $\tilde{d}(0)$ goes to zero faster than $\varepsilon$. We
now proceed by induction. Let us assume that for some $k\geq 0$, the vectors
$\begin{bmatrix}\tilde{x}(k)\\\ \tilde{u}(k)\end{bmatrix}$ and $\tilde{d}(k)$
converge to zero at least as fast as $\varepsilon$ and faster than
$\varepsilon$, respectively, and that
$\lim_{\varepsilon\to{0}}R(\tilde{x}(k),\tilde{u}(k))={0}$. Since

$\tilde{x}(k+1)=\begin{bmatrix}A&B\end{bmatrix}\begin{bmatrix}\tilde{x}(k)\\\
\tilde{u}(k)\end{bmatrix}+\tilde{d}(k)$

then $\tilde{x}(k+1)$ converges to zero at least as fast as $\varepsilon$ and
the same property holds for $\begin{bmatrix}\tilde{x}(k+1)\\\
\tilde{u}(k+1)\end{bmatrix}$ by design, since $\tilde{u}(k+1)=\varepsilon
u(k)$. The identity (31) written at time $k+1$ implies that
$\lim_{\varepsilon\to{0}}R(\tilde{x}(k+1),\tilde{u}(k+1))={0}$. By (30)
written for $k+1$, we conclude that $\tilde{d}(k+1)$ converges to zero faster
than $\varepsilon$ as claimed. $\blacksquare$

Similar arguments can be used to show Assumption 2. To this end, recall from
(27) that $\tilde{X}_{1,T}=\varepsilon X^{l}_{1,T}+\tilde{\Psi}$, where
$X^{l}_{1,T}$ has full row rank. Thus, by applying the Young’s inequality we
obtain

$\displaystyle\tilde{X}_{1,T}\tilde{X}_{1,T}^{\top}$ $\displaystyle=$
$\displaystyle(\varepsilon X^{l}_{1,T}+\tilde{\Psi})(\varepsilon
X^{l}_{1,T}+\tilde{\Psi})^{\top}$ (32) $\displaystyle\succeq$
$\displaystyle\frac{\varepsilon^{2}}{2}X^{l}_{1,T}(X^{l}_{1,T})^{\top}-\tilde{\Psi}\tilde{\Psi}^{\top}$

Hence, a sufficient condition for Assumption 2 to hold is that

$\frac{\gamma\varepsilon^{2}}{2}X^{l}_{1,T}(X^{l}_{1,T})^{\top}\succeq\gamma\tilde{\Psi}\tilde{\Psi}^{\top}+\tilde{D}_{0,T}\tilde{D}_{0,T}^{\top}$

The result follows by noting that $X^{l}_{1,T}$ does not depend on
$\varepsilon$ and is full row rank, while the entries of $\tilde{\Psi}$ and
$\tilde{D}_{0,T}$ are functions of $\tilde{d}$ which goes to zero faster than
$\varepsilon$ (meaning that $\tilde{\Psi}\tilde{\Psi}^{\top}$ and
$\tilde{D}_{0,T}\tilde{D}_{0,T}^{\top}$ converge to zero faster than
$\varepsilon^{2}$).

We summarize the results in the following theorem. Note that, by an abuse of
notation, in the statement below as well as in the subsequent Theorem 5, when
we refer to Assumptions 1 and 2 and to problem (16), we let the matrices
$U_{0,1,T},X_{0,T},X_{1,T},D_{0,T}$ therein to be replaced by the matrices
$\tilde{U}_{0,1,T},\tilde{X}_{0,T},\tilde{X}_{1,T},\tilde{D}_{0,T}$.

###### Theorem 4

Consider a nonlinear system as in (1), with state-space dimension $n$ and with
equilibrium $(\overline{x},\overline{u})=(0,0)$, and let the corresponding
linearized system be controllable. Consider any experiment
$(x_{0},u_{[0,T-1]})$ with persistently exciting input $u_{[0,T-1]}$ of order
$n+1$. Then, for any $\gamma>0$, there exists $\overline{\varepsilon}$ such
that, for all $\varepsilon\in(0,\overline{\varepsilon})$, the experiment
$(\varepsilon x_{0},\varepsilon u_{[0,T-1]})$ satisfies Assumptions 1 and 2.
$\blacksquare$

Theorem 4 gives a principled method to satisfy Assumptions 1 and 2. In fact,
if one can perform multiple experiments on the system at the equilibrium
($x_{0}=0$), one can apply scaled versions of the input sequence $u_{[0,T-1]}$
until Assumptions 1 and 2 are satisfied. In connection with the example of
Section 3.1, an intuitive explanation of Theorem 4 is given in Figure 1.

Figure 1: Pictorial representation of Theorem 4 for the example of Section
3.1. The blue curve represents the set of points satisfying
$[u(0)\,\,\,u(1)]=[\theta\,\,\,\theta+\theta^{2}]$ for which Assumption 1
fails to hold (for ease of illustration, we only report the first two entries
of the input). For each experiment (circles in the blue curve), the scaled
versions $\varepsilon[u(0)\,\,\,u(1)]$ are represented by the coloured
straight lines.

We can further strengthen the result.

###### Theorem 5

Consider a nonlinear system as in (1), with state-space dimension $n$ and with
equilibrium $(\overline{x},\overline{u})=(0,0)$, and let the corresponding
linearized system be controllable. Consider any experiment $(\varepsilon
x_{0},\varepsilon u_{[0,T-1]})$ with persistently exciting input sequence
$u_{[0,T-1]}$ of order $n+1$. Then there exists $\overline{\varepsilon}$ such
that, for all $\varepsilon\in(0,\overline{\varepsilon})$:

  1. (i)

Assumptions 1 and 2 are satisfied;

  2. (ii)

problem (16) is feasible;

  3. (ii)

the solution returns a stabilizing controller.

$\blacksquare$

_Proof._ Point (i) has been previously shown.

(ii). Take any experiment $(x_{0},u_{[0,T-1]})$ with persistently exciting
input $u_{[0,T-1]}$. If this experiment is applied to the linearized system,
namely $(x^{l}_{0},u^{l}_{[0,T-1]})=(x_{0},u_{[0,T-1]})$, we obtain from (De
Persis and Tesi, 2019a, Theorem 3) that

$\displaystyle\begin{bmatrix}X^{l}_{0,T}\,\overline{Q}&X^{l}_{1,T}\overline{Q}\\\
\overline{Q}^{\top}(X^{l}_{1,T})^{\top}&X^{l}_{0,T}\,\overline{Q}\end{bmatrix}\succ
0$ (33)

for some matrix $\overline{Q}$. (We resort to Theorem 3 of De Persis and Tesi
(2019a) since it involves strict inequalities).

Recall now that for any scaled experiment $(\varepsilon x_{0},\varepsilon
u_{[0,T-1]})$ we have

$\displaystyle\tilde{X}_{0,T}=\varepsilon
X^{l}_{0,T}+\tilde{\Xi},\,\,\tilde{X}_{1,T}=\varepsilon
X^{l}_{1,T}+\tilde{\Psi}$ (34)

where $\tilde{\Xi}$ and $\tilde{\Psi}$ converge to zero faster than
$\varepsilon$. Thus, since $\overline{Q}$ is fixed, by taking $\varepsilon$
sufficiently small we obtain

$\begin{array}[]{rl}\begin{bmatrix}\tilde{X}_{0,T}\,\overline{Q}&\tilde{X}_{1,T}\overline{Q}\\\
\overline{Q}^{\top}\tilde{X}_{1,T}^{\top}&\tilde{X}_{0,T}\,\overline{Q}\end{bmatrix}=&\varepsilon\begin{bmatrix}{X}_{0,T}^{l}\,\overline{Q}&{X}_{1,T}^{l}\overline{Q}\\\
\overline{Q}^{\top}({X}_{1,T}^{l})^{\top}&{X}_{0,T}^{l}\,\overline{Q}\end{bmatrix}\\\
&+\begin{bmatrix}\tilde{\Xi}\,\overline{Q}&\tilde{\Psi}\overline{Q}\\\
\overline{Q}^{\top}\tilde{\Psi}^{\top}&\tilde{\Xi}\,\overline{Q}\end{bmatrix}\succ
0\end{array}$

This also implies that

$\displaystyle\begin{bmatrix}\tilde{X}_{0,T}\,\overline{Q}-\alpha\tilde{X}_{1,T}\tilde{X}_{1,T}^{\top}&\tilde{X}_{1,T}\overline{Q}\\\
\overline{Q}^{\top}\tilde{X}_{1,T}^{\top}&\tilde{X}_{0,T}\,\overline{Q}\end{bmatrix}\succ
0$ (35)

for sufficiently small $\alpha$. This shows the feasibility of the first of
(16).

Regarding the second of (16), first notice that the above LMI remains
satisfied if we scale it by $\delta>0$, that is

$\displaystyle\begin{bmatrix}\tilde{X}_{0,T}\,\hat{Q}-\delta\alpha\tilde{X}_{1,T}\tilde{X}_{1,T}^{\top}&\tilde{X}_{1,T}\hat{Q}\\\
\hat{Q}^{\top}\tilde{X}_{1,T}^{\top}&\tilde{X}_{0,T}\,\hat{Q}\end{bmatrix}\succ
0,\quad\hat{Q}:=\delta\overline{Q}$ (36)

for all $\delta>0$. Consider now the matrix

$I_{T}-\hat{Q}(\tilde{X}_{0,T}\,\hat{Q})^{-1}\hat{Q}^{\top}$

which is well defined since $\tilde{X}_{0,T}\,\hat{Q}\succ 0$. Exploiting the
relation $\hat{Q}:=\delta\overline{Q}$, we obtain

$I_{T}-\hat{Q}(\tilde{X}_{0,T}\,\hat{Q})^{-1}\hat{Q}^{\top}=I_{T}-\delta\overline{Q}(\tilde{X}_{0,T}\,\overline{Q})^{-1}\overline{Q}^{\top}\succ
0$

for sufficiently small $\delta$. By the Schur complement, this shows the
feasibility of the second of (16), which together with (36) concludes the
proof of (ii).

(iii). Here we exploit the fact that in the design formulation we search for
the solution maximizing $\alpha$. Let $(\overline{Q},\overline{\alpha})$ be
the solution to (16) associated to the experiment $(x_{0},u_{[0,T-1]})$
carried out on the linearized system, with corresponding matrices
$X^{l}_{0,T}$ and $X^{l}_{1,T}$. Consider the scaled experiment $(\varepsilon
x_{0},\varepsilon u_{[0,T-1]})$ still on the linearized system, resulting in
state matrices $\varepsilon X^{l}_{0,T}$ and $\varepsilon X^{l}_{1,T}$. Denote
its solution by $(Q,\alpha)$. It follows that $\alpha=\overline{\alpha}$, that
is the optimal value of $\alpha$ does not change. To see this fact, we first
show that $(Q,\alpha)=(\varepsilon\overline{Q},\overline{\alpha})$ is
feasible. For the first of (16), the scaled experiment $(\varepsilon
x_{0},\varepsilon u_{[0,T-1]})$ gives

$\displaystyle\begin{bmatrix}\varepsilon X^{l}_{0,T}\,Q-\varepsilon^{2}\alpha
X_{1,T}^{l}(X_{1,T}^{l})^{\top}&\varepsilon X_{1,T}^{l}Q\\\ \varepsilon
Q^{\top}(X_{1,T}^{l})^{\top}&\varepsilon X_{0,T}^{l}\,Q\end{bmatrix}$
$\displaystyle=$
$\displaystyle\begin{bmatrix}\varepsilon^{2}X^{l}_{0,T}\,\overline{Q}-\varepsilon^{2}\overline{\alpha}X_{1,T}^{l}(X_{1,T}^{l})^{\top}&\varepsilon^{2}X_{1,T}^{l}\overline{Q}\\\
\varepsilon^{2}\overline{Q}^{\top}(X_{1,T}^{l})^{\top}&\varepsilon^{2}X_{0,T}^{l}\,\overline{Q}\end{bmatrix}\succeq
0$ (37)

since $(\overline{Q},\overline{\alpha})$ is a feasible solution to (16)
associated with the experiment $(x_{0},u_{[0,T-1]})$.

Concerning the second of (16), the scaled experiment $(\varepsilon
x_{0},\varepsilon u_{[0,T-1]})$ gives

$\displaystyle\begin{bmatrix}I_{T}&Q\\\ Q^{\top}&\varepsilon
X_{0,T}^{l}\,Q\end{bmatrix}=\begin{bmatrix}I_{T}&\varepsilon\overline{Q}\\\
\varepsilon\overline{Q}^{\top}&\varepsilon^{2}X_{0,T}^{l}\,\overline{Q}\end{bmatrix}\succeq
0$ (38)

since $(\overline{Q},\overline{\alpha})$ is a feasible solution to (16)
associated with $(x_{0},u_{[0,T-1]})$. In fact, by the Schur complement, (38)
holds if and only if
$I_{T}-\overline{Q}(X_{0,T}^{l}\,\overline{Q})^{-1}\overline{Q}^{\top}\succeq
0$.

We finally show that for the scaled experiment we cannot have a solution with
$\alpha>\overline{\alpha}$. Consider the solution $(Q,\alpha)$ associated with
the scaled experiment $(\varepsilon x_{0},\varepsilon u_{[0,T-1]})$ and
suppose by contradiction that $\alpha>\overline{\alpha}$. We show that in this
case $(\hat{Q},\hat{\alpha}):=(Q/\varepsilon,\alpha)$ is a solution associated
with the experiment $(x_{0},u_{[0,T-1]})$, contradicting the fact that
$\overline{\alpha}$ is the optimal value. For the first of (16) the solution
$(Q,\alpha)$ associated with the scaled experiment $(\varepsilon
x_{0},\varepsilon u_{[0,T-1]})$ gives

$\displaystyle\begin{bmatrix}\varepsilon
X^{l}_{0,T}\,Q-\alpha\varepsilon^{2}X_{1,T}^{l}(X_{1,T}^{l})^{\top}&\varepsilon
X_{1,T}^{l}Q\\\ \varepsilon Q^{\top}(X_{1,T}^{l})^{\top}&\varepsilon
X_{0,T}^{l}\,Q\end{bmatrix}\succeq 0$ (39)

From the relation $(\hat{Q},\hat{\alpha})=(Q/\varepsilon,\alpha)$ we would
thus get

$\displaystyle\begin{bmatrix}\varepsilon^{2}X^{l}_{0,T}\,\hat{Q}-\varepsilon^{2}\hat{\alpha}X_{1,T}^{l}(X_{1,T}^{l})^{\top}&\varepsilon^{2}X_{1,T}^{l}\hat{Q}\\\
\varepsilon^{2}\hat{Q}^{\top}(X_{1,T}^{l})^{\top}&\varepsilon^{2}X_{0,T}^{l}\,\hat{Q}\end{bmatrix}\succeq
0$ (40)

showing that $(\hat{Q},\hat{\alpha})$ with $\hat{\alpha}>\alpha$ is a solution
to the first of (16) for the experiment $(x_{0},u_{[0,T-1]})$. Similarly, for
the second of (16) the solution $(Q,\alpha)$ associated with the scaled
experiment $(\varepsilon x_{0},\varepsilon u_{[0,T-1]})$ gives

$\displaystyle\begin{bmatrix}I_{T}&Q\\\ Q^{\top}&\varepsilon
X_{0,T}^{l}\,Q\end{bmatrix}\succeq 0$ (41)

From the relation $(\hat{Q},\hat{\alpha})=(Q/\varepsilon,\alpha)$ we would
thus get

$\displaystyle\begin{bmatrix}I_{T}&Q\\\ Q^{\top}&\varepsilon
X_{0,T}^{l}\,Q\end{bmatrix}=\begin{bmatrix}I_{T}&\varepsilon\hat{Q}\\\
\varepsilon\hat{Q}^{\top}&\varepsilon^{2}X_{0,T}^{l}\,\hat{Q}\end{bmatrix}\succeq
0$ (42)

To summarize, for the linearized system all the experiments $(\varepsilon
x_{0},\varepsilon u_{[0,T-1]})$ result in the same optimal value for $\alpha$.
The proof is concluded noting that for the nonlinear system, the state
matrices $X_{0,T}$ and $X_{1,T}$ associated with the experiment $(\varepsilon
x_{0},\varepsilon u_{[0,T-1]})$ converge to $\varepsilon X^{l}_{0,T}$ and
$\varepsilon X^{l}_{1,T}$, respectively, as $\varepsilon$ converges to zero,
and this implies that the optimal value of $\alpha$ for the nonlinear system
converges to $\overline{\alpha}$. In turn, this implies that for $\varepsilon$
sufficiently small the solution returns a stabilizing controller since: (a)
$\overline{\alpha}$ is fixed; and, as stated in Theorem 4, (b) for any given
$\gamma$, Assumption 2 is satisfied for $\varepsilon$ sufficiently small.
$\blacksquare$

###  3.3 Practical considerations

The foregoing analysis also gives a simple method to assess whether the
solution returns a stabilizing controller. Specifically, if we have some prior
knowledge on $D_{0,T}$, that is on the type of nonlinearity, then we can
directly assess the fulfilment of the condition
$\gamma<\alpha^{2}/(4+2\alpha)$ by computing the smallest value of $\gamma$
that satisfies Assumption 2. If instead we have no prior knowledge on
$D_{0,T}$, one can perform multiple experiments on the system, each one being
an $\varepsilon$-scaled version of the experiment $(\varepsilon
x_{0},\varepsilon u_{[0,T-1]})$. In this case, it follows from the foregoing
analysis that $\alpha$ converges (faster than linearly) to the optimal _fixed_
value $\overline{\alpha}$ associated with solution for the linearized system.

| $\|K-\overline{K}\|$ | $|\alpha-\overline{\alpha}|$ | Stability achieved | $\gamma<\alpha^{2}/(4+2\alpha)$
---|---|---|---|---
|  |  |  | fulfilled
linear system | 0 | 0 | YES | YES
($\varepsilon$-invariant) |  |  |  |
nonlinear system | 0.0082 | 2e-7 | YES | YES
($\varepsilon=1$) |  |  |  |
nonlinear system | 0.0020 | 5e-8 | YES | YES
($\varepsilon=0.5$) |  |  |  |
nonlinear system | 0.0001 | 2e-9 | YES | YES
($\varepsilon=0.1$) |  |  |  |
nonlinear system | 1e-5 | 3e-10 | YES | YES
($\varepsilon=0.01$) |  |  |  |
Table 1: Some tests for the inverted pendulum system.

We illustrate these considerations through a numerical example. Consider the
Euler discretization of an inverted pendulum

$\begin{array}[]{rl}x_{1}(k+1)=&x_{1}(k)+\Delta x_{2}(k)\\\\[14.22636pt]
x_{2}(k+1)=&\displaystyle\frac{\Delta g}{\ell}\sin
x_{1}(k)+\left(1-\frac{\Delta\mu}{m\ell^{2}}\right)x_{2}(k)+\frac{\Delta}{m\ell^{2}}u(k)\\\
\end{array}$

where we simplified the sampled times $k\Delta$ in $k$, with $\Delta$ the
sampling time. In the model $m$ is the mass to be balanced, $\ell$ is the
distance from the base to the center of mass of the balanced body, $\gamma$ is
the coefficient of rotational friction, $g$ is the acceleration due to
gravity. The states $x_{1},x_{2}$ are the angular position and velocity,
respectively, $u$ is the applied torque. The system has an unstable
equilibrium in $(\overline{x},\overline{u})=(0,0)$ corresponding to the
pendulum upright position.

The first experiment is generated starting at the equilibrium ($x_{0}=0$),
with random input sequence within the interval $[-5,5]$, which causes a
displacement of about $28$ degrees from the equilibrium position. All the
other experiments are generated starting at the equilibrium and by scaling the
input sequence of a factor $\varepsilon$. The solution for the linearized
system results in $\overline{\alpha}=0.01136$ and
$\overline{K}=[-56.9163\,\,\,-12.6186]$. (Notice that we allow for such large
input sequences since we start from $x_{0}=0$.)

As the examples shows, $\alpha$ converges to $\overline{\alpha}$ very quickly.
In this respect, an interesting point observed in simulations is that also $K$
converges to $\overline{K}$. In fact, simulations indicate that the solution
$Q$ for the nonlinear system with experiment $(\varepsilon x_{0},\varepsilon
u_{[0,T-1]})$ converges to $\varepsilon\overline{Q}$, where $\overline{Q}$ is
the solution for the linearized system with experiment $(x_{0},u_{[0,T-1]})$.
This property immediately gives that $K=U_{0,1,T}Q(X_{0,T}Q)^{-1}$ converges
to $\varepsilon U_{0,1,T}\varepsilon\overline{Q}(\varepsilon
X^{l}_{0,T}\varepsilon\overline{Q})^{-1}=\overline{K}$.

##  4 Conclusions

We have studied the conditions under which a controller that stabilizes
unknown nonlinear systems in the first approximation directly from
experimental data can be designed. We have shown via an example that these
conditions may not hold even when the data have been obtained applying
persistently exciting inputs. This is in sharp contrast with the results
derived for linear systems. Nevertheless, for those experiments for which the
conditions for the existence of a stabilizer of the nonlinear system do not
hold, we have shown that a suitable scaling of the initial conditions and the
applied input leads to the fulfilment of these conditions, and hence they
allow the design of a data-based stabilizer of the unknown nonlinear system.
We regard these results as a principled method to design stabilizers in the
first approximation and will also have an impact on the nonlocal design of
stabilizers for nonlinear systems, which we are currently investigating.
Important practical aspects have been neglected. The scaling of the
experiments may make the obtained data more sensitive to noise. Although the
approach of De Persis and Tesi (2019a) allows for a robust analysis of the
effect of noise on the controller design, this effect has not been studied in
the current paper and deserves attention.

## References

  * Baggio et al. (2019) Baggio, G., Katewa, V., and Pasqualetti, F. (2019).  Data-driven minimum-energy controls for linear systems.  _IEEE Control Systems Letters_ , 3(3), 589–594.
  * Berberich et al. (2019) Berberich, J., Romer, A., Scherer, C.W., and Allgöwer, F. (2019).  Robust data-driven state-feedback design.  _arXiv preprint arXiv:1909.04314_.
  * Bisoffi et al. (2019) Bisoffi, A., De Persis, C., and Tesi, P. (2019).  Data-based guarantees of set invariance properties.  _arXiv preprint arXiv:1911.12293_.
  * Campi et al. (2002) Campi, M., Lecchini, A., and Savaresi, S. (2002).  Virtual reference feedback tuning: a direct method for the design of feedback controllers.  _Automatica_ , 38(8), 1337–1346.
  * Campi and Savaresi (2006) Campi, M. and Savaresi, S. (2006).  Direct nonlinear control design: the virtual reference feedback tuning (vrft) approach.  _IEEE Transactions on Automatic Control_ , 51(1), 14–27.
  * Chiuso (2016) Chiuso, A. (2016).  Regularization and bayesian learning in dynamical systems: Past, present and future.  _Annual Reviews in Control_ , 41, 24–38.
  * Coulson et al. (2019a) Coulson, J., Lygeros, J., and Dörfler, F. (2019a).  Data-enabled predictive control: In the shallows of the deepc.  In _2019 18th European Control Conference (ECC)_.
  * Coulson et al. (2019b) Coulson, J., Lygeros, J., and Dörfler, F. (2019b).  Regularized and distributionally robust data-enabled predictive control.  In _2019 IEEE 58th Conference on Decision and Control (CDC)_.
  * Dahleh et al. (2003) Dahleh, M., Dahleh, M., and Verghese, G. (2003).  Lectures notes on dynamic systems & control. Available at shorturl.at/cocz3.
  * Dai and Sznaier (2018) Dai, T. and Sznaier, M. (2018).  A moments based approach to designing MIMO data driven controllers for switched systems.  In _57th IEEE Conference on Decision and Control (CDC)_.
  * De Persis and Tesi (2019a) De Persis, C. and Tesi, P. (2019a).  Formulas for data-driven control: Stabilization, optimality and robustness.  _IEEE Transactions on Automatic Control_ , 65(3), 909–924.
  * De Persis and Tesi (2019b) De Persis, C. and Tesi, P. (2019b).  On persistency of excitation and formulas for data-driven control.  In _IEEE 58th Conference on Decision and Control (CDC)_.
  * Formentin et al. (2013) Formentin, S., Karimi, A., and Savaresi, S. (2013).  Optimal input design for direct data-driven tuning of model-reference controllers.  _Automatica_ , 49(6), 1874–1882.
  * Gonçalves da Silva et al. (2019) Gonçalves da Silva, R., Bazanella, A., Lorenzini, C., and Campestrini, L. (2019).  Data-driven LQR control design.  _IEEE Control Systems Letters_ , 3(1), 180–185.
  * Hjalmarsson and Gevers (1998) Hjalmarsson, H. and Gevers, M. (1998).  Iterative feedback tuning: Theory and applications.  _IEEE Control Systems Magazine_ , 18(4), 26–41.
  * Huang et al. (2019) Huang, L., Coulson, J., Lygeros, J., and Dörfler, F. (2019).  Data-enabled predictive control for grid-connected power converters.  In _2019 IEEE 58th Conference on Decision and Control (CDC)_.
  * Karimi et al. (2007) Karimi, A., Van Heusden, K., and Bonvin, D. (2007).  Non-iterative data-driven controller tuning using the correlation approach.  In _2007 European Control Conference (ECC)_.
  * Markovsky and Rapisarda (2008) Markovsky, I. and Rapisarda, P. (2008).  Data-driven simulation and control.  _International Journal of Control_ , 81(12), 1946–1959.
  * Markovsky et al. (2005) Markovsky, I., Willems, J., Rapisarda, P., and De Moor, B. (2005).  Algorithms for deterministic balanced subspace identification.  _Automatica_ , 41(5), 755–766.
  * Novara et al. (2016) Novara, C., Formentin, S., Savaresi, S., and Milanese, M. (2016).  Data-driven design of two degree-of-freedom nonlinear controllers: The D2-IBC approach.  _Automatica_ , 72, 19–27.
  * Recht (2019) Recht, B. (2019).  A tour of reinforcement learning: The view from continuous control.  _Annual Review of Control, Robotics, and Autonomous Systems_ , 2, 253–279.
  * Tanaskovic et al. (2017) Tanaskovic, M., Fagiano, L., Novara, C., and Morari, M. (2017).  Data-driven control of nonlinear systems: An on-line direct approach.  _Automatica_ , 75, 1–10.
  * van Waarde et al. (2020) van Waarde, H.J., Eising, J., Trentelman, H.L., and Camlibel, M.K. (2020).  Data informativity: a new perspective on data-driven analysis and control.  _IEEE Transactions on Automatic Control_ , 65(11), 4753–4768.
  * Wabersich and Zeilinger (2018) Wabersich, K.P. and Zeilinger, M.N. (2018).  Scalable synthesis of safety certificates from data with application to learning-based control.  In _2018 European Control Conference (ECC)_.
  * Willems et al. (2005) Willems, J., Rapisarda, P., Markovsky, I., and De Moor, B. (2005).  A note on persistency of excitation.  _Systems & Control Letters_, 54(4), 325–329.


supersonic relative velocities ($v_{\rm rel}/c_{\rm s,hot}$). Numerical
simulations of such high Mach number cloud-wind interactions are challenging,
so it remains unclear how the mass, momentum, and energy transfer depends on
the Mach number. There is, however, an indication that the characteristic
cloud destruction time scale should be lengthened by a factor of
$\sqrt{1+\mathcal{M}}$ (Scannapieco & Brüggen, 2015).

_Distributed cloud sweep-up in the subsonic region_ : The simple volumetric
injection of hot material at small radii in our model drives the wind which
immediately runs into a thin shell of cold clouds outside of the injection
region. This abrupt introduction of the cold clouds into the hot wind neglects
the fact that many (if not all) of the cold clouds are initially co-spatial
with the supernova remnants that are powering the wind. A more physical
approach may be to introduce the clouds distributed uniformly within the same
(subsonic) region that the hot material is added.

_Unexplored dependencies_ : In addition to the long list of missing
ingredients and model extensions there are also parameter dependencies that we
have yet to fully explore. In particular, the structure of multiphase galactic
winds has interesting albeit less significant dependence on $v_{\rm c}$, $Z$,
and $r_{\star}$, that we have not presented for the sake of brevity. We
encourage interested readers to explore these dependencies using our code that
is available here: https://github.com/dfielding14/MultiphaseGalacticWind/.

##  6 Summary

This work presents a new analytic model that describes the steady state
structure of a spherically symmetric multiphase galactic wind that is
comprised of cold clouds and a volume-filling hot component. The novel aspect
of our model, which is also the fundamental driver of much of the interesting
behavior of the winds, is our treatment of the exchange of mass, momentum, and
energy between the hot and cold phase. This exchange is regulated in the
turbulent radiative mixing layers that develop between the phases. The
transfer rates are set by the competition between turbulence, which shreds the
clouds, and cooling, which grows the clouds, with the crucial parameter being
the ratio of the mixing time to the cooling time $\xi=t_{\rm mix}/t_{\rm
cool}$. When $\xi>1$ the cloud is in the growth regime, and when $\xi<1$ the
cloud is in the shredding regime. Simulations of turbulent radiative mixing
layers and cloud-wind interactions have, recently, illuminated the nature of
these competing processes and enabled us to write down intuitive analytic
approximations that are ideally suited for use in a galactic wind model. By
including this essential and inherently small scale phenomena into a global
context we have been able to explore not only the impact of the hot wind on
the cold clouds, but the previously under-appreciated back reaction of the
cold clouds on the hot wind itself. A schematic summary of our model is shown
in Fig. 1.

In § 2 we lay the foundation for our model by exploring the general behavior
of a hot wind under the influence of arbitrary sources/sinks of mass,
momentum, and energy. Limiting cases reproduce well studied cases of galactic
wind solutions, driven by centrally injected mass and energy with or without
gravity and cooling (see § 2.1). We then refine this generic hot wind radial
structure formulation to account for the specific case in which material can
be either gained from, or lost to, another phase that is moving with a
velocity relative to the hot phase (see § 2.3). In addition to the change in
mass, momentum, and energy from this transfer, our formulation also accounts
for additional forces, such as drag/ram pressure, and additional energy source
terms, such as radiative cooling (see Eq. 15). This leads us to the
fundamental equations that govern the radial evolution of the hot phase of a
galactic wind (see Eq. 16 to Eq. 19). Although these equations for the hot
phase velocity $v$, density $\rho$, pressure $P$, and entropy $K$, have up to
ten separate terms, each term has a well-defined physical interpretation,
which we use to understand what is driving the properties of the wind (see
Fig. 6).

The source terms for the hot phase come from the interaction of the hot phase
with the cold clouds that are embedded within it. In § 3, we present a
framework for modeling the cloud-wind interaction that, although simple,
encapsulates the essential behaviors. Specifically, our framework describes
the mass transfer, accelerations, and energetics of this interaction. We
present a physically motivated parameterization of the mass evolution of an
embedded cold cloud that is inspired by and consistent with extensive recent
work on the topic (see § 3.2 and Fig. 3). In developing this parameterization
we stress where our specific choices are underconstrained by simulations or
oversimplified and can be refined in future work.

Finally, we present the full model for the coevolution of the cold clouds and
the hot wind in § 4. This model allows us to self-consistently describe the
radial structure of the hot wind and the properties of the cold clouds for a
given set of wind launching parameters. The primary parameters that determine
the wind properties are the SFR, hot phase mass loading $\eta_{\rm M,hot}$,
hot phase energy loading $\eta_{\rm E}$, cold phase mass loading $\eta_{\rm
M,cold}$, and the initial mass of cold clouds $M_{\rm cl}$.888Additional
parameters include the size of the star forming region $r_{\star}$, hot phase
metallicity, and cold phase metallicity, we, however, keep these fixed at 300
pc, 2 $Z_{\odot}$, and $Z_{\odot}$, respectively. We demonstrate how the
presence of cold clouds can dramatically impact the structure of the hot wind
and how the clouds themselves evolve as they flow out into the surrounding
medium. Some of the notable features of the multiphase wind solutions are as
follows.

  1. i.

We find that the multiphase model presented here naturally leads to galactic
winds in which the energy from supernovae is contained within a low-density,
high-specific-energy fluid, while cold dense clouds are accelerated and
entrained by the turbulent radiative mixing layers that form at their
boundaries. This allows the winds to be effective agents of feedback in the
circumgalactic medium, providing a mechanism by which star formation in the
galactic disk can slow cooling and accretion as part of a global preventive
self-regulation mechanism. In addition, the winds can transport cold clouds to
large radii, in agreement with observations999We note that the amount of mass
in the hot phase is relatively small, meaning that it is very challenging for
cosmological galaxy-scale simulations (almost all of which adopt fixed mass
resolution) to resolve. Such simulations tend to mix the two phases, resulting
in a wind which cools much more rapidly than the multiphase structure
envisioned here, requiring the winds to be essentially momentum-driven, and
therefore highly mass loaded. We demonstrate in the current paper that is not
required and that, in real galactic systems, the high-specific energy hot wind
does not (necessarily) lose significant amounts of energy so that low-mass
loaded winds can effectively transmit energy into the circumgalactic medium..

  2. ii.

The bi-directional interaction between the hot and cold phases is of critical
importance in the model. The acceleration of the cold clouds, which is almost
always dominated by the accretion of cooled wind material onto the clouds
during the turbulent radiative mixing layer entrainment process (see Fig. 6),
comes at the expense of the hot phase velocity, which can be reduced by a
factor $>3$ relative to an otherwise identical wind without the cold clouds.
Furthermore, the fact that cold clouds are predominantly accelerated by the
exchange of mass with the hot phase explains the close correspondence between
the degree of metal mixing and the clouds velocity (see Fig. 4 and Fig. 5).

  3. iii.

The thermalization of the relative kinetic energy as cloud material is mixed
into the hot wind leads to a significant heating of the wind. When cooling is
unimportant in the hot phase, which happens when the hot phase mass loading is
low ($\eta_{\rm M,hot}\sim 0.1$), this heating increases the hot wind entropy,
and results in a flatter temperature profile (see Fig. 4 and Fig. 6). Flatter
temperature profiles have been found in high resolution galactic wind
simulations (Fielding et al., 2017a; Schneider et al., 2020), as well as in
observations (Lopez et al., 2020).

  4. iv.

In galactic winds with a low hot phase density, which occurs when either
$\eta_{\rm M,hot}$ and/or SFR is small, cold clouds tend to lose mass (see
Fig. 4). This is true even with the most massive cold clouds $M_{\rm
cl}=10^{6}\ M_{\odot}$, which start out in the growing regime ($\xi>1$) but
quickly transition to being shredded ($\xi<1$). This is because the cooling
rate in the turbulent radiative mixing layers drops precipitously as the
pressure drops. As the clouds are shredded the mass flux in the cold phase
drops while the hot phase becomes increasingly mass loaded. In cases in which
there is a large initial cold mass loading factor and the clouds are rapidly
shredded ($\xi\lesssim 0.1$) the rapid addition of material into the hot phase
can cause the wind to _fail_ by driving the Mach number down below unity.
Winds with small $\eta_{\rm M,hot}$ and/or low SFRs are able to sustain higher
total mass loading ($\eta_{\rm M,tot}=\eta_{\rm M,hot}+\eta_{\rm M,cold}$),
but are likely to become increasingly dominated by the hot phase at larger
distances (see Fig. 7, Fig. 8, and Fig. 9).

  5. v.

On the other hand, in galactic winds with a higher hot phase density, which
occurs when either $\eta_{\rm M,hot}$ and/or the SFR is large, cold clouds
tend to gain mass because of the short cooling times in their turbulent
radiative mixing layers (see Fig. 5). As a result of cloud growth, the mass
flux in the hot phase drops and the density profile steepens. In these rapid
cloud-growth systems the hot phase energy flux is dramatically reduced. In
some cases the hot wind loses all of its thermal energy and cools down to the
clouds temperature ($10^{4}$ K) through the combination of cooling in the
mixing layer and in the volume filling phase. These winds have in effect been
_poisoned_ by sweeping up too many cold clouds. This processes limits the
maximum cold phase mass flux, and as a result winds with large $\eta_{\rm
M,hot}$ and/or high SFRs are not able to sustain as high of a total mass
loading factor ($\eta_{\rm M,tot}=\eta_{\rm M,hot}+\eta_{\rm M,cold}$), but
are likely to remain multiphase (or even predominantly cold) to large
distances (see Fig. 7, Fig. 8, and Fig. 9).

The model we have presented provides a flexible, intuitive, and self-
consistent framework to understand the nature of galactic winds. Our analytic
approach is an ideal complement to simulations because it reproduces many of
the key processes in both small and large scale simulations while being easily
extensible and incurring essentially no computational cost, making it well-
suited for comparison to diverse observations. This first incarnation of the
model sacrifices some physical realism for the sake of simplicity, but there
is a clear path forward to include extensions. These extensions could include
improved parameterizations, the inclusion of additional processes, a full
accounting for the phase structure necessary for detailed observational
comparisons, and the model’s implementation as a cosmological simulation
subgrid model. We suspect that these extensions will only strengthen the key
finding of this work that accounting for the full galactic wind phase
structure is crucial for understanding the flow of mass, momentum, and energy
into and out of galaxies which regulates their formation over cosmic time.

We have derived pleasure and profit from discussing the ideas in this paper
with many colleagues. In particular we would like to thank Matthew Smith, Eve
Ostriker, Alberto Bolatto, Eliot Quataert, Chang-Goo Kim, John Forbes, Matthew
Abruzzo, Anthony Chow, Zirui Chen, Evan Schneider, Max Gronke, Peng Oh, and
the SMAUG team for many useful discussions and suggestions. We thank Lucy
Reading-Ikkanda for help making Fig. 1. DBF is supported by the Simons
Foundation through the Flatiron Institute. GLB acknowledges financial support
from the Simons Foundation, NSF (grant OAC-1835509), and NASA (grant
NNX15AB20G). We thank KITP for hosting the Fundamentals of Gaseous Halos
workshop (supported by the NSF under Grant No. NSF PHY-1748958) during which
some of the key ideas of this work were formulated. _Code availability:_ A
basic implementation of the code used to solve for the steady state structure
of a multiphase wind that was used to create all of the figures in this paper
is available here: https://github.com/dfielding14/MultiphaseGalacticWind/.
_Software:_ matplotlib (Hunter, 2007), scipy (Virtanen et al., 2020), numpy
(Harris et al., 2020), ipython (Perez & Granger, 2007), CMasher (van der
Velden, 2020).

## References

  * Abruzzo et al. (2021) Abruzzo, M. W., Bryan, G. L., & Fielding, D. B. 2021, arXiv e-prints, arXiv:2101.10344.  https://arxiv.org/abs/2101.10344
  * Afruni et al. (2019) Afruni, A., Fraternali, F., & Pezzulli, G. 2019, A&A, 625, A11, doi: 10.1051/0004-6361/201835002
  * Afruni et al. (2021) —. 2021, MNRAS, 501, 5575, doi: 10.1093/mnras/staa3759
  * Anglés-Alcázar et al. (2017) Anglés-Alcázar, D., Faucher-Giguère, C.-A., Kereš, D., et al. 2017, MNRAS, 470, 4698, doi: 10.1093/mnras/stx1517
  * Armillotta et al. (2016) Armillotta, L., Fraternali, F., & Marinacci, F. 2016, MNRAS, 462, 4157, doi: 10.1093/mnras/stw1930
  * Ashley et al. (2020) Ashley, T., Fox, A. J., Jenkins, E. B., et al. 2020, ApJ, 898, 128, doi: 10.3847/1538-4357/ab9ff8
  * Banda-Barragán et al. (2020) Banda-Barragán, W., Brüggen, M., Heesen, V., et al. 2020, arXiv e-prints, arXiv:2011.05240.  https://arxiv.org/abs/2011.05240
  * Begelman & Fabian (1990) Begelman, M. C., & Fabian, A. C. 1990, MNRAS, 244, 26P
  * Begelman & McKee (1990) Begelman, M. C., & McKee, C. F. 1990, ApJ, 358, 375, doi: 10.1086/168994
  * Behroozi et al. (2013) Behroozi, P. S., Wechsler, R. H., & Conroy, C. 2013, ApJ, 770, 57, doi: 10.1088/0004-637X/770/1/57
  * Benson et al. (2003) Benson, A. J., Bower, R. G., Frenk, C. S., et al. 2003, ApJ, 599, 38, doi: 10.1086/379160
  * Bolatto et al. (2013) Bolatto, A. D., Warren, S. R., Leroy, A. K., et al. 2013, Nature, 499, 450, doi: 10.1038/nature12351
  * Brüggen & Scannapieco (2016) Brüggen, M., & Scannapieco, E. 2016, ApJ, 822, 31, doi: 10.3847/0004-637X/822/1/31
  * Brüggen & Scannapieco (2020) —. 2020, ApJ, 905, 19, doi: 10.3847/1538-4357/abc00f
  * Burchett et al. (2021) Burchett, J. N., Rubin, K. H. R., Prochaska, J. X., et al. 2021, ApJ, 909, 151, doi: 10.3847/1538-4357/abd4e0
  * Bustard et al. (2016) Bustard, C., Zweibel, E. G., & D’Onghia, E. 2016, ApJ, 819, 29, doi: 10.3847/0004-637X/819/1/29
  * Butsky & Quinn (2018) Butsky, I., & Quinn, T. R. 2018, ArXiv e-prints.  https://arxiv.org/abs/1803.06345
  * Chen et al. (2010a) Chen, H.-W., Wild, V., Tinker, J. L., et al. 2010a, ApJ, 724, L176, doi: 10.1088/2041-8205/724/2/L176
  * Chen et al. (2010b) Chen, Y.-M., Tremonti, C. A., Heckman, T. M., et al. 2010b, AJ, 140, 445, doi: 10.1088/0004-6256/140/2/445
  * Chevalier & Clegg (1985) Chevalier, R. A., & Clegg, A. W. 1985, Nature, 317, 44, doi: 10.1038/317044a0
  * Chisholm et al. (2018) Chisholm, J., Bordoloi, R., Rigby, J. R., & Bayliss, M. 2018, MNRAS, 474, 1688, doi: 10.1093/mnras/stx2848
  * Chisholm et al. (2017) Chisholm, J., Tremonti, C. A., Leitherer, C., & Chen, Y. 2017, MNRAS, 469, 4831, doi: 10.1093/mnras/stx1164
  * Chisholm et al. (2015) Chisholm, J., Tremonti, C. A., Leitherer, C., et al. 2015, ApJ, 811, 149, doi: 10.1088/0004-637X/811/2/149
  * Cicone et al. (2014) Cicone, C., Maiolino, R., Sturm, E., et al. 2014, A&A, 562, A21, doi: 10.1051/0004-6361/201322464
  * Cooper et al. (2009) Cooper, J. L., Bicknell, G. V., Sutherland, R. S., & Bland-Hawthorn, J. 2009, ApJ, 703, 330, doi: 10.1088/0004-637X/703/1/330
  * Cottle et al. (2020) Cottle, J., Scannapieco, E., Brüggen, M., Banda-Barragán, W., & Federrath, C. 2020, ApJ, 892, 59, doi: 10.3847/1538-4357/ab76d1
  * Cowie & McKee (1977) Cowie, L. L., & McKee, C. F. 1977, ApJ, 211, 135, doi: 10.1086/154911
  * Cowie et al. (1981) Cowie, L. L., McKee, C. F., & Ostriker, J. P. 1981, ApJ, 247, 908, doi: 10.1086/159100
  * Crain et al. (2015) Crain, R. A., Schaye, J., Bower, R. G., et al. 2015, MNRAS, 450, 1937, doi: 10.1093/mnras/stv725
  * Creasey et al. (2013) Creasey, P., Theuns, T., & Bower, R. G. 2013, MNRAS, 429, 1922, doi: 10.1093/mnras/sts439
  * Davé et al. (2019) Davé, R., Anglés-Alcázar, D., Narayanan, D., et al. 2019, MNRAS, 486, 2827, doi: 10.1093/mnras/stz937
  * de Avillez (2000) de Avillez, M. A. 2000, MNRAS, 315, 479, doi: 10.1046/j.1365-8711.2000.03464.x
  * Dekel & Silk (1986) Dekel, A., & Silk, J. 1986, ApJ, 303, 39, doi: 10.1086/164050
  * Drake et al. (2020) Drake, J. F., Pfrommer, C., Reynolds, C. S., et al. 2020, arXiv e-prints, arXiv:2007.07931.  https://arxiv.org/abs/2007.07931
  * Dursi & Pfrommer (2008) Dursi, L. J., & Pfrommer, C. 2008, ApJ, 677, 993, doi: 10.1086/529371
  * Emerick et al. (2018) Emerick, A., Bryan, G. L., & Mac Low, M.-M. 2018, ApJ, 865, L22, doi: 10.3847/2041-8213/aae315
  * Farber & Gronke (2021) Farber, R. J., & Gronke, M. 2021, arXiv e-prints, arXiv:2107.07991.  https://arxiv.org/abs/2107.07991
  * Fielding et al. (2018) Fielding, D., Quataert, E., & Martizzi, D. 2018, MNRAS, 481, 3325, doi: 10.1093/mnras/sty2466
  * Fielding et al. (2017a) Fielding, D., Quataert, E., Martizzi, D., & Faucher-Giguère, C.-A. 2017a, MNRAS, 470, L39, doi: 10.1093/mnrasl/slx072
  * Fielding et al. (2017b) Fielding, D., Quataert, E., McCourt, M., & Thompson, T. A. 2017b, MNRAS, 466, 3810, doi: 10.1093/mnras/stw3326
  * Fielding et al. (2020a) Fielding, D. B., Ostriker, E. C., Bryan, G. L., & Jermyn, A. S. 2020a, ApJ, 894, L24, doi: 10.3847/2041-8213/ab8d2c
  * Fielding et al. (2020b) Fielding, D. B., Tonnesen, S., DeFelippis, D., et al. 2020b, ApJ, 903, 32, doi: 10.3847/1538-4357/abbc6d
  * Gatto et al. (2017) Gatto, A., Walch, S., Naab, T., et al. 2017, MNRAS, 466, 1903, doi: 10.1093/mnras/stw3209
  * Girichidis et al. (2021) Girichidis, P., Naab, T., Walch, S., & Berlok, T. 2021, MNRAS, doi: 10.1093/mnras/stab1203
  * Girichidis et al. (2016) Girichidis, P., Naab, T., Walch, S., et al. 2016, ApJ, 816, L19, doi: 10.3847/2041-8205/816/2/L19
  * Gray et al. (2019) Gray, W. J., Scannapieco, E., & Lehnert, M. D. 2019, ApJ, 875, 110, doi: 10.3847/1538-4357/ab1004
  * Gronke & Oh (2018) Gronke, M., & Oh, S. P. 2018, ArXiv e-prints.  https://arxiv.org/abs/1806.02728
  * Gronke & Oh (2020) —. 2020, MNRAS, 492, 1970, doi: 10.1093/mnras/stz3332
  * Gronke et al. (2021) Gronke, M., Oh, S. P., Ji, S., & Norman, C. 2021, arXiv e-prints, arXiv:2107.13012.  https://arxiv.org/abs/2107.13012
  * Grønnow et al. (2018) Grønnow, A., Tepper-García, T., & Bland-Hawthorn, J. 2018, ArXiv e-prints.  https://arxiv.org/abs/1805.03903
  * Gutcke et al. (2017) Gutcke, T. A., Stinson, G. S., Macciò, A. V., Wang, L., & Dutton, A. A. 2017, MNRAS, 464, 2796, doi: 10.1093/mnras/stw2539
  * Hafen et al. (2019) Hafen, Z., Faucher-Giguère, C.-A., Anglés-Alcázar, D., et al. 2019, MNRAS, 488, 1248, doi: 10.1093/mnras/stz1773
  * Hafen et al. (2020) —. 2020, MNRAS, 494, 3581, doi: 10.1093/mnras/staa902
  * Harris et al. (2020) Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Nature, 585, 357, doi: 10.1038/s41586-020-2649-2
  * Hayes et al. (2016) Hayes, M., Melinder, J., Östlin, G., et al. 2016, ApJ, 828, 49, doi: 10.3847/0004-637X/828/1/49
  * Heckman et al. (2017) Heckman, T., Borthakur, S., Wild, V., Schiminovich, D., & Bordoloi, R. 2017, ApJ, 846, 151, doi: 10.3847/1538-4357/aa80dc
  * Heckman et al. (2015) Heckman, T. M., Alexandroff, R. M., Borthakur, S., Overzier, R., & Leitherer, C. 2015, ApJ, 809, 147, doi: 10.1088/0004-637X/809/2/147
  * Heckman et al. (1990) Heckman, T. M., Armus, L., & Miley, G. K. 1990, ApJS, 74, 833, doi: 10.1086/191522
  * Heckman et al. (2000) Heckman, T. M., Lehnert, M. D., Strickland, D. K., & Armus, L. 2000, ApJs, 129, 493, doi: 10.1086/313421
  * Hill et al. (2012) Hill, A. S., Joung, M. R., Mac Low, M.-M., et al. 2012, ApJ, 750, 104, doi: 10.1088/0004-637X/750/2/104
  * Hodges-Kluck et al. (2020) Hodges-Kluck, E. J., Yukita, M., Tanner, R., et al. 2020, ApJ, 903, 35, doi: 10.3847/1538-4357/abb884
  * Hopkins (2012) Hopkins, P. F. 2012, MNRAS, 423, 2016, doi: 10.1111/j.1365-2966.2012.20730.x
  * Hopkins et al. (2021a) Hopkins, P. F., Chan, T. K., Ji, S., et al. 2021a, MNRAS, 501, 3640, doi: 10.1093/mnras/staa3690
  * Hopkins et al. (2021b) Hopkins, P. F., Chan, T. K., Squire, J., et al. 2021b, MNRAS, 501, 3663, doi: 10.1093/mnras/staa3692
  * Hopkins et al. (2014) Hopkins, P. F., Kereš, D., Oñorbe, J., et al. 2014, MNRAS, 445, 581, doi: 10.1093/mnras/stu1738
  * Hu (2019) Hu, C.-Y. 2019, MNRAS, 483, 3363, doi: 10.1093/mnras/sty3252
  * Huang et al. (2021) Huang, S., Katz, N., Cottle, J., et al. 2021, arXiv e-prints, arXiv:2106.01511.  https://arxiv.org/abs/2106.01511
  * Huang et al. (2020) Huang, S., Katz, N., Scannapieco, E., et al. 2020, MNRAS, 497, 2586, doi: 10.1093/mnras/staa1978
  * Hummels et al. (2013) Hummels, C. B., Bryan, G. L., Smith, B. D., & Turk, M. J. 2013, MNRAS, 430, 1548, doi: 10.1093/mnras/sts702
  * Hunter (2007) Hunter, J. D. 2007, Computing in Science and Engineering, 9, 90, doi: 10.1109/MCSE.2007.55
  * Jana et al. (2020) Jana, R., Gupta, S., & Nath, B. B. 2020, MNRAS, 497, 2623, doi: 10.1093/mnras/staa2025
  * Ji et al. (2019) Ji, S., Oh, S. P., & Masterson, P. 2019, MNRAS, 487, 737, doi: 10.1093/mnras/stz1248
  * Johnson et al. (2017) Johnson, S. D., Chen, H.-W., Mulchaey, J. S., Schaye, J., & Straka, L. A. 2017, ApJ, 850, L10, doi: 10.3847/2041-8213/aa9370
  * Joung & Mac Low (2006) Joung, M. K. R., & Mac Low, M.-M. 2006, ApJ, 653, 1266, doi: 10.1086/508795
  * Kacprzak et al. (2015) Kacprzak, G. G., Muzahid, S., Churchill, C. W., Nielsen, N. M., & Charlton, J. C. 2015, ApJ, 815, 22, doi: 10.1088/0004-637X/815/1/22
  * Kanjilal et al. (2021) Kanjilal, V., Dutta, A., & Sharma, P. 2021, MNRAS, 501, 1143, doi: 10.1093/mnras/staa3610
  * Kereš et al. (2005) Kereš, D., Katz, N., Weinberg, D. H., & Davé, R. 2005, MNRAS, 363, 2, doi: 10.1111/j.1365-2966.2005.09451.x
  * Kim & Ostriker (2018) Kim, C.-G., & Ostriker, E. C. 2018, ApJ, 853, 173, doi: 10.3847/1538-4357/aaa5ff
  * Kim et al. (2017) Kim, C.-G., Ostriker, E. C., & Raileanu, R. 2017, ApJ, 834, 25, doi: 10.3847/1538-4357/834/1/25
  * Kim et al. (2020a) Kim, C.-G., Ostriker, E. C., Somerville, R. S., et al. 2020a, ApJ, 900, 61, doi: 10.3847/1538-4357/aba962
  * Kim et al. (2020b) Kim, C.-G., Ostriker, E. C., Fielding, D. B., et al. 2020b, ApJ, 903, L34, doi: 10.3847/2041-8213/abc252
  * Klein et al. (1994) Klein, R. I., McKee, C. F., & Colella, P. 1994, ApJ, 420, 213, doi: 10.1086/173554
  * Komarov et al. (2018) Komarov, S., Schekochihin, A. A., Churazov, E., & Spitkovsky, A. 2018, Journal of Plasma Physics, 84, 905840305, doi: 10.1017/S0022377818000399
  * Koo & McKee (1992) Koo, B.-C., & McKee, C. F. 1992, ApJ, 388, 93, doi: 10.1086/171132
  * Krieger et al. (2021) Krieger, N., Walter, F., Bolatto, A. D., et al. 2021, ApJ, 915, L3, doi: 10.3847/2041-8213/ac01e9
  * Lehnert & Heckman (1996) Lehnert, M. D., & Heckman, T. M. 1996, ApJ, 462, 651, doi: 10.1086/177180
  * Leroy et al. (2015) Leroy, A. K., Walter, F., Martini, P., et al. 2015, ApJ, 814, 83, doi: 10.1088/0004-637X/814/2/83
  * Li & Bryan (2020) Li, M., & Bryan, G. L. 2020, ApJ, 890, L30, doi: 10.3847/2041-8213/ab7304
  * Li et al. (2017) Li, M., Bryan, G. L., & Ostriker, J. P. 2017, ApJ, 841, 101, doi: 10.3847/1538-4357/aa7263
  * Li et al. (2020) Li, Z., Hopkins, P. F., Squire, J., & Hummels, C. 2020, MNRAS, 492, 1841, doi: 10.1093/mnras/stz3567
  * Lochhaas et al. (2020) Lochhaas, C., Bryan, G. L., Li, Y., Li, M., & Fielding, D. 2020, MNRAS, 493, 1461, doi: 10.1093/mnras/staa358
  * Lochhaas et al. (2018) Lochhaas, C., Thompson, T. A., Quataert, E., & Weinberg, D. H. 2018, ArXiv e-prints.  https://arxiv.org/abs/1804.00741
  * Lopez et al. (2020) Lopez, L. A., Mathur, S., Nguyen, D. D., Thompson, T. A., & Olivier, G. M. 2020, ApJ, 904, 152, doi: 10.3847/1538-4357/abc010
  * Mac Low & McCray (1988) Mac Low, M.-M., & McCray, R. 1988, ApJ, 324, 776, doi: 10.1086/165936
  * Mac Low et al. (1994) Mac Low, M.-M., McKee, C. F., Klein, R. I., Stone, J. M., & Norman, M. L. 1994, ApJ, 433, 757, doi: 10.1086/174685
  * Mandelker et al. (2020) Mandelker, N., Nagai, D., Aung, H., et al. 2020, MNRAS, 494, 2641, doi: 10.1093/mnras/staa812
  * Mao & Ostriker (2018) Mao, S. A., & Ostriker, E. C. 2018, ApJ, 854, 89, doi: 10.3847/1538-4357/aaa88e
  * Marcolini et al. (2005) Marcolini, A., Strickland, D. K., D’Ercole, A., Heckman, T. M., & Hoopes, C. G. 2005, MNRAS, 362, 626, doi: 10.1111/j.1365-2966.2005.09343.x
  * Marinacci et al. (2010) Marinacci, F., Binney, J., Fraternali, F., et al. 2010, MNRAS, 404, 1464, doi: 10.1111/j.1365-2966.2010.16352.x
  * Martin (1999) Martin, C. L. 1999, ApJ, 513, 156, doi: 10.1086/306863
  * Martin (2005) —. 2005, ApJ, 621, 227, doi: 10.1086/427277
  * Martin & Bouché (2009) Martin, C. L., & Bouché, N. 2009, ApJ, 703, 1394, doi: 10.1088/0004-637X/703/2/1394
  * Martini et al. (2018) Martini, P., Leroy, A. K., Mangum, J. G., et al. 2018, ApJ, 856, 61, doi: 10.3847/1538-4357/aab08e
  * Martizzi et al. (2016) Martizzi, D., Fielding, D., Faucher-Giguère, C.-A., & Quataert, E. 2016, MNRAS, 459, 2311, doi: 10.1093/mnras/stw745
  * Mathews & Baker (1971) Mathews, W. G., & Baker, J. C. 1971, ApJ, 170, 241, doi: 10.1086/151208
  * McCourt et al. (2015) McCourt, M., O’Leary, R. M., Madigan, A.-M., & Quataert, E. 2015, MNRAS, 449, 2, doi: 10.1093/mnras/stv355
  * McKeith et al. (1995) McKeith, C. D., Greve, A., Downes, D., & Prada, F. 1995, A&A, 293, 703
  * Meinecke et al. (2021) Meinecke, J., Tzeferacos, P., Ross, J. S., et al. 2021, arXiv e-prints, arXiv:2105.08461.  https://arxiv.org/abs/2105.08461
  * Mellema et al. (2003) Mellema, G., Kurk, J., & Röttgering, H. 2003, New A Rev., 47, 243, doi: 10.1016/S1387-6473(03)00033-2
  * Melso et al. (2019) Melso, N., Bryan, G. L., & Li, M. 2019, ApJ, 872, 47, doi: 10.3847/1538-4357/aafaf5
  * Mitchell et al. (2020) Mitchell, P. D., Schaye, J., Bower, R. G., & Crain, R. A. 2020, MNRAS, 494, 3971, doi: 10.1093/mnras/staa938
  * Morita (1982) Morita, K. 1982, PASJ, 34, 65
  * Muratov et al. (2015) Muratov, A. L., Kereš, D., Faucher-Giguère, C.-A., et al. 2015, MNRAS, 454, 2691, doi: 10.1093/mnras/stv2126
  * Naab & Ostriker (2017) Naab, T., & Ostriker, J. P. 2017, ARA&A, 55, 59, doi: 10.1146/annurev-astro-081913-040019
  * Nelson et al. (2021) Nelson, D., Byrohl, C., Peroux, C., Rubin, K. H. R., & Burchett, J. N. 2021, arXiv e-prints, arXiv:2106.09023.  https://arxiv.org/abs/2106.09023
  * Nelson et al. (2013) Nelson, D., Vogelsberger, M., Genel, S., et al. 2013, MNRAS, 429, 3353, doi: 10.1093/mnras/sts595
  * Nelson et al. (2019) Nelson, D., Pillepich, A., Springel, V., et al. 2019, MNRAS, 2010, doi: 10.1093/mnras/stz2306
  * Nelson et al. (2020) Nelson, D., Sharma, P., Pillepich, A., et al. 2020, arXiv e-prints, arXiv:2005.09654.  https://arxiv.org/abs/2005.09654
  * Newman et al. (2012) Newman, S. F., Genzel, R., Förster-Schreiber, N. M., et al. 2012, ApJ, 761, 43, doi: 10.1088/0004-637X/761/1/43
  * Nielsen et al. (2015) Nielsen, N. M., Churchill, C. W., Kacprzak, G. G., Murphy, M. T., & Evans, J. L. 2015, ApJ, 812, 83, doi: 10.1088/0004-637X/812/1/83
  * Oppenheimer & Davé (2008) Oppenheimer, B. D., & Davé, R. 2008, MNRAS, 387, 577, doi: 10.1111/j.1365-2966.2008.13280.x
  * Orlando et al. (2005) Orlando, S., Peres, G., Reale, F., et al. 2005, A&A, 444, 505, doi: 10.1051/0004-6361:20052896
  * Osterbrock (1989) Osterbrock, D. E. 1989, Astrophysics of gaseous nebulae and active galactic nuclei
  * Pandya et al. (2021) Pandya, V., Fielding, D., Anglés-Alcázar, D., et al. 2021, arXiv e-prints, arXiv:2103.06891.  https://arxiv.org/abs/2103.06891
  * Perez & Granger (2007) Perez, F., & Granger, B. E. 2007, Computing in Science and Engineering, 9, 21, doi: 10.1109/MCSE.2007.53
  * Péroux et al. (2020) Péroux, C., Nelson, D., van de Voort, F., et al. 2020, MNRAS, 499, 2462, doi: 10.1093/mnras/staa2888
  * Pettini et al. (2001) Pettini, M., Shapley, A. E., Steidel, C. C., et al. 2001, ApJ, 554, 981, doi: 10.1086/321403
  * Pillepich et al. (2018) Pillepich, A., Nelson, D., Hernquist, L., et al. 2018, MNRAS, 475, 648, doi: 10.1093/mnras/stx3112
  * Prochaska et al. (2011) Prochaska, J. X., Weiner, B., Chen, H.-W., Mulchaey, J., & Cooksey, K. 2011, ApJ, 740, 91, doi: 10.1088/0004-637X/740/2/91
  * Prochaska et al. (2017) Prochaska, J. X., Werk, J. K., Worseck, G., et al. 2017, ApJ, 837, 169, doi: 10.3847/1538-4357/aa6007
  * Ptak et al. (1997) Ptak, A., Serlemitsos, P., Yaqoob, T., Mushotzky, R., & Tsuru, T. 1997, AJ, 113, 1286, doi: 10.1086/118342
  * Quataert et al. (2021a) Quataert, E., Jiang, Y.-F., & Thompson, T. A. 2021a, arXiv e-prints, arXiv:2106.08404.  https://arxiv.org/abs/2106.08404
  * Quataert et al. (2021b) Quataert, E., Thompson, T. A., & Jiang, Y.-F. 2021b, arXiv e-prints, arXiv:2102.05696.  https://arxiv.org/abs/2102.05696
  * Roberg-Clark et al. (2016) Roberg-Clark, G. T., Drake, J. F., Reynolds, C. S., & Swisdak, M. 2016, ApJ, 830, L9, doi: 10.3847/2041-8205/830/1/L9
  * Rubin et al. (2014) Rubin, K. H. R., Prochaska, J. X., Koo, D. C., et al. 2014, ApJ, 794, 156, doi: 10.1088/0004-637X/794/2/156
  * Rubin et al. (2011) Rubin, K. H. R., Prochaska, J. X., Ménard, B., et al. 2011, ApJ, 728, 55, doi: 10.1088/0004-637X/728/1/55
  * Rudie et al. (2019) Rudie, G. C., Steidel, C. C., Pettini, M., et al. 2019, ApJ, 885, 61, doi: 10.3847/1538-4357/ab4255
  * Rupke et al. (2005a) Rupke, D. S., Veilleux, S., & Sanders, D. B. 2005a, ApJS, 160, 87, doi: 10.1086/432886
  * Rupke et al. (2005b) —. 2005b, ApJ, 632, 751, doi: 10.1086/444451
  * Rupke et al. (2019) Rupke, D. S. N., Coil, A., Geach, J. E., et al. 2019, Nature, 574, 643, doi: 10.1038/s41586-019-1686-1
  * Salem & Bryan (2014) Salem, M., & Bryan, G. L. 2014, MNRAS, 437, 3312, doi: 10.1093/mnras/stt2121
  * Scannapieco & Brüggen (2015) Scannapieco, E., & Brüggen, M. 2015, ApJ, 805, 158, doi: 10.1088/0004-637X/805/2/158
  * Schneider et al. (2020) Schneider, E. E., Ostriker, E. C., Robertson, B. E., & Thompson, T. A. 2020, ApJ, 895, 43, doi: 10.3847/1538-4357/ab8ae8
  * Schneider & Robertson (2017) Schneider, E. E., & Robertson, B. E. 2017, ApJ, 834, 144, doi: 10.3847/1538-4357/834/2/144
  * Schneider & Robertson (2018) —. 2018, ApJ, 860, 135, doi: 10.3847/1538-4357/aac329
  * Shapley et al. (2003) Shapley, A. E., Steidel, C. C., Pettini, M., & Adelberger, K. L. 2003, ApJ, 588, 65, doi: 10.1086/373922
  * Sharma & Nath (2013) Sharma, M., & Nath, B. B. 2013, ApJ, 763, 17, doi: 10.1088/0004-637X/763/1/17
  * Silich et al. (2011) Silich, S., Bisnovatyi-Kogan, G., Tenorio-Tagle, G., & Martínez-González, S. 2011, ApJ, 743, 120, doi: 10.1088/0004-637X/743/2/120
  * Silich et al. (2003) Silich, S., Tenorio-Tagle, G., & Muñoz-Tuñón, C. 2003, ApJ, 590, 791, doi: 10.1086/375133
  * Smith et al. (2018) Smith, M. C., Sijacki, D., & Shen, S. 2018, MNRAS, 478, 302, doi: 10.1093/mnras/sty994
  * Socrates et al. (2008) Socrates, A., Davis, S. W., & Ramirez-Ruiz, E. 2008, ApJ, 687, 202, doi: 10.1086/590046
  * Somerville & Davé (2015) Somerville, R. S., & Davé, R. 2015, ARA&A, 53, 51, doi: 10.1146/annurev-astro-082812-140951
  * Sparre et al. (2020) Sparre, M., Pfrommer, C., & Ehlert, K. 2020, MNRAS, 499, 4261, doi: 10.1093/mnras/staa3177
  * Spitzer (1956) Spitzer, L. 1956, Physics of Fully Ionized Gases
  * Springel & Hernquist (2003) Springel, V., & Hernquist, L. 2003, MNRAS, 339, 312, doi: 10.1046/j.1365-8711.2003.06207.x
  * Steidel et al. (2010) Steidel, C. C., Erb, D. K., Shapley, A. E., et al. 2010, ApJ, 717, 289, doi: 10.1088/0004-637X/717/1/289
  * Stern et al. (2020a) Stern, J., Fielding, D., Faucher-Giguère, C.-A., & Quataert, E. 2020a, MNRAS, 492, 6042, doi: 10.1093/mnras/staa198
  * Stern et al. (2020b) Stern, J., Faucher-Giguère, C.-A., Fielding, D., et al. 2020b, arXiv e-prints, arXiv:2006.13976.  https://arxiv.org/abs/2006.13976
  * Strickland & Heckman (2009) Strickland, D. K., & Heckman, T. M. 2009, ApJ, 697, 2030, doi: 10.1088/0004-637X/697/2/2030
  * Suchkov et al. (1994) Suchkov, A. A., Balsara, D. S., Heckman, T. M., & Leitherer, C. 1994, ApJ, 430, 511, doi: 10.1086/174427
  * Suchkov et al. (1996) Suchkov, A. A., Berman, V. G., Heckman, T. M., & Balsara, D. S. 1996, ApJ, 463, 528, doi: 10.1086/177267
  * Sutherland & Dopita (1993) Sutherland, R. S., & Dopita, M. A. 1993, ApJs, 88, 253, doi: 10.1086/191823
  * Tan & Oh (2021) Tan, B., & Oh, S. P. 2021, arXiv e-prints, arXiv:2105.11496.  https://arxiv.org/abs/2105.11496
  * Tan et al. (2021) Tan, B., Oh, S. P., & Gronke, M. 2021, MNRAS, doi: 10.1093/mnras/stab053
  * Tanner et al. (2016) Tanner, R., Cecil, G., & Heitsch, F. 2016, ApJ, 821, 7, doi: 10.3847/0004-637X/821/1/7
  * Thompson et al. (2016) Thompson, T. A., Quataert, E., Zhang, D., & Weinberg, D. H. 2016, MNRAS, 455, 1830, doi: 10.1093/mnras/stv2428
  * Tonnesen & Bryan (2021) Tonnesen, S., & Bryan, G. L. 2021, ApJ, 911, 68, doi: 10.3847/1538-4357/abe7e2
  * van de Voort et al. (2021) van de Voort, F., Bieri, R., Pakmor, R., et al. 2021, MNRAS, 501, 4888, doi: 10.1093/mnras/staa3938
  * van der Velden (2020) van der Velden, E. 2020, The Journal of Open Source Software, 5, 2004, doi: 10.21105/joss.02004
  * Vijayan et al. (2018) Vijayan, A., Sarkar, K. C., Nath, B. B., Sharma, P., & Shchekinov, Y. 2018, MNRAS, 475, 5513, doi: 10.1093/mnras/sty198
  * Virtanen et al. (2020) Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nature Methods, 17, 261, doi: 10.1038/s41592-019-0686-2
  * Vogelsberger et al. (2013) Vogelsberger, M., Genel, S., Sijacki, D., et al. 2013, MNRAS, 436, 3031, doi: 10.1093/mnras/stt1789
  * Voit et al. (2017) Voit, G. M., Meece, G., Li, Y., et al. 2017, ApJ, 845, 80, doi: 10.3847/1538-4357/aa7d04
  * Walter et al. (2002) Walter, F., Weiss, A., & Scoville, N. 2002, ApJ, 580, L21, doi: 10.1086/345287
  * Wang (1995) Wang, B. 1995, ApJ, 444, 590, doi: 10.1086/175633
  * Weaver et al. (1977) Weaver, R., McCray, R., Castor, J., Shapiro, P., & Moore, R. 1977, ApJ, 218, 377, doi: 10.1086/155692
  * Weinberger et al. (2020) Weinberger, R., Springel, V., & Pakmor, R. 2020, ApJS, 248, 32, doi: 10.3847/1538-4365/ab908c
  * Weiner et al. (2009) Weiner, B. J., Coil, A. L., Prochaska, J. X., et al. 2009, ApJ, 692, 187, doi: 10.1088/0004-637X/692/1/187
  * Werk et al. (2014) Werk, J. K., Prochaska, J. X., Tumlinson, J., et al. 2014, ApJ, 792, 8, doi: 10.1088/0004-637X/792/1/8
  * Westmoquette et al. (2011) Westmoquette, M. S., Smith, L. J., & Gallagher, J. S., I. 2011, MNRAS, 414, 3719, doi: 10.1111/j.1365-2966.2011.18675.x
  * Westmoquette et al. (2007) Westmoquette, M. S., Smith, L. J., Gallagher, J. S., I., et al. 2007, ApJ, 671, 358, doi: 10.1086/522693
  * Westmoquette et al. (2009) —. 2009, ApJ, 696, 192, doi: 10.1088/0004-637X/696/1/192
  * White & Frenk (1991) White, S. D. M., & Frenk, C. S. 1991, ApJ, 379, 52, doi: 10.1086/170483
  * White & Rees (1978) White, S. D. M., & Rees, M. J. 1978, MNRAS, 183, 341, doi: 10.1093/mnras/183.3.341
  * Wiener et al. (2019) Wiener, J., Zweibel, E. G., & Ruszkowski, M. 2019, MNRAS, 489, 205, doi: 10.1093/mnras/stz2007
  * Wiersma et al. (2009) Wiersma, R. P. C., Schaye, J., & Smith, B. D. 2009, MNRAS, 393, 99, doi: 10.1111/j.1365-2966.2008.14191.x
  * Williams & McKee (1997) Williams, J. P., & McKee, C. F. 1997, ApJ, 476, 166, doi: 10.1086/303588
  * Xu & Stone (1995) Xu, J., & Stone, J. M. 1995, ApJ, 454, 172, doi: 10.1086/176475
  * Zhang et al. (2017) Zhang, D., Thompson, T. A., Quataert, E., & Murray, N. 2017, MNRAS, 468, 4801, doi: 10.1093/mnras/stx822