# Safe Online Dynamics Learning with Initially Unknown Models and Infeasible
Safety Certificates

Alexandre Capone1, Ryan K. Cosner2, Aaron D. Ames2, and Sandra Hirche1
1Alexandre Capone and Sandra Hirche are with the Department of Electrical and
Computer Engineering, Technical University of Munich, 80333 Munich, Germany.
Emails: $\\{\texttt{alexandre.capone, hirche}\\}$@tum.de).2 Ryan K. Cosner and
Aaron D. Ames are with the Department of Mechanical and Civil Engineering,
California Institute of Technology, Pasadena, CA 91125, USA. Emails:
$\\{\texttt{rkcosner<EMAIL_ADDRESS>
###### Abstract

Control tasks with high levels of uncertainty and safety requirements are
increasingly common. Typically, techniques that guarantee safety during
learning and control utilize constraint-based safety certificates, which can
be leveraged to compute safe control inputs. However, if model uncertainty is
very high, the corresponding certificates are potentially invalid, meaning no
control input satisfies the constraints imposed by the safety certificate.
This paper considers a learning-based setting with a safety certificate based
on a control barrier function second-order cone program. If the control
barrier function certificate is valid, our approach leverages the control
barrier function to guarantee safety. Otherwise, our method explores the
system to recover the feasibility of the control barrier function constraint
as fast as possible. To this end, we employ a method inspired by well-
established tools from Bayesian optimization. We show that if the sampling
frequency is high enough, we recover the feasibility of a control barrier
function second-order cone program, guaranteeing safety. To the best of our
knowledge, this corresponds to the first algorithm that guarantees safety
through online learning without requiring a prior model or backup safe non-
learning-based controller.

##  I INTRODUCTION

With the growing proliferation of robotics in safety-critical fields, e.g.,
autonomous vehicles, medical robotics, and aerospace systems, the need for
methods which ensure the safety of systems and their users has become
paramount. In control theory, guaranteeing safety has become synonymous with
guaranteeing the forward invariance of some user-defined safe set [1]. To
achieve this, several methods have been developed, including Model Predictive
Control (MPC) with state constraints [2], Reachability-based methods [3], and
Control Barrier Functions (CBFs) [1], all of which are capable of providing
rigorous mathematical guarantees of safety.

Unfortunately, the safety guarantees of the previously mentioned methods
typically rely on an assumed perfect knowledge of the system dynamics, which
is not necessarily available in practice. In general, dynamics models have
some associated errors that must be accounted for to achieve safety. This
error is often modeled using data-driven frameworks, such as neural networks
[4, 5] or Gaussian processes (GPs) [6, 7]. In order to still be able to
guarantee safety whenever a residual model is involved, accurate model error
bounds are generally required and need to be carefully considered in the
control design. This has been carried out, e.g., for MPC [8, 9], reachability-
based techniques [10], and control-barrier functions [11, 12]. However, a
considerable drawback of the aforementioned techniques is that they all rely
on the existence of a (potentially conservative) backup safe controller that
they can fall back on in case the system reaches a critical state that the
main algorithm cannot address.

In this work, we present an online learning-based algorithm that is able to
guarantee feasibility without requiring an a-priori safe controller. Our
approach learns the time-derivative of a Control Barrier function online and
attempts to solve an associated quadratic program (QP) that takes the model
error into account. Whenever the QP becomes infeasible, our approach leverages
techniques from Bayesian optimization to ensure that the time-derivative of
the Control Barrier Function is learned fast enough in order to recover
feasibility of the associated SOCP before the system becomes unsafe. Our
method is applicable for a rich class of potential systems and only requires
the knowledge of a CBF for the system. To the best of our knowledge, this
represents the first approach that allows for safe control without any prior
system model or backup safe control law.

The remainder of this paper is structured as follows111Notation:For scalars
$a_{1},...a_{N}$, we use the notation $\text{diag}(a_{1},\ldots,a_{N})$ to
refer to the $N$-dimensional diagonal matrix with $a_{1},...a_{N}$ as diagonal
entries. We employ $|\cdot|$ to denote the determinant operator, and
$\lVert\cdot\rVert_{2}$ to refer to the $2$-norm of both vectors and matrices.
For any set $\mathcal{S}$, we employ $\mathcal{P}(\mathcal{S})$ to refer to
the power set of $\mathcal{S}$.. In Section II we describe the problem setting
considered in this paper, together with some background on control barrier
functions and Gaussian processes. Our main contribution, which includes the
recovery of feasibility of a CBF-SOCP by means of online learning, is
presented in Section III. In Section V, we provide a numerical validation of
our approach on using a cruise control and a quadcopter model. We finalize the
paper with some conclusions, in Section VII.

##  II Background and Problem Setting

Consider the control affine system

$\displaystyle\dot{x}=f(x)+g(x)u$ (1)

where $x\in\mathcal{X}\subseteq\mathbb{R}^{n}$ and
$u\in\mathcal{U}\subseteq\mathbb{R}^{m}$, and
$f:\mathbb{R}^{n}\to\mathbb{R}^{n}$ and
$g:\mathbb{R}^{n}\to\mathbb{R}^{n\times m}$ are (partially) unknown locally
Lipschitz continuous functions that represent the drift dynamics and the input
matrix, respectively. We further assume that $\mathcal{X}$ is an open and
connected set and that $\mathcal{U}$ is compact. When a locally Lipschitz
controller $\pi:\mathbb{R}^{n}\to\mathbb{R}^{m}$ is used, we can define the
closed-loop system:

$\displaystyle\dot{x}=f(x)+g(x)\pi(x)$ (2)

whose solutions we assume are forward complete (i.e., exist for all $t\geq
0$). We also assume to know an upper bound $M_{\dot{x}}\in\mathbb{R}_{+}$ for
the time-derivative of the dynamics, as specified by the following assumption.
i.e.,

######  Assumption 1.

There exists a known positive constant $M_{\dot{x}}$, such that
$\|\dot{x}\|_{2}\leq M_{\dot{x}}$ holds for all $x\in\mathcal{X}$ and
$u\in\mathcal{U}$.

To define safety we consider a safe set $\mathcal{C}$, which we aim to render
forward invariant, as specified in the following.

######  Definition 1 (Forward Invariance (Safety)).

A set $\mathcal{C}\subset\mathcal{X}$ is forward invariant for system (2) if
$x(0)\in\mathcal{C}$ implies that $x(t)\in\mathcal{C}$ for all $t\geq 0$. The
closed loop system (2) is said to be safe with respect to $\mathcal{C}$ if
$\mathcal{C}$ is forward invariant.

Our goal is to design a control law $\pi$ that renders the closed loop system
(2) safe. To this end, we consider the case where the safe set $\mathcal{C}$
corresponds to the superlevel set of some known continuously differentiable
function $h:\mathbb{R}^{n}\to\mathbb{R}$ with 0 a regular value 222A function
$h$ has 0 as a regular value if $h(x)=0\implies\frac{\partial h}{\partial
x}(x)\neq 0$. :

$\displaystyle\mathcal{C}$
$\displaystyle\triangleq\\{x\in\mathbb{R}^{n}~{}|~{}h(x)\geq 0\\}.$

This specification of $\mathcal{C}$ allows us to achieve safety through the
control barrier function framework, which we discuss in the following section.

###  II-A Safety through Control Barrier Functions

Before introducing control barrier functions, we briefly introduce the
concepts of class $\mathcal{K}_{\infty}$ ($\mathcal{K}_{\infty}^{e})$ function
and barrier functions.

We call a continuous function $\alpha:\mathbb{R}\to\mathbb{R}$ an extended
class $\mathcal{K}_{\infty}$ ($\mathcal{K}_{\infty}^{e})$ if it is strictly
monotonically increasing and satisfies $\alpha(0)=0$ and is radially
unbounded. Barrier Functions (BFs) can be used to synthesize controllers
ensuring the safety of the closed-loop system (2) with respect to a given set
$\mathcal{C}$.

######  Definition 2 (Barrier Function (BF) [1]).

Let $\mathcal{C}\subset\mathcal{X}$ be the 0-superlevel set of a continuously
differentiable function $h:\mathbb{R}^{n}\to\mathbb{R}$ with zero a regular
value. The function $h$ is a barrier function (BF) for (2) if there exists
some $\alpha\in\mathcal{K}_{\infty}^{e}$ such that

$\displaystyle\frac{dh}{dt}(x)=\frac{\partial h}{\partial
x}(x)(f(x)+g(x)\pi(x))\geq-\alpha(h(x)),\quad\forall x\in\mathcal{X}$ (3)

Barrier functions have the following associated safety guarantee:

######  Lemma 1 ([1]).

Let $\mathcal{C}\subset\mathbb{R}^{n}$ be the 0-superlevel set of a
continuously differentiable function $h:\mathbb{R}^{n}\to\mathbb{R}$ with
$\frac{\partial h}{\partial x}(x)\neq 0$ when $h(x)=0$. If $h$ is a BF for (2)
on $\mathcal{C}$, then the system (2) is safe with respect to $\mathcal{C}$.

Furthermore, control barrier functions (CBFs) serve a similar role for system
(1) without a pre-defined controller and can be used to generate safe control
inputs.

######  Definition 3 (Control Barrier Function (CBF) [1]).

Let $\mathcal{C}\subset\mathcal{X}$ be the 0-superlevel set of a continuously
differentiable function $h:\mathbb{R}^{n}\to\mathbb{R}$ with 0 a regular
value. The function $h$ is a control barrier function (CBF) for (1) on
$\mathcal{C}$ if there exists an $\alpha\in\mathcal{K}_{\infty}^{e}$ such that
for all $x\in\mathcal{X}$:

$\displaystyle\begin{split}\sup_{u\in\mathcal{U}}\dot{h}(x,u)\triangleq\sup_{u\in\mathcal{U}}\frac{\partial
h}{\partial x}(x)\left(f(x)+g(x)u\right)>-\alpha(h(x)).\end{split}$ (4)

In this paper, we consider Lipschitz continuous functions $h$ which yield
compact safe sets, as stated in the following assumption.

######  Assumption 2.

The safe set $\mathcal{C}$ is compact and the function $h$ admits a Lipschitz
constant $L_{h}$, i.e., $|h(x)-h(x^{\prime})|\leq L_{h}\|x-x^{\prime}\|_{2}$
holds for all $x,x^{\prime}\in\mathcal{X}$.

Given a CBF $h$ and corresponding $\alpha\in\mathcal{K}_{\infty}^{e}$, there
exist a corresponding non-empty set of safe inputs for each $x\in\mathcal{X}$
given as $\Pi_{\textrm{CBF}}:\mathcal{X}\to\mathcal{P}(\mathcal{U})$:

$\displaystyle\Pi_{\textrm{CBF}}(x)\triangleq\left\\{u\in\mathcal{U}~{}\Big{|}~{}\dot{h}(x,u)\geq-\alpha(h(x))\right\\}$
(5)

which can be used to synthesize controllers which guarantee safety.

######  Lemma 2 ([1]).

Let $\mathcal{C}\subset\mathcal{X}$ be a 0-superlevel set of a continuously
differentiable function $h:\mathbb{R}^{n}\to\mathbb{R}$ with 0 a regular
value. If $h$ is a CBF for (1), then the set $\Pi_{\textrm{CBF}}(x)$ is non-
empty for all $x\in\mathcal{X}$ and additionally if
$\pi:\mathcal{X}\to\mathcal{U}$ is a locally Lipschitz controller with
$\pi(x)\in\Pi_{\textrm{CBF}}(x)$ for all $x\in\mathcal{X}$, the closed loop
system (2) is safe with respect to $\mathcal{C}$.

One common method for synthesizing controllers which satisfy
$\pi(x)\in\Pi_{\textrm{CBF}}(x)$ for all $x\in\mathcal{X}$ is by employing the
CBF condition (4) as a constraint in a safety-filter which enforces safety
while achieving minimal deviation from a nominal controller
$\pi_{\textrm{nom}}:\mathcal{X}\to\mathcal{U}$:

$\displaystyle\pi_{\textrm{CBF-
QP}}(x)=\operatorname*{argmin}_{u\in\mathcal{U}}$
$\displaystyle\|u-\pi_{\textrm{nom}}(x)\|_{2}$ (CBF-QP) s.t.
$\displaystyle\frac{\partial h}{\partial
x}(x)\left(f(x)+g(x)u\right)\geq-\alpha(h(x)).$

Given the affine form of this constraint, this safety filter is a quadratic
program (QP) which has a closed form solution [1] and can be solved rapidly
enough for online safe controller synthesis [13] for robotic systems.

Since we seek to consider systems with unknown models, we assume an additional
degree of robustness to the CBF constraint which will allow us to find “safer”
control actions to compensate for this uncertainty:

######  Assumption 3 (Robust CBF Feasibility).

There exists some $\alpha\in\mathcal{K}_{\infty}^{e}$ and a positive scalar
$\epsilon>0$, such that

$\displaystyle\sup_{u\in\mathcal{U}}\frac{\partial h}{\partial
x}(x)\left(f(x)+g(x)u\right)\geq-\alpha(h(x))+\epsilon$ (6)

holds for all $x\in\mathcal{X}$.

We note that given the existence of a CBF, this assumption is not highly
conservative. This can be be seen through the lens of Input-to-State Safety,
as shown in the following example.

######  Example 1.

If the CBF condition (4) holds on $\mathcal{X}$ with $\alpha(r)=ar$ for some
$a>0$, then the Assumption 3 holds for any function $h^{\prime}(x)\triangleq
h(x)+\frac{\epsilon}{a}$ with $\epsilon>0$ such that
$\mathcal{C}_{\epsilon}=\\{x\in\mathbb{R}^{n}~{}|~{}h^{\prime}(x)\geq
0\\}\subset\mathcal{X}$. This can be seen by taking

$\displaystyle\sup_{u\in\mathcal{U}}\dot{h}^{\prime}(x,u)$
$\displaystyle=\sup_{u\in\mathcal{U}}\dot{h}(x,u)$
$\displaystyle=\sup_{u\in\mathcal{U}}\frac{\partial h}{\partial
x}(x)\left(f(x)+g(x)u\right)$ $\displaystyle\geq-ah(x)+\epsilon-\epsilon$
$\displaystyle=-ah^{\prime}(x)+\epsilon$

Thus we recover (6) for $h^{\prime}$. We note that
$\mathcal{C}\subset\mathcal{C}_{\epsilon}\subset\mathcal{X}$, so to attain
this robust CBF feasibility, the set $\mathcal{C}$ should be chosen with some
allowable margin such that, $h(x)\in[-\frac{\epsilon}{a},0]$ does not imply
catastrophic failure.

###  II-B Gaussian Processes and Reproducing Kernel Hilbert Spaces

Consider a prior model of the system $\hat{f}$ and $\hat{g}$, which can be set
to zero whenever no prior knowledge is available. We then model the residual
$f-\hat{f}$ and $g-\hat{g}$ of the unknown functions $f$ and $g$ jointly by
using a Gaussian process (GP) model. In order to then obtain theoretical
guarantees on the growth of the model error as new data is added, we assume
that the residual model belongs to a reproducing kernel Hilbert space (RKHS).
This will be codified later in this section. In the following, we review both
GPs and RKHSs, and provide some preliminary theoretical results.

######  Definition 4 (Gaussian Process (GP) [14]).

A Gaussian Process (GP) is a collection of random variables, any finite number
of which have (consistent) joint Gaussian distributions.

GPs are fully specified by a prior mean, which we set to zero without loss of
generality, and kernel $k:\mathcal{X}\times\mathcal{X}\rightarrow\mathbb{R}$.
Though GPs are typically defined for scalar functions, a GP model for the
vector-valued functions $f$ and $g$ can easily be obtained by defining a
separate GP prior for each dimension. We then employ composite kernels

$\displaystyle\begin{split}k_{i}(\mathbf{z},\mathbf{z}^{\prime})\triangleq
k_{f_{i}}(x,x^{\prime})+\sum_{j=1}^{m}u_{j}k_{g_{(i,j)}}(x,x^{\prime})u_{j}^{\prime},\end{split}$
(7)

where $\mathbf{z}\triangleq(x^{\top},u^{\top})^{\top}$, to model the entries
of $f(x)+g(x)u$ jointly. Here $k_{f_{i}}$ and $k_{g_{(i,j)}}$ are kernels that
capture the behavior of the individual entries of $f$ and $g$, respectively.
In practice, if little is known about the system, then so-called “universal
kernels” (e.g. Matérn or squared-exponential kernel) are employed, which can
approximate continuous functions arbitrarily accurately [15]. In this paper,
we consider the practically relevant squared-exponential kernel

$\displaystyle k_{f_{i}}(\mathbf{z},\mathbf{z}^{\prime})$
$\displaystyle=\sigma_{f_{i}}^{2}\exp\left(-\frac{\|\mathbf{z}-\mathbf{z}^{\prime}\|_{2}^{2}}{2l_{f_{i}}^{2}}\right),$
(8) $\displaystyle k_{g_{(i,j)}}(\mathbf{z},\mathbf{z}^{\prime})$
$\displaystyle=\sigma_{g_{(i,j)}}^{2}\exp\left(-\frac{\|\mathbf{z}-\mathbf{z}^{\prime}\|_{2}^{2}}{2l_{g_{(i,j)}}^{2}}\right),$
(9)

although our approach also is applicable with other kernels. The reason why we
employ a composite kernel (7) to model the entries of $f$ and $g$ jointly is
because this allows us to leverage measurements of $\dot{x}$ to improve our
learned model.

Consider $N$ noisy measurements
$\mathcal{D}_{i,N}\triangleq\\{\mathbf{z}^{(q)},y_{i}^{(q)}\\}_{q=1,...,N}$ of
the i-th entry of the time-derivative of the state, where

$\displaystyle\begin{split}y_{i}^{(q)}=&f_{i}(x^{(q)})-\hat{f}_{i}(x^{(q)})\\\
&+\sum_{j=1}^{m}\left(g_{i,j}(x^{(q)})-\hat{g}_{i,j}(x^{(q)})\right)u_{j}^{(q)}+\xi^{(q)}_{i}\end{split}$
(10)

and the measurement noise $\xi^{(q)}$ satisfies the following assumption.

######  Assumption 4.

The measurement noise $\xi^{(q)}$ is iid zero-mean Gaussian noise with
covariance $\sigma_{i,\text{ns}}^{2}$ for every $i\in\\{1,...,n\\}$.

The posterior of the model at an arbitrary state $\mathbf{z}^{*}$ is then
given by

$\displaystyle\begin{split}&\mu_{i,N}(\mathbf{z}^{*})\triangleq\\\
&\qquad\hat{f}_{i}(x^{*})+\sum_{j=1}^{m}\hat{g}_{i,j}(x^{*})u_{j}^{*}+\mathbf{k}_{i,*}^{\top}\left(\mathbf{K}_{i,N}+\sigma_{i,\text{ns}}^{2}\mathbf{I}_{N}\right)^{-1}\mathbf{y}_{i},\\\
&\sigma^{2}_{i,N}(\mathbf{z}^{*})\triangleq
k_{i}(\mathbf{z}^{*},\mathbf{z}^{*})-\mathbf{k}_{i,*}^{\top}\left(\mathbf{K}_{i,N}+\sigma_{i,\text{ns}}^{2}\mathbf{I}_{N}\right)^{-1}\mathbf{k}_{i,*},\end{split}$

where
$\mathbf{k}_{i,*}\triangleq(k_{i}(\mathbf{z}^{*},\mathbf{z}^{(1)}),...,k_{i}(\mathbf{z}^{*},\mathbf{z}^{(N)}))^{\top}$,
and the entries of the covariance matrix are given by
$[\mathbf{K}_{i,N}]_{pq}=k_{i}(\mathbf{z}^{(p)},\mathbf{z}^{(q)})$. In our
approach, we employ $\mu_{i,N}$ as a model of our system, and $\sigma_{i,N}$
as a measure of uncertainty, which we employ to inform the data-collection
trigger.

In the following, we refer to the mean and covariance matrix of the full
multivariate Gaussian process model as

$\displaystyle\begin{split}\boldsymbol{\mu}_{N}(\cdot)&\triangleq(\mu_{{1,N}}(\cdot),\ldots,\mu_{{n,N}}(\cdot))^{\top},\\\
\boldsymbol{\Sigma}_{N}^{2}(\cdot)&\triangleq\text{diag}(\sigma^{2}_{1,N}(\cdot),\ldots,\sigma^{2}_{n,N}(\cdot)).\end{split}$
(11)

Furthermore, we employ $\mathbf{z}=(x^{\top},u^{\top})^{\top}$ and
$\mathcal{D}_{N}=\\{\mathbf{z}^{(q)},\mathbf{y}^{(q)}\\}_{q=1,...,N}$, where
$\mathbf{y}^{(q)}\triangleq(y_{1}^{(q)},...,y_{n}^{(q)})^{\top}$, to refer to
concatenated state and input, and the full data set corresponding to $N$
system measurements, respectively.

In the following, we assume that the entries
$f_{i}(x)-\hat{f}_{i}(x)+\sum_{j=1}^{m}(g_{i,j}(x)-\hat{g}_{i,j}(x))u_{j}$ of
the residual system dynamics belong to the reproducing kernel Hilbert space
generated by the composite kernel (7).

######  Assumption 5.

For every $i=1,..,n$, the function
$f_{i}(x)-\hat{f}_{i}(x)+\sum_{j=1}^{m}(g_{i,j}(x)-\hat{g}_{i,j}(x))u_{j}$
belongs to the RKHS with reproducing kernel $k_{i}$, and the corresponding
RKHS norm is bounded by a known positive scalar $B_{i}\in(0,\infty)$.

Assumption 5 implies that the entries of the system dynamics can be expressed
as a weighted sum of evaluations of the composite kernel (7). This is not very
restrictive, as the corresponding function space is considerably rich.

We can then leverage Assumption 5 to obtain a bound on the model error.

######  Lemma 3 ([16, Lemma 1]).

Let Assumptions 4 and 5 hold. Then, for any $i=1,...,n$, with probability at
least $1-\delta$,

$\displaystyle\bigg{|}f_{i}(x^{*})+\sum_{j=1}^{m}g(x^{*})u_{j}^{*}-\mu_{i,N}(\mathbf{z}^{*})\bigg{|}\leq\beta_{i,N}\sigma_{i,N}(\mathbf{z}^{*})$
(12)

holds for all $\mathbf{z}\in\mathcal{X}\times\mathcal{U}$ and all
$N\in\mathbb{N}$, where

$\displaystyle\beta_{i,N}\triangleq
B_{i}+\sigma_{i,\text{ns}}\sqrt{2\left(\gamma_{i,N}+1+\ln{(n\delta^{-1})}\right)}$
(13)

and

$\displaystyle\gamma_{i,N}\triangleq\max_{\mathbf{z}^{(1)},...,\mathbf{z}^{(N)}}\frac{1}{2}\left|\mathbf{I}_{N}+\sigma_{i,\text{ns}}^{-1}\mathbf{K}_{i,N}\right|$
(14)

corresponds to the maximal information gain after $N$ rounds of data
collection.

While Lemma 3 allows us to bound the GP model error at any given point
$\mathbf{z}^{*}$ given the data with high probability, to recover feasibility
we need also to understand how $\gamma_{i,N}$ and $\sigma^{2}_{i,N}$ change as
we add more data. To achieve this, we will also employ the following results.

######  Lemma 4 ([17, Lemma 5.4]).

The posterior variance satisfies

$\displaystyle\sum_{q=1}^{N}\sigma^{2}_{i,q}\left(\mathbf{z}^{(q-1)}\right)\leq\frac{2}{\ln{\left(1+\sigma_{i,\text{ns}}^{-2}\right)}}\gamma_{i,N}.$
(15)

Lemma 4 allows us to bound the mean of the error as new data is added.

In order to put everything together, we now only need to determine how quickly
$\gamma_{i,N}$ grows. To this end, we employ the two following results, which
imply that the growth of the error bound is logarithmic at worst.

######  Lemma 5 ([17, Theorem 5]).

Let $k_{f_{i}}$, $k_{g_{i,j}}$, be squared-exponential kernels. Then, there
exist finite positive constants $C_{f_{i}},C_{g_{i}}\in(0,\infty)$, such that,
for $N\in\mathbb{N}$ large enough,

$\displaystyle\max_{\mathbf{z}^{(1)},...\mathbf{z}^{(N)}}\frac{1}{2}\left|\mathbf{I}_{N}+\sigma_{i,\text{ns}}^{-1}\mathbf{K}_{f_{i},N}\right|$
$\displaystyle\leq C_{f_{i}}\left(\log\left(N\right)\right)^{n+m+1}$ (16)
$\displaystyle\max_{\mathbf{z}^{(1)},...\mathbf{z}^{(N)}}\frac{1}{2}\left|\mathbf{I}_{N}+\sigma_{i,\text{ns}}^{-1}\mathbf{K}_{g_{i,j},N}\right|$
$\displaystyle\leq C_{g_{i}}\left(\log\left(N\right)\right)^{n+m+1},$ (17)

where the entries of $\mathbf{K}_{f_{i},N}$ and $\mathbf{K}_{g_{i,j},N}$ are
given by

$\left[\mathbf{K}_{f_{i},N}\right]_{q,l}=k_{f_{i}}(\mathbf{z}_{q},\mathbf{z}_{l}),\qquad\left[\mathbf{K}_{g_{i,j},N}\right]_{q,l}=k_{g_{i,j}}(\mathbf{z}_{q},\mathbf{z}_{l}).$

######  Lemma 6.

Let $k_{i}$ be given as in (7). Then there exists a scalar $C_{i}$, such that

$\displaystyle\gamma_{i,N}\leq C_{i}\left(\log\left(N\right)\right)^{n+m+1}.$
(18)

###### Proof.

The proof follows directly from Lemma 5, together with Theorem 2 and Theorem 3
in [18]. ∎

######  Lemma 7.

Let $\beta_{N}\triangleq\max_{i}\beta_{i,N}$. Then there exists a scalar
$C_{\beta}\in(0,\infty)$, such that, for all $N>1$,

$\displaystyle\beta_{N}\leq
C_{\beta}\left(\log\left(N\right)\right)^{\frac{n+m+1}{2}}.$

###### Proof.

The proof follows directly from the definition of $\beta_{i,N}$ and Lemma 6. ∎

Lemmas 3, 4, 5, 6 and 7 imply that the GP error bound grows logarithmically,
i.e., sublinearly. We will leverage this to show that, if we learn our system
fast enough, then the model error becomes small enough to recover feasibility
before the boundary of the safe set is reached. We will need the following
result to determine the required learning rate, which provides an upper bound
on the point when linear growth overtakes logarithmic growth.

######  Lemma 8.

For $a,b>0$, and any $\lambda$ with

$\displaystyle\lambda\geq\frac{2}{a}\left(\log(b)-\log\left(\frac{a}{2}\right)\right)$

it holds that

$\displaystyle\exp(a\lambda)\geq b(1+\lambda).$

###### Proof.

For any $c>0$, it holds that

$\displaystyle
b\exp(\lambda\exp(-c))=b\sum_{j=0}^{\infty}\exp(-cj)\frac{\lambda^{j}}{j!}$
$\displaystyle=$ $\displaystyle
b\sum_{j=0}^{\infty}\frac{(\lambda\exp(-c))^{j}}{j!}\geq
b\exp(-c)\sum_{j=0}^{\infty}\frac{\lambda^{j}}{j!}$ $\displaystyle\geq$
$\displaystyle b\exp(-c)(1+\lambda).$

Hence, for $c=-\log(a/2)$, we have

$\displaystyle b(1+\lambda)\leq$
$\displaystyle\exp\left(\log(b)+\lambda\frac{a}{2}-\log\left(\frac{a}{2}\right)\right)$
$\displaystyle\leq$
$\displaystyle\exp\left(\lambda\frac{a}{2}+\lambda\frac{a}{2}\right)=\exp(a\lambda).$

∎

##  III Recovering Feasibility with GPs

In this section, we present an algorithm geared towards guaranteeing safety by
exploring the state and input space efficiently, such that the feasibility of
a CBF-based criterion is guaranteed before the boundary of the safe set is
reached, which in turn implies safety. We start by showing that safety is
guaranteed whenever a CBF-based criterion is satisfied, then present our
algorithm, which ensures safety.

For all $N\in\mathbb{N}$, let
$\Pi_{N}:\mathcal{X}\to\mathcal{P}(\mathbb{R}^{m})$ be defined as

$\displaystyle\Pi_{N}(x)\triangleq\left\\{\begin{array}[]{l|c}&\frac{\partial
h}{\partial x}(x)\boldsymbol{\mu}_{N}(x,u)\\\ u\in\mathcal{U}&\geq\\\
&-\alpha(h(x))+L_{h}\beta_{N}\sqrt{\text{tr}\left(\mathbf{\Sigma}_{N}^{2}(x,u)\right)}\end{array}\right\\}$

where $L_{h}>0$ is the global Lipschitz constant of $h$ on $\mathcal{X}$
specified in Assumption 2.

Note that, due to the nature of the composite kernels (7),
$\boldsymbol{\mu}_{N}(x,u)$ is a linear function of $u$ and
$\mathbf{\Sigma}_{N}^{2}(x,u)$ is a diagonal matrix whose entries are positive
definite quadratic functions of $u$. Hence, $\Pi_{N}$ is defined using a
second-order cone (SOC) constraint.

The specification of $\Pi_{N}$ is motivated by the fact that if there exists a
locally Lipschitz $\pi:\mathcal{X}\to\mathcal{U}$ such that
$\pi(x)\in\Pi_{N}(x)$ for all $x\in\mathcal{X}$, then $\mathcal{C}$ can be
rendered forward invariant, i.e., safety is robustly guaranteed for the
uncertain system. We show this formally in the following.

######  Lemma 9.

Let Assumptions 1, 2, 3, 4 and 5 hold, let $\Pi_{N}(x)$ be non-empty for all
$x\in\mathcal{X}$ and an arbitrary $N\in\mathbb{N}$, and let
$\pi:\mathcal{X}\to\mathcal{U}$ be a locally Lipschitz controller such that
$\pi(x)\in\Pi_{N}(x)$ for all $x\in\mathcal{X}$. Then the closed loop system
is safe with respect to $\mathcal{C}$ with probability at least $1-\delta$.

###### Proof.

By Lemma 3, with probability at least $1-\delta$,

$\displaystyle\frac{dh}{dt}(x)=\frac{\partial h}{\partial
x}(x)\left(f(x)+g(x)\pi(x)\right)$ $\displaystyle=$
$\displaystyle\frac{\partial h}{\partial
x}(x)\left(\boldsymbol{\mu}_{N}(x,\pi(x))\right)$
$\displaystyle+\frac{\partial h}{\partial
x}(x)\left(f(x)+g(x)\pi(x)-\boldsymbol{\mu}_{N}(x,\pi(x))\right)$
$\displaystyle\geq$ $\displaystyle\frac{\partial h}{\partial
x}(x)\left(\boldsymbol{\mu}_{N}(x,\pi(x))\right)$
$\displaystyle-\left\|\frac{\partial h}{\partial
x}(x)\right\|_{2}\left\|\left(f(x)+g(x)\pi(x)-\boldsymbol{\mu}_{N}(x,\pi(x))\right)\right\|_{2}$
$\displaystyle\geq$ $\displaystyle\frac{\partial h}{\partial
x}(x)\left(\boldsymbol{\mu}_{N}(x,\pi(x))\right)-\left\|\frac{\partial
h}{\partial
x}(x)\right\|_{2}\beta_{N}\sqrt{\text{tr}\left(\mathbf{\Sigma}_{N}^{2}(x,u)\right)}$
$\displaystyle\geq$ $\displaystyle-\alpha(h(x)),$

holds for all $x\in\mathcal{X}$, where the last inequality comes from the fact
that $\pi(x)\in\Pi_{N}(x)$ for all $x$. Since $\pi$ satisfies the the
properties of Lemma 2 with probability $1-\delta$, the set $\mathcal{C}$ is
safe with probability $1-\delta$. ∎

(a) System state exits the region where (GP-CBF-SOCP) is infeasible.

(b) Algorithm starts navigating and collecting data by applying (19).

(c) Theorem 1 guarantees that the feasibility of (GP-CBF-SOCP) is recovered
before the boundary of the safe set is reached, which enables us to compute
safe inputs.

Figure 1: Illustration of Algorithm 1 and Theorem 1.

Given a nominal controller $\pi_{\textrm{nom}}:\mathcal{X}\to\mathcal{U}$, we
now formulate the robust CBF-SOCP that takes the GP model error into account:

$\displaystyle\begin{matrix}\pi_{\textrm{safe}}(x)=\operatorname*{argmin}_{u\in\mathcal{U}}\|u-\pi_{\textrm{nom}}(x)\|_{2}\\\
\qquad\hskip 20.0pt\textrm{s.t. }u\in\Pi_{N}(x).\end{matrix}$

We now provide data sampling schemes that recover feasibility of the GP-CBF-
SOCP before reaching the boundary of the safe set $\mathcal{C}$. The concept
and theoretical guarantees of our approach borrow from the well-known GP upper
confidence bound algorithm [17], which aims to find the maximum of an unknown
function by balancing exploration and exploitation. Whenever the SOC
constraint is infeasible (i.e., $\Pi_{N}(x)$ is empty) given the present data
set $\mathcal{D}_{N}$, we propose choosing the input $u$ by solving the
optimization problem

$\displaystyle\begin{split}u^{(N+1)}=&\arg\sup_{u\in\mathcal{U}}\Bigg{[}\frac{\partial
h}{\partial x}\left(x^{(N+1)}\right)\mu_{N}\left(x^{(N+1)},u\right)\\\
&+L_{h}\beta_{N}\sqrt{\text{tr}\left(\mathbf{\Sigma}_{N}^{2}\left(x^{(N+1)},u\right)\right)}\Bigg{]}.\end{split}$
(19)

The idea behind solving (19) is to try and find the control input $u$ that
maximizes the unknown time derivative $\frac{dh}{dt}(x)$ of $h(x)$.

######  Remark 1.

Note that we are adding the model uncertainty to the model in (19), as opposed
to subtracting it as we do in the formulation of the feasible set of the SOCP
$\Pi_{N}(x)$. This is because we want to emphasize exploration over safety
whenever the (GP-CBF-SOCP) is infeasible.

It can be shown that, under the assumption that the state and input space can
be sampled arbitrarily, i.e., without constraints imposed by the dynamical
system, solving (19) recursively eventually yields a result that is
arbitrarily close to the true maximum $\sup_{u\in\mathcal{U}}\frac{\partial
h}{\partial x}(x)(f(x)+g(x)u)$, eventually recovering feasibility of the (GP-
CBF-SOCP) [17]. However, in our setting, we need to recover feasibility before
the system reaches the boundary of the safe set $\mathcal{C}$, i.e., we need
to sample sufficiently frequently so as to recover feasibility before the
boundary of the safe set is reached. To this end, we present an algorithm that
samples new data points at pre-specified time intervals in a way that
guarantees that $\Delta N$ data points are sampled before the boundary of the
safe set $\mathcal{C}$ is reached. This is achieved by leveraging the upper
bound on the time-derivative of the system to compute the maximal sampling
frequency required to collect at least $\Delta N$ data points before reaching
the boundary. This is presented in Algorithm 1. Note that we also require the
control policy $\bar{\pi}(x^{(q)})$ that is applied to the system directly
after sampling to be Lipschitz continuous and to correspond to $u^{(q)}$ at
$x^{(q)}$. This corresponds to Line 10 in Algorithm 1.

We can then show that, if we choose the number of data points $\Delta N$ to be
collected before reaching the boundary $\partial\mathcal{C}$ high enough,
then, with high probability, Algorithm 1 always recovers feasibility before
reaching the boundary $\partial\mathcal{C}$, guaranteeing safety at all times,
and stop collecting data after at most $\Delta N$ points have been collected,
meaning that the Algorithm 1 only requires a finite amount of memory.

Algorithm 1 CBF-Control and Sampling Scheme with Guaranteed Feasibility
Recovery

0: Size of data set to be collected $\Delta N$

1: if $\Pi_{N}(x)$ is not empty then

2: Solve (GP-CBF-SOCP) and apply $\pi_{\text{safe}}(x)$

3: else

4: Set $\Delta t={h(x)}\left({M_{\|\frac{\partial h}{\partial
x}\|}M_{\dot{x}}\Delta N}\right)^{-1}$

5: Set $q=0$, $t_{\text{infeasible}}=t$.

6: while $\Pi_{N+q}(x)$ is empty do

7: Set $q=q+1$.

8: Set $x^{(q)}=x$ and compute $u^{(q)}$ by solving (19).

9: while $t\leq t_{\text{infeasible}}+q\Delta t$ do

10: Apply any admissible locally Lipschitz controller $\bar{\pi}$ with
$\bar{\pi}(x^{(q)})=u^{(q)}$ to the system.

11: end while

12: Set $\mathbf{y}^{(q)}=\dot{x}^{(q)}+\mathbf{\xi}^{(q)}$.

13: Set $\mathcal{D}_{\Delta N}=\mathcal{D}_{\Delta
N}\bigcap\\{\mathbf{z}^{(q)},\mathbf{y}^{(q)}\\}$.

14: Update GP and $\Pi_{N+q}$.

15: end while

16: Set $N=N+q$ and $\Delta N=\Delta N-q$.

17: end if

######  Theorem 1.

Let Assumptions 1, 2, 3, 4 and 5 hold, and let $x(0)$ be within the interior
of the safe set $\mathcal{C}$, i.e.,
$x(0)\in\mathcal{C}\backslash\partial\mathcal{C}$. Furthermore, let $N$ denote
the initial amount of data points used to train the GP, and choose

$\eta\triangleq\left(\sum_{i=1}^{n}\frac{4L_{h}C_{\beta}C_{i}}{\epsilon\ln{\left(1+\sigma_{i,\text{ns}}^{-2}\right)}}\right)^{-\frac{2}{3(n+m+1)}}\frac{1}{3(n+m+1)}$

where $C_{i}$ and $C_{\beta}$ are chosen as in Lemma 6 and Lemma 7,
respectively. Then, if we employ Algorithm 1 with

$\displaystyle\Delta N\geq$
$\displaystyle\eta^{-1}\left(\log\left(N\right)-\log\left(\eta\right)\right),$

with probability at least $1-\delta$, feasibility is recovered using a finite
number of points before the state $x$ reaches $\partial\mathcal{C}$ and only
samples the state space at most $\Delta N$ times, thus for any piecewise
locally Lipschitz controller $\pi(x)\in\Pi_{N}(x)$ the system (2) safe with
respect to $\mathcal{C}$ with probability at least $1-\delta$.

###### Proof.

Note that Algorithm 1 collects at least $\Delta N$ data points before leaving
the safe set $\mathcal{C}$. Hence, it is sufficient to show that $\Pi_{N}(x)$
is non-empty for all $x\in\mathcal{X}$ after $\Delta N$ points have been
collected. Safety with respect to $\mathcal{C}$ with probability at least
$1-\delta$ then follows from Lemma 9.

We then show that $\Pi_{N}(x)$ is non-empty for all $x\in\mathcal{X}$ after at
most $\Delta N$ points have been collected by contradiction. Let
$\tilde{N}\triangleq N+q$. Note that Algorithm 1 only collects data points
whenever $\Pi_{\tilde{N}}(x)$ is empty, implying

$\displaystyle\frac{\partial h}{\partial
x}(x^{(\tilde{N}+1)})\mu_{\tilde{N}}(x^{(\tilde{N}+1)},u)$ $\displaystyle<$
$\displaystyle
L_{h}\beta_{\tilde{N}}\sqrt{\text{tr}\left(\mathbf{\Sigma}_{\tilde{N}}^{2}(x^{(\tilde{N}+1)},u)\right)}-\alpha(h(x^{(\tilde{N}+1)}))$

holds for all $u\in\mathcal{U}$. Now, let

$\displaystyle u^{**}_{\tilde{N}+1}\triangleq$
$\displaystyle\arg\sup_{u\in\mathcal{U}}\frac{\partial h}{\partial
x}\left(x^{(\tilde{N}+1)}\right)\left(f\left(x^{(\tilde{N}+1)}\right)+g\left(x^{(\tilde{N}+1)}\right)u\right),$

and note that, with probability at least $1-\delta$,

$\displaystyle\begin{split}&\frac{\partial h}{\partial
x}\left(x^{(\tilde{N}+1)}\right)\left(f\left(x^{(\tilde{N}+1)}\right)+g\left(x^{(\tilde{N}+1)}\right)u^{**}_{\tilde{N}+1}\right)\\\
\leq&\frac{\partial h}{\partial
x}\left(x^{(\tilde{N}+1)}\right)\mu_{\tilde{N}}\left(x^{(\tilde{N}+1)},u^{**}_{\tilde{N}+1}\right)\\\
&+L_{h}\beta_{\tilde{N}}\sqrt{\text{tr}\left(\mathbf{\Sigma}_{\tilde{N}}^{2}\left(x^{(\tilde{N}+1)},u^{**}_{\tilde{N}+1}\right)\right)}\\\
\leq&\frac{\partial h}{\partial
x}\left(x^{(\tilde{N}+1)}\right)\mu_{\tilde{N}}\left(x^{(\tilde{N}+1)},u^{(\tilde{N}+1)}\right)\\\
&+L_{h}\beta_{\tilde{N}}\sqrt{\text{tr}\left(\mathbf{\Sigma}_{\tilde{N}}^{2}\left(x^{(\tilde{N}+1)},u^{(\tilde{N}+1)}\right)\right)}.\end{split}$
(20)

We then obtain

$\displaystyle-\alpha(h(x^{(\tilde{N}+1)}))$ $\displaystyle>$
$\displaystyle\frac{\partial h}{\partial
x}(x^{(\tilde{N}+1)})\mu_{\tilde{N}}\left(x^{(\tilde{N}+1)},u^{(\tilde{N}+1)}\right)$
$\displaystyle-
L_{h}\beta_{\tilde{N}}\sqrt{\text{tr}\left(\mathbf{\Sigma}_{\tilde{N}}^{2}(x^{(\tilde{N}+1)},u^{(\tilde{N}+1)})\right)}$
$\displaystyle\geq$ $\displaystyle\frac{\partial h}{\partial
x}\left(x^{(\tilde{N}+1)}\right)\left(f\left(x^{(\tilde{N}+1)}\right)+g\left(x^{(\tilde{N}+1)}\right)u^{**}_{\tilde{N}+1}\right)$
$\displaystyle-2L_{h}\beta_{\tilde{N}}\sqrt{\text{tr}\left(\mathbf{\Sigma}_{\tilde{N}}^{2}(x^{(\tilde{N}+1)},u^{(\tilde{N}+1)})\right)}$
$\displaystyle\geq$
$\displaystyle-\alpha\left(h(x^{(\tilde{N}+1)})\right)+\epsilon$
$\displaystyle-2L_{h}\beta_{\tilde{N}}\sqrt{\text{tr}\left(\mathbf{\Sigma}_{\tilde{N}}^{2}(x^{(\tilde{N}+1)},u^{(\tilde{N}+1)})\right)}$

i.e.,

$\displaystyle
2L_{h}\beta_{\tilde{N}}\sqrt{\text{tr}\left(\mathbf{\Sigma}_{\tilde{N}}^{2}\left(x^{(\tilde{N}+1)},u^{(\tilde{N}+1})\right)\right)}\geq\epsilon.$

For simplicity of exposition, we introduce

$\psi\triangleq\left(\sum_{i=1}^{n}\frac{4L_{h}}{\ln{\left(1+\sigma_{i,\text{ns}}^{-2}\right)}}C_{\beta}C_{i}\right).$

By summing up the posterior covariance terms over the collected data points
and employing Lemmas 6 and 7, we then obtain

$\displaystyle-\Delta N\epsilon\geq$ $\displaystyle\sum_{q=1}^{\Delta
N}-2L_{h}\beta_{N+q}\sqrt{\text{tr}\left(\mathbf{\Sigma}_{N+q}^{2}(x^{({N}+q)},u^{(N+q)}_{{N}+q})\right)}$
$\displaystyle\geq$ $\displaystyle-\sum_{q=1}^{\Delta
N}\left(2L_{h}\beta_{N+q}\sum_{i=1}^{n}\frac{2}{\ln{\left(1+\sigma_{i,\text{ns}}^{-2}\right)}}\gamma_{i,N+q}\right)$
$\displaystyle\geq$ $\displaystyle-\psi\log(N+\Delta N)^{\frac{3}{2}(n+m+1)}.$

Hence,

$\displaystyle\exp\left(\left(\frac{\Delta
N\epsilon}{\psi}\right)^{\frac{2}{3(n+m+1)}}\right)\leq N+\Delta N,$

holds for all $\Delta N>0$. Without loss of generality, we can substitute
$\Delta N$ with $(\Delta N)^{\frac{3(n+m+1)}{2}}$, which yields

$\displaystyle\exp\left(\Delta
N\left(\frac{\epsilon}{\psi}\right)^{\frac{2}{3(n+m+1)}}\right)\leq N+(\Delta
N)^{\frac{3(n+m+1)}{2}}$ $\displaystyle\leq$ $\displaystyle(N+\Delta
N)^{\frac{3(n+m+1)}{2}},$

i.e.,

$\displaystyle\exp\left(2\Delta N\eta\right)\leq N+\Delta N.$

Since we have

$\displaystyle\Delta N\geq$
$\displaystyle\eta^{-1}\left(\left(\log\left(N\right)-\log\left(\eta\right)\right)\right),$

this is a contradiction by Lemma 8. By the same argument, Algorithm 1 collects
at most $\Delta N$ points in the state space where $\Pi_{N+q}(x)$ is empty.

∎

######  Remark 2.

Although Algorithm 1 employs a temporal trigger to stipulate when exploration
takes place, a trigger based on the distance to the boundary or the value of
the CBF can also be employed, e.g., by dividing the CBF into $\Delta N$
segments and sampling whenever a new segment is reached.

##  IV Discussion

###  IV-A Choice of $\Delta N$.

Though the constants $\mathcal{C}_{i}$, $C_{\beta}$ and $\beta_{i}$ required
by Theorem 1 can be computed explicitly [17, 19], they can result in very
conservative values for $\Delta N$. However, this is only a necessary
criterion, as opposed to a sufficient one, and practice we can choose lower
values of $\Delta N$. This is because $\Delta N$ corresponds to the maximal
amount of data that is required in order for $\Pi_{\tilde{N}}(x)$ to be non-
empty for all $x\in\mathcal{C}$, and often only a subset of $\mathcal{C}$ is
visited during control, i.e., it is sufficient for $\Pi_{\tilde{N}}(x)$ to be
non-empty only for a subset of the safe set $\mathcal{C}$.

###  IV-B Practicability of Online Exploration.

A practical concern that may generally arise with control inputs $u$ that are
geared towards exploring the state and input spaces is that the corresponding
inputs may place a significant strain on the system and lead to undesirable
behavior, e.g., if the input exhibits a high frequency and amplitude. However,
computing the input by solving (19) corresponds to choosing the optimistically
safest input under uncertainty. Hence it is reasonable to expect the
corresponding input to be acceptable for the system, i.e., not damaging.

Figure 2: Value of control barrier function $h(\mathbf{x})$ for cruise control
example. Since no prior model is available for control, safety is obtained by
efficiently learning a model online. Figure 3: Estimated (worst-case) time
derivative of control barrier function for cruise control example. Spikes are
due to training data set updates, leading to decreased model uncertainty.
Positive values indicate infeasibility of the (GP-CBF-SOCP).

##  V Numerical Validation

We now showcase how our approach performs using two numerical simulations. We
start with a cruise control system, then present results for a quadrotor with
ground dynamics, which is more complex. Note that, in the following, we assume
to have either no prior model (cruise control example) or to know only the
model component corresponding to the time-derivatives (quadrotor), which is
insufficient to implement any state-of-the-art approach. We finish this
section by comparing our approach to persistence of excitation.

###  V-A Cruise Control

Our approach is employed to learn the road vehicle model presented in [11]
while simultaneously applying an adaptive cruise control system. In the
following, we omit phyisical dimensions when describing the system model. The
state space model is as in (1), with unknown state-dependent functions

$\displaystyle
f(x)=\begin{bmatrix}-\frac{1}{m}(\zeta_{0}+\zeta_{1}v+\zeta_{2}v^{2})\\\
v_{0}-v\end{bmatrix},\quad g(x)=\begin{bmatrix}0\\\ \frac{1}{m}\end{bmatrix}$
(21)

and state $x=[v\ z]^{\top}$. Here $z$ denotes the distance between the ego
vehicle and the target vehicle in front, $v$ denotes the ego vehicle speed,
$m=1650$ its mass, and $\zeta_{0}=0.2$, $\zeta_{1}=10$, $\zeta_{2}=0.5$ are
parameters that specify the rolling resistance. As a control barrier function,
we employ $h(\mathbf{x})=z-T_{h}v$, where $T_{h}=1.8$, which aims to maintain
a safe distance between the ego vehicle and the vehicle in front. The nominal
controller $\pi_{\textrm{nom}}(x)$ used for the (GP-CBF-SOCP) is a
P-controller $\pi_{\textrm{nom}}=-10(v-v_{d})$, where $v_{d}=24$ corresponds
to the desired velocity. We assume to have $N=10$ data points at the start of
the simulation, which we employ exclusively to learn the lengthscales and
signal variances of the covariance kernels by minimizing the posterior
likelihood [14]. Though the amount of data $\Delta N$ stipulated by Theorem 1
yields strict theoretical guarantees, this represents only a sufficient, and
not a necessary condition, and can be conservative in practice. Hence, to
additionally showcase the practical applicability of our approach, we set
$\Delta t=10^{-5}$, which corresponds to a high rate for many practical
applications.

We simulate the system for $100$ seconds. The CBF value can be seen in Fig. 2,
the worst-case estimated value of the CBF constraint is depicted in Fig. 3.
The CBF value $h(x)$ is always above zero, meaning that safety is always kept
at all times. This is to be expected from Theorem 1. The (CBF-SOCP) is
infeasible during many simulation instances, particularly in the beginning,
which leads to a high rate of exploratory inputs, obtained by solving (19).
However, feasibility is recovered after approximately $6$ seconds, after which
a safe input can be obtained without further exploration. Note that the
boundary of the safe set is not reached. This is because the proposed
algorithm stops exploring when a safe controller can be computed, i.e., the
(GP-CBF-SOCP) becomes feasible. Although closer proximity to the boundary of
the safe set can be potentially achieved through exploration, we leave this to
future work.

Figure 4: Value of control barrier functions for quadrotor example.. Figure 5:
Estimated (worst-case) time derivative of control barrier function for
quadrotor example.

(a) Cruise control.

(b) Quadrotor.

Figure 6: Rate of failure using our approach and random inputs for the cruise
control (left) and quadrotor (right) settings.

###  V-B Quadrotor

We now showcase our approach using a numerical simulation of a quadrotor with
ground dynamics. The quadrotor dynamics are specified by the functions

$\displaystyle\dot{p}$ $\displaystyle=v,\quad$ $\displaystyle\dot{v}$
$\displaystyle=g_{\textup{gr}}e_{z}+\zeta(p_{z})Re_{z}T,$ (22a)
$\displaystyle\dot{R}$ $\displaystyle=R[\omega]_{\times}$ (22b)

where $p\in\mathbb{R}^{3}$ is the global position, $v\in\mathbb{R}^{3}$ the
global velocity, and $R\in\text{SO}(3)$ the system orientation. The parameter
$g_{\textup{gr}}=9.81$ denotes gravity,
$[\;\cdot\;]_{\times}:\mathbb{R}^{3}\rightarrow\mathbb{R}^{3\times 3}$ is the
skew-symmetric mapping, and $e_{z}=[0\ 0\ 1]^{\top}$ is the unit-$z$ vector.
The function $\zeta:\mathbb{R}_{+}\rightarrow[0,1]$ models ground effects and
takes the quadrotor height $p_{z}$ as an input variable. It is computed as
[20]

$\displaystyle\zeta(p_{z})=1-\rho\left(\frac{r_{\text{rot}}}{4p_{z}}\right)^{2},$
(23)

where $\rho=5$ and $r_{\text{rot}}=0.09$ is the rotor radious. The control
inputs are $T$, the thrust acceleration, and $\omega\in SO(3)$, the body-frame
angle rates:

$u=[T\ \omega]^{\top}\in\mathbb{R}\times SO(3).$

In this example, we employ two CBFs. The first is
$h(\mathbf{x})=10(p_{z}-T_{z}v_{z})$, where $T_{z}=0.1$, and is geared toward
keeping altitude higher than zero. The second CBF utilizes the differential
flatness of the quadrotor to synthesize a viable safe set. In particular we
are concerned with restricting the system to safe positions defined as the
0-superlevel set of $h_{p}(p)=r^{2}-\|p\|^{2}$. We extend $h_{p}$ to include
velocities as $h_{e}(p,v)=\dot{h}_{p}(p,v)+\alpha h_{p}(p)$ for some
$\alpha>0$ to produce a relative degree 1 CBF for a double integrator system
as in [21]. To include orientation we add a rotation term to produce a CBF for
the drone, $h(p,v,R)=h_{e}(p,v)-\lambda(1-\frac{1}{2r}\frac{\partial
h}{\partial p}Re_{z})$, with some $\lambda\in(0,r^{2}/2)$ which shrinks the
safe set to ensure that the thrust vector is pointing inwards whenever
$h_{e}(p,v)=0$. When computing the exploring input (19), we alternate between
CBFs. The nominal controller $\pi_{\textrm{nom}}(x)$ used for the (GP-CBF-
SOCP) corresponds to a differentially flat controller, computed as in [22],
and we consider bounded thrust, with $|T|\leq 15000$. Note that $\omega\in
SO(3)$ is already bounded. Similarly to the cruise control setting, we assume
to have $N=10$ data points at the start of the simulation to learn the kernel
hyperparameters and set $\Delta t=10^{-5}$.

We simulate the system for $50$ seconds. The CBF value can be seen in Fig. 4,
the worst-case estimated value of the CBF constraint is depicted in Fig. 5.
Similarly to the cruise control case, instances, when the (GP-CBF-SOCP) is
infeasible, are initially frequent, leading to a high data collection rate.
After two seconds, feasibility is recovered, and safety is guaranteed by
solving the (GP-CBF-SOCP).

##  VI Comparison with Persistence of Excitation

Theorem 1 states that it is sufficient to collect data by applying (19) to the
system in order to guarantee safety. However, other types of control inputs
may also satisfy this requirement. In the following, we investigate how our
approach performs compared to persistence of excitation. More specifically, we
apply our method to the system with the following difference: instead of
computing exploratory inputs by solving (19), we sample the control inputs
from a uniform distribution on $\mathcal{U}$, which corresponds to a
persistently exciting signal [23]. We again consider the adaptive cruise
control and quadrotor settings and investigate how our approach and the random
input-based one perform if the maximal sampling frequency $\frac{1}{\Delta t}$
is bounded. This is relevant, as many practical settings do not allow for
arbitrarily high sampling frequencies.

We perform $100$ simulations with different initial conditions, uniformly
sampled from a region within the safe set. We report how often each method
fails, i.e., leads to a positive value for the CBF during the simulation. The
average number of failures for both settings is shown in Figure 6. A non-zero
failure rate is expected at low sampling frequencies since we cannot
efficiently learn the system model if too little data is collected. However,
as can be seen, our approach nonetheless performs better than the random
control input-based approach. This is because the inputs applied to the system
are geared towards recovering the feasibility of the (GP-CBF-SOCP), whereas
random inputs are not. This also results in a higher data efficiency, as
reflected in the average collected data, which is lower for our approach. At
higher sampling rates, data efficiency is higher, leading to less collected
data. This is because a higher sampling rate means that the elapsed time
between data collection and model update is smaller, i.e., the model captures
the true system more faithfully immediately after an update. This suggests
that, if the maximal sampling frequency is low, then particular effort should
be put into updating the GP as fast as possible.

##  VII Conclusion

We have presented an online learning-based approach to recovering feasibility
of a CBF-QP before reaching the boundary of a safe set, thus guaranteeing
safety with high probability. This is achieved by leveraging tools commonly
used in Bayesian optimization to devise an exploration approach that makes
progress towards learning safe inputs. In future, we aim to apply the proposed
approach to real-life and more complex systems.

## References

  * [1] A. D. Ames, X. Xu, J. W. Grizzle, and P. Tabuada, “Control Barrier Function Based Quadratic Programs for Safety Critical Systems,” _IEEE Transactions on Automatic Control_ , vol. 62, no. 8, pp. 3861–3876, Aug. 2017, arXiv: 1609.06408. [Online]. Available: http://arxiv.org/abs/1609.06408
  * [2] K. P. Wabersich and M. N. Zeilinger, “Linear Model Predictive Safety Certification for Learning-Based Control,” in _2018 IEEE Conference on Decision and Control (CDC)_ , Dec. 2018, pp. 7130–7135, iSSN: 2576-2370.
  * [3] S. Bansal, M. Chen, S. Herbert, and C. J. Tomlin, “Hamilton-jacobi reachability: A brief overview and recent advances,” in _2017 IEEE 56th Annual Conference on Decision and Control (CDC)_. IEEE, 2017, pp. 2242–2253.
  * [4] A. J. Taylor, A. Singletary, Y. Yue, and A. D. Ames, “Learning for Safety-Critical Control with Control Barrier Functions,” p. 10.
  * [5] N. Csomay-Shanklin, R. K. Cosner, M. Dai, A. J. Taylor, and A. D. Ames, “Episodic Learning for Safe Bipedal Locomotion with Control Barrier Functions and Projection-to-State Safety,” p. 13.
  * [6] A. Lederer, A. Begzadić, N. Das, and S. Hirche, “Safe learning-based control of elastic joint robots via control barrier functions,” _arXiv preprint arXiv:2212.00478_ , 2022.
  * [7] I. D. J. Rodriguez, U. Rosolia, A. D. Ames, and Y. Yue, “Learning unstable dynamics with one minute of data: A differentiation-based gaussian process approach,” _arXiv preprint_ , 2021.
  * [8] L. Hewing, K. P. Wabersich, M. Menner, and M. N. Zeilinger, “Learning-based model predictive control: Toward safe learning in control,” _Annual Review of Control, Robotics, and Autonomous Systems_ , vol. 3, pp. 269–296, 2020\.
  * [9] T. Koller, F. Berkenkamp, M. Turchetta, and A. Krause, “Learning-based model predictive control for safe exploration,” in _IEEE Conference on Decision and Control_ , 2018, pp. 6059–6066.
  * [10] A. K. Akametalu, J. F. Fisac, J. H. Gillula, S. Kaynama, M. N. Zeilinger, and C. J. Tomlin, “Reachability-based safe learning with gaussian processes,” in _53rd IEEE Conference on Decision and Control_ , 2014, pp. 1424–1431.
  * [11] F. Castañeda, J. J. Choi, B. Zhang, C. J. Tomlin, and K. Sreenath, “Pointwise feasibility of gaussian process-based safety-critical control under model uncertainty,” in _2021 60th IEEE Conference on Decision and Control (CDC)_ , 2021, pp. 6762–6769.
  * [12] P. Jagtap, G. J. Pappas, and M. Zamani, “Control barrier functions for unknown nonlinear systems using gaussian processes,” in _2020 59th IEEE Conference on Decision and Control (CDC)_. IEEE, 2020, pp. 3699–3704.
  * [13] T. Gurriet, A. Singletary, J. Reher, L. Ciarletta, E. Feron, and A. Ames, “Towards a framework for realizable safety critical control through active set invariance,” in _2018 ACM/IEEE 9th International Conference on Cyber-Physical Systems (ICCPS)_. IEEE, 2018, pp. 98–106.
  * [14] C. E. Rasmussen and C. K. Williams, “Gaussian processes for machine learning. 2006,” _The MIT Press, Cambridge, MA, USA_ , 2006.
  * [15] C. A. Micchelli, Y. Xu, and H. Zhang, “Universal kernels,” _Journal of Machine Learning Research_ , vol. 7, no. Dec, pp. 2651–2667, 2006.
  * [16] A. Capone and S. Hirche, “Backstepping for partially unknown nonlinear systems using Gaussian processes,” _IEEE Control Systems Letters_ , vol. 3, pp. 416–421, 2019.
  * [17] N. Srinivas, A. Krause, S. M. Kakade, and M. W. Seeger, “Information-theoretic regret bounds for Gaussian process optimization in the bandit setting,” _IEEE Transactions on Information Theory_ , vol. 58, no. 5, pp. 3250–3265, 2012.
  * [18] A. Krause and C. Ong, “Contextual gaussian process bandit optimization,” in _Advances in Neural Information Processing Systems_ , J. Shawe-Taylor, R. Zemel, P. Bartlett, F. Pereira, and K. Weinberger, Eds., vol. 24. Curran Associates, Inc., 2011. [Online]. Available: https://proceedings.neurips.cc/paper/2011/file/f3f1b7fc5a8779a9e618e1f23a7b7860-Paper.pdf
  * [19] S. R. Chowdhury and A. Gopalan, “On kernelized multi-armed bandits,” _arXiv preprint arXiv:1704.00445_ , 2017.
  * [20] L. Danjun, Z. Yan, S. Zongying, and L. Geng, “Autonomous landing of quadrotor based on ground effect modelling,” in _2015 34th Chinese Control Conference (CCC)_ , 2015, pp. 5647–5652.
  * [21] Q. Nguyen and K. Sreenath, “Exponential control barrier functions for enforcing high relative-degree safety-critical constraints,” in _2016 American Control Conference (ACC)_. IEEE, 2016, pp. 322–328.
  * [22] M. Faessler, A. Franchi, and D. Scaramuzza, “Differential flatness of quadrotor dynamics subject to rotor drag for accurate tracking of high-speed trajectories,” _IEEE Robotics and Automation Letters_ , vol. 3, no. 2, pp. 620–626, 2018.
  * [23] B. D. Anderson, “Adaptive systems, lack of persistency of excitation and bursting phenomena,” _Automatica_ , vol. 21, no. 3, pp. 247–258, 1985. [Online]. Available: https://www.sciencedirect.com/science/article/pii/0005109885900585