# Fluids in art: “The water’s language was a wondrous one, some narrative on a
recurrent subject …”
Unexpected discoveries are made when classical paintings are analysed
on the subject of the faithful portrayal of ubiquitous fluid phenomena.

R. Krechetnikov  University of Alberta, Edmonton, Alberta, Canada T6G 2E1

###### Abstract

Artists spent a great deal of time studying anatomy for precise rendering of
the human body as well as light, shadows, and perspective for convincing
representation of the three-dimensional world. But in many paintings, they
also had to depict fluids in their static and dynamic states – a subject they
could not study thoroughly, which led to a number of glaring
misrepresentations or deliberate deceits.

##  I “Painters make their defeat certain by attempting to draw running water,
which is a lustrous object in rapid motion.” (John Ruskin 1)

Since ancient times, artists studied proportions – the Greek sculptor
Polykleitos (c. 450-420 B.C.), known for his ideal bronze Doryphoros and an
influential Canon describing the proportions to be followed in sculpture – as
part of human anatomy for realistic and convincing representations of their
subject matter. In the treatise Della pittura (1435; “On Painting”) Leon
Battista Alberti urged painters to construct the human figure as it exists in
nature, with musculature built on the skeleton, and only then draped in skin.
For the same reason, it is believed that attempts to develop a system of
perspective began around the fifth century B.C. in Greece, to achieve
illusionism in theatrical scenery. While first intuitive, these studies were
also driven by a systematic look at the underpinning science – the best
examples including Leonardo da Vinci who performed dissections to uncover the
structure of the human body resulting in the famous anatomical drawings, which
are among the most significant achievements of the Renaissance. Altogether, on
the artist side, there has been an effort to convince the connoisseur of art.
Can we say the same about rendering omnipresent fluid phenomena?

Similarly, artists had to represent still and flowing water since days of
yore, e.g. the Ezekiel Mural “The well of the wilderness: Moses gives water to
the tribes” at Dura Europos synagogue in III century or wall painting in The
Tombs of Beni Hasan around 1900 B.C. with drawings of fishing and even
humorously rescuing a human who fell in water 2. Many centuries later, J.
Ruskin – an English writer, philosopher, art critic and polymath – writes 3:
“We want not to see water anatomized … After the entire failure of all artists
that ever lived before Turner in land and skies, we are prepared to find that
they had not the least idea of water. When they thought they painted water, in
fact, they were like ‘those happier children, sliding on dry ground,’ and had
not the chance of wetting a foot.” In his opinion, J.M.W. Turner is “the real
Triton of the sea, as he was Titan of the earth”.

Figure 1: “The Ninth Wave” (1850) by Ivan Aivazovsky (1817-1900), The State
Russian Museum, Mikhailovsky Palace.

Seascapes, in fact, serve as a good example where artists are prone to
mistakes and illustrate the fundamental reason for that: we will use one of
them to introduce the reader to the subject that this article is dedicated to.
In particular, breaking waves have always been the focus of maritime painters
as they reflect the destructive power of nature with one of the most famous
and analyzed 4; 5; 6 paintings of this sort being “The Great Wave of Kanagawa”
by Hokusai. Another example illustrated in Fig. 1 is Aivazovsky’s depiction of
“The Ninth wave” – a colossal wave larger than all the others and the subject
of legend among sailors – formed in the middle of the sea about to crash on
sailors survived after a shipwreck; however, the artist painted from the
shore! Such a depicted wave can form only close to the coast, but not in the
open sea, and is called a shore-break. As we now know from water wave theory,
when a wave travels across the open sea, it gains speed, but when it reaches a
shallow coastline, it begins to slow down due to the friction with the shallow
bottom causing the wave to break. Later in his life, after painting about 3000
seascapes, in his autobiography Aivazovsky wrote: “A painter, who just copies
nature, becomes its slave, with his limbs tied up … The living elements (of
nature) are elusive to artist’s brush.” One of his latest paintings “The Wave”
(1898) amazes with the details of the portrayed complexity of sea waves.

##  II “Imagine a child playing beside a pool of water, testing to see how
close he can place his finger to the water surface” (Cortat & Miklavcic7)

Let us start with depiction of still water, which presumably should be easier
to render. For accuracy, we turn to academism, a favorite representative of
which is William-Adolphe Bouguereau, who achieved a level of technical skill
that was virtually unparalleled by his colleagues.

Figure 2: “L’Aurore” (Dawn) (1881) by William-Adolphe Bouguereau, Birmingham
Museum of Art, Birmingham, Alabama, US.

In one of his most acclaimed and reverential works, “L’Aurore” in Fig. 2, the
water is like a mirror lightly perspired from a morning dew; as we can see,
the reflection of a twin toe is rising to the surface to meet the other, and
hence to touch the water surface. The artist’s idea behind this part of the
painting was to convey the sensitive, unbearable lightness of being.
Allegorically, for scientists this painting illustrates probing nature of
water by carefully approaching it. Leaving aside that L’Aurore defies gravity,
while water does not, the first look at the canvas brings up the question 7
formulated in this section title and, under closer examination, shows
impossibility of water to behave this way as there is no capillary rise where
her toe seems to be touching the water surface. Obviously, the answer is
dictated by the van der Waals interaction between finger/toe and water
separated by the distance $h$ leading to the force per area $S$ characterized
by the Hamaker constant $A_{H}\approx 10^{-20}\,\mathrm{J}$:

$\displaystyle\frac{F_{\mathrm{vdW}}}{S}=\frac{A_{H}}{6\pi h^{3}},$

which should be balanced by the gravity force $F_{g}$ acting on the bulge of
water of height $h$ and of base area $S$ as well as the surface tension force
$F_{\gamma}=\gamma/R$ tending to flatten the bulge, where the radius of
curvature $R\approx S/(2\pi h)$. As a result,

$\displaystyle\frac{F_{\mathrm{vdW}}}{S}=\frac{F_{g}}{S}+\frac{F_{\gamma}}{S}\
\Rightarrow\ \frac{A_{H}}{6\pi h^{3}}=\rho gh+\frac{2\pi\gamma}{S}h.$

Taking the area $S=10^{-4}\,\mathrm{m^{2}}$ to be that of the finger/toe tip,
cf. Fig. 9 in Appendix, we can see that $F_{g}$ and $F_{\gamma}$ are of the
same order and hence the minimal separation distance $h$ can be estimated via
$h\approx\left(A_{H}/6\pi\rho g\right)^{1/4}=1\,\mathrm{\mu m}$. However, this
distance is much smaller than the normal tremor amplitude of human limbs 8,
which starts at $O(100)\,\mathrm{\mu m}$ and increases to $O(1)\,\mathrm{cm}$
with time of continuous extension of a limb ($\sim 1\,\mathrm{hr}$);
presumably, the model had to pose much longer than that and hence touching
water surface is unavoidable leading to capillary wetting. This vicious circle
only means that Bouguereau’s representation is physically impossible despite
his realistic genre. In fact, his other painting “The Wave” suffers from a
similar issue: the female figure appears to be completely dry despite the
intense sea waves. Speaking of tremor in limbs, we can, albeit humorously,
mention that the object of his painting is supposedly healthy, though science
is capable of spotting neurological disorders in paintings as in the famous
“Christina’s World” by A. Wyeth 9. Same as how osteoarthritis, known to
increase tremor, in the hands of Michelangelo Buonarroti was spotted from his
portrait (c. 1535) by Jacopino del Conte 10.

##  III “Bacchus hath drowned more men than Neptune” (Thomas Fuller)

Figure 3: “Bacchus” (c. 1596) by Michelangelo Merisi da Caravaggio
(1571-1610), Uffizi, Florence.

First look at the stark, dramatically lit super-realistic painting by
Caravaggio, known for his uncunny level of detail, cf. Fig. 3, suggests that
the concentric waves in the wine glass – likely concocted by the artist with
the idea to add some life to the painting as the flat wine surface is
otherwise unappealing – could be due to shaking hands of Bacchus, possibly
caused by excessive drinking this deity is known for. Caravaggio himself was a
heavy drinker, a womaniser, and of violent character which led him to
committing three murders. This, nevertheless, did not prevent him from being
one of the greatest artists of all times able to master an unprecedented
accuracy. Thanks to the latter we can see that the wavelength $\lambda\approx
1\,\mathrm{cm}$ is shorter than the depth of wine in the glass. Hence, we can
determine the wave period from the standard relation for the gravity-driven
waves, $\omega=(2\pi g/\lambda)^{1/2}$, which gives circular frequency
$\omega=\mathcal{O}(10^{2})\,\mathrm{rad\,s^{-1}}$. Should the waves be
generated from a point source – e.g. by dipping and withdrawing a finger,
which would be awkward to achieve as both hands of Bacchus are occupied – then
we know that the wavelength must be changing with the distance from the
origin. Namely, according to the following expression for the wavelength
$\lambda=8\pi r^{2}/(9t^{2})$ for a given instant of time $t$ the wavelength
must increase as $r^{2}$, which contradicts the painting, in which $\lambda$
does not change with the distance from the origin! Also, the waves are
obviously shown for the time when they were able to propagate to the glass
wall of radius $R$ at the wine level, that is from the wave propagation limit
$R\sim gt^{2}/4$ we find $t>\left(4R/g\right)^{1/2}$.

Hence, we conclude that the depicted waves, if true, must have a different
origin rather than from a finger dipped in the glass. Another potential source
of wave generation could be due to vibration of the wine glass itself, induced
for example by hitting it against the table or some other object (or rubbing a
moistened finger around the wineglass rim). Given the glass radius $R\approx
10\,\mathrm{cm}$, thickness $a\approx 2\,\mathrm{mm}$, density per unit area
$\rho_{g}\approx 5\,\mathrm{kg\,m^{-2}}$, Young modulus $E\approx
48\,\mathrm{GPa}$, Poisson modulus $\nu=0.22$, we can estimate the circular
frequency of axisymmetric waves (formula is for a circular plate as the glass
is almost flat):

$\displaystyle\omega_{g}=\frac{\lambda^{2}}{R^{2}}\left(\frac{D}{\rho_{g}}\right)^{1/2},\
D=\frac{Ea^{3}}{12(1-\nu^{2})},$

which yields for $\lambda^{2}=10-10^{3}$:
$\omega_{g}=10^{3}-10^{5}\,\mathrm{rad\,s^{-1}}$. Obviously, this frequency is
incommensurate with the frequency of the observed water waves
$\omega=\mathcal{O}(10^{2})\,\mathrm{rad\,s^{-1}}$ calculated above, not to
mention that these waves would have to propagate from the glass perimeter
towards the center and hence wavelength must change and becomes shorter with a
spike in the middle. Hence, this is not an option as well.

Finally, sloshing wine in the glass would not be able to induce symmetric
waves of a uniform length, but would be a more convincing way to add life to
the painting. Altogether, one can conclude that the depicted waves are
inconsistent with the physics of water waves.

##  IV “The Vikings’ voyage was long and dangerous” (Belikov & Knyazeva 11)

Sidebar: His name is attached to Roerich Pact, which is an inter-American
treaty protecting cultural objects from destruction for military purposes.
Notably, Roerich also designed religious art for churches throughout both
Russia and Ukraine.

Laminar water waves on a larger scale have been depicted numerous times, but a
notable artwork in this regard is by Nikolai Roerich on semi-fantastic motifs
about the past Russia, shown in Fig. 4. The canvas portrays a caravan of
Viking ships sailed by merchants from distant countries. Obviously, they are
not defenseless as on the rook sides there are large shields hanging, while
inside the rooks people are shackled in iron armor.

The painting offers a number of phenomena suitable for basic estimates; the
analysis here follows that by A. Stasenko 12, though the conclusions are
somewhat different. First, the waves in front of the boat are depicted as
‘supersonic’ relative to the boat speed as they propagate ahead of the boat.
By comparing with the head of Varangians, one can estimate their wavelength as
$\lambda=35\,\mathrm{cm}$. Since, clearly, the wave propagation happens in the
deep water regime, their speed can be calculated from the corresponding
gravity-driven limit, $u_{g}=(g\lambda/2\pi)^{1/2}=0.7\,\mathrm{m\,s^{-1}}$.
Second, the boat speed $v$ can be evaluated from the height $h$ of the bow
wave based on the energy conservation for a mass of water $m$ moving with
velocity $v_{\perp}=v\sin{\alpha}$ towards to the boat wall at the bow of
angle $2\alpha$: $mgh\leq mv_{\perp}^{2}/2$. Estimating the bow wave height
from the painting to be $h\approx 25\,\mathrm{cm}$ and the bow angle
$\alpha\approx\frac{\pi}{2}$, we conclude that
$v\geq\left(2gh\right)^{1/2}/\sin{\alpha}\approx 2.3\,\mathrm{m\,s^{-1}}$.
Obviously, there is a contradiction to the observation made earlier that the
water waves are ‘supersonic’. Given that the reflections on the waves ahead of
the boat seem to be quite realistic suggesting that they were copied from a
real model, one might guess that it is the wave height at the bow, which was
misrepresented for an artistic purpose. The conundrum resolution depends on
which part of the painting you trust. The artist’s idea with waves propagating
ahead of the boat was to project forceful arrival of Vikings. However, Roerich
cannot be excused for making the phenomena unrealistic since linear water wave
theory has already existed by that time due to George Biddell Airy (1841).

Figure 4: “Merchants from Overseas” (1901) is one of the most recognizable
paintings by Nicholas Roerich (1874-1947), Tretyakov’s gallery in Moscow.

##  V “What care I that the virtue of some sixteen-year-old maid was the price
for Ingres’ La Source? That the model died of drink and disease in the
hospital is nothing when compared with the essential that I should have La
Source, that exquisite dream of innocence.” (George Moore)

Figure 5: “La Source” (1820-1856) by Jean-Auguste-Dominique Ingres
(1780-1867), Musée d’Orsay, Paris.

Waves are ubiquitous manifestations of fluid dynamics not only on horizontal
water surfaces discussed above, but also on flowing jets and streams, as we
can see in Ingres’ painting, cf. Fig. 5. A reader might say: “Enough nudes!”
Gustav Klimt would have objected “All art is erotic”. It just happen that
nudes constitute a substantial part of art, especially from that period. And
there is a good reason for that, i.e. to show people outside of their
historical context. The work on “La Source” was begun in Florence around 1820
and completed in Paris around 1856, when he was seventy-six years old, already
being the famous president of the École des Beaux-Arts. A water spring, which
in French is la source, is represented here in the guise of a girl pouring
water from a jug; in classical literature, a source is also sacred to the
Muses of poetic inspiration. The young daughter of Ingres’ concierge served as
a model for the painting.

Despite the long period of work on this painting, there are a few issues which
become apparent due to meticulous attention to details provided by the artist.
First, by noting that there is a well-defined circle – an expanding wave on
the water pond surface – from general theory of gravity-driven waves we know
that its radius is $r=g(t_{d}-t_{i})^{2}/4$ where $t_{d}$ is the movement of
depiction and $t_{i}$ the moment of the first jet impact on the pond, which
gives us $t_{d}-t_{i}\approx 0.2\,\mathrm{s}$, i.e. $t_{d}\approx
0.6\,\mathrm{s}$. The moment $t_{i}$ follows from the basic motion under
gravity giving $t_{i}\approx 0.4\,\mathrm{s}$; the corresponding jet velocity
is $V_{i}\approx 5\,\mathrm{m/s}$ after it travels $\approx 1\,\mathrm{m}$.
Finally, the water exit velocity $V_{0}$ from the pitcher at the moment of
depiction follows from Bernoulli’s equation $V_{0}=(2gh)^{1/2}$, where for $h$
we estimate $h=0.2\,\mathrm{m}$ thus giving $V_{0}\approx 2\,\mathrm{m/s}$.

A close look reveals that there are three jets: one on the left of the nude’s
hand (viewer’s point of view) reaches the pond, while the other two enveloping
her hand traveled only about $0.5\,\mathrm{m}$. The only explanation for the
different lengths of these jets is that the pitcher was tilted more some time
after initial pouring started in order to overflow the obstruction due to the
hand. The jet dynamics raises even more troubling issues. For that notice that
the left jet diameter $d_{1}\approx 1.5\,\mathrm{cm}$ based on the comparison
with the model’s iris. At the exit from the pitcher the same jet is
approximately elliptic with two major axes $a$ and $b$, so that the equivalent
radius $d_{e}=(ab)^{1/2}\approx 2\,\mathrm{cm}$ at the top of the jet. The
axis-switching wavelength varies approximately from $\lambda_{1}\approx
5\,\mathrm{cm}$ at the top of the jet to $\lambda_{1}\approx
8-10\,\mathrm{cm}$ at the bottom. The axis-switching of elliptic jets studied
by Lord Rayleigh 13 and recently improved to account for gravity 14 yields for
the circular frequency of axis-switching:

$\displaystyle\omega=\omega_{0}\left(1+\frac{2gz}{V_{0}^{2}}\right)^{1/4},\
\omega_{0}=2\sqrt{3}\left(\frac{\sigma}{\rho}\right)^{1/2}d_{e}^{-3/2},$

where $d_{e}$ is an equivalent diameter of the jet, $z$ is the distance from
the exit from pitcher, and $V_{0}$ the jet exit velocity. The wavelength is
then found from $\lambda(z)=2\pi V(z)/\omega(z)$. Plugging in numbers shows
that the instability wavelength varies from $\lambda\approx 0.15\,\mathrm{m}$
at the top to $\approx 0.5\,\mathrm{m}$ at the bottom, which agrees with the
laboratory measurements 14. Thus, while Ingres correctly captured the
qualitative feature that wavelength of axis-switching must increase along the
jet travel direction, the depicted wavelength is several times shorter! This
indicates that he likely used a scale-down experiment as a model for his
painting. More importantly, the Weber $We=\frac{\rho V^{2}d}{\sigma}$ and
Reynolds $Re=\frac{\rho Vd}{\mu}$ numbers for the jet on the left for $d$
varying from $d_{1}$ to $d_{e}$:

$\displaystyle We\approx(1-4)\cdot 10^{5},\ Re\approx 5\cdot 10^{6},$

imply that the jet is far in the unstable regime accompanied by jet breakup
and fish-bone structures! Presumably, Ingres knew that the jet breaks up into
droplets, but likely found it unappealing. Notably, its namesake by Gustave
Courbet (1862) does not suffer from incorrect depiction of water jet.

##  VI “We are meaner than flies; flies have their virtues, we are nothing but
bubbles.” (Satyricon, Petronius, 1st century C.E.)

The symbolism of Homo bulla (man is a bubble) was already proverbial by the
first century B.C. (“De Re Rustica”, Varro, 116-27 B.C.). In the baroque and
rococo eras soap bubbles conveyed allegorically the fleeting nature of earthly
pleasures and to remind the viewers of the fragility of the human existence
15; 16. As they have been a source of inspiration to artists, so for
scientists as well: “Make a soap bubble and observe it; you could spend a
whole life studying it,” (Lord Kelvin).

Figure 6: “Soap bubbles” (Les Bulles de Savon) (1733-1734) by Jean-Baptiste-
Siméon Chardin, Metropolitan Museum of Art, New York City. A version of it was
exhibited at the Paris Salon of 1739, where it was widely acclaimed.

Jean-Siméon Chardin111Besides him, soap bubbles appear in works of many other
artists in various countries including contemporaries Gerrit Dou, Karel
Dujardin, Charles Van Loo, Jan Steen, Edward Wilson, John Everett Millais,
Thomas Couture, Edouard Manet, Hans Heyerdahl, to name a few. Not so long time
ago, the art world was shaken with the stolen painting “Child with a Soap
Bubble” in 1999, which is believed to belong to Rembrandt. The painting was
returned to the municipal museum of Draguignan in March, 2014. transformed the
genre of still-life painting and contributed to making France the origin of a
new dimension of critical understanding in art, a new feeling of life, in
which individual intelligence was opposed to wealth and power. Compared to
still life painters from the previous century such as Frans Hals, he added
something more profound such a smile shedding light on hypocrisy, the mood of
quiet absorption, so familiar to scientists, that can lead to new knowledge.
Chardin’s picture in Fig. 6 definitely hints at the transient nature of life.
Such images resonate in the context of the period’s growing interest in the
experience of childhood, culminating in Jean-Jacques Rousseau’s influential
treatise “Émile”, or “on Education” (1762) about the right way to raise
children. Rousseau himself, however, betrayed his own and sent them to
orphanage immediately upon birth without ever knowing or seeing them. The
price to pay for genius?

First, from the painting in figure 6 one learns that split straw was used to
create soap bubbles. One may ask why is the bubble conspicuously elongated and
why is it not hanging much down due to gravity compared to paintings of other
artists? Careful examination of the canvas suggests that with the average
radius $R\approx 10\,\mathrm{cm}$ the main axes of the ellipsoid bubble differ
by about $1\,\mathrm{cm}$, i.e. $\delta R\equiv
R_{\mathrm{max}}-R_{\mathrm{min}}\approx 1\,\mathrm{cm}$. Also, the bubble
main axis is at about $\theta_{0}=\frac{\pi}{6}$ from the vertical axis and
deviates from that of the straw by about $\delta\theta\approx
0.1\,\mathrm{rad}$. Hence, the bubble volume $V_{b}\approx 4\cdot
10^{-3}\,\mathrm{m^{3}}$. We also need the soap film thickness. Since we
cannot see colors on the bubble, it must be above the visible light
wavelength. At the same time, given inefficiency of the old method to produce
soap bubbles, it is likely that the depicted soap bubble is close to bursting,
so we can take for its thickness $h=1\,\mathrm{\mu m}$, which gives for the
volume of water in the soap film $V_{\mathrm{SF}}=4\pi R^{2}h\approx
10^{-7}\,\mathrm{m^{3}}$ and mass $m_{\mathrm{SF}}\approx 0.1\,\mathrm{g}$.
Also, the soap was not as efficient as the one produced commercially nowadays
such as sodium dodecyl sulfate, which lowers surface tension of water from
$0.07\,\mathrm{N\,m^{-1}}$ to about $0.03\,\mathrm{N\,m^{-1}}$; hence we take
for the surface tension of the soap solution in the painting
$\sigma=0.05\,\mathrm{N\,m^{-1}}$.

Next, applying Bernouilli’s equation to the streamline impinging the front of
the soap bubble, cf. Fig. 10 in Appendix, we can relate the stagnation
pressure $p_{s}$ at the umbilic of the bubble with that at the exit of the
straw $p_{e}$, i.e. $p_{s}=p_{e}+\frac{1}{2}\rho u_{e}^{2}$. Also, from the
Young-Laplace law we know the pressure jump across the bubble interface at the
stagnation point $s$:

$\displaystyle p_{s}-p_{0}=\frac{4\sigma}{R-\frac{3}{2}\delta
R}\approx\frac{4\sigma}{R}\left(1+\frac{3}{2}\frac{\delta R}{R}\right).$

Taking into account that should there be no flow, $u=0$, the bubble would be
round, $p_{e}-p_{0}=\frac{4\sigma}{R}$. Altogether, this gives the expression
for the ejection velocity $u^{2}=\frac{12\sigma}{\rho}\frac{\delta R}{R^{2}}$,
which yields $u\approx 1\,\mathrm{m\,s^{-1}}$.

Finally, considering the horizontal ($x$) and vertical ($y$) components of the
force acting on the bubble,

$\displaystyle F_{x}\approx p_{s}S\sin{\theta_{0}},\ F_{y}\approx
p_{s}S\cos{\theta_{0}}+\delta m\,g,$

where $S=\pi R^{2}$ is the cross-section of the bubble and $\delta
m\,g=\delta\rho\,V_{b}\,g$ is the buoyancy force due to the density difference
between air inside and outside the bubble, we get

$\displaystyle\tan{(\theta_{0}+\delta\theta)}=\frac{F_{x}}{F_{y}}\
\Rightarrow\ \delta\rho=\frac{2\,\delta\theta\,p_{s}S}{gV_{b}}\approx
10^{-2}\,\mathrm{kg\,m^{-3}},$

i.e. the air is a bit heavier inside the bubble, which is due to exhaled
microdroplets of water (evaporated water, not in the form of microdroplets, in
fact makes air lighter!) and $\mathrm{CO_{2}}$. Of course, the bubble
equilibrium could be easily swayed by the convecting air, which was not taken
into account. The consistency of the calculations with the artist depiction
suggests that there was no significant flow of air at the moment of
observations, which is perhaps true for the sunny afternoon in French Riviera.
Chardin produced quite a number of paintings with soap bubbles, which means
that he inevitably studied them well and hence depicted accurately.

##  VII Summary

One might suppose that the difference between the artist and the scientist is
that the artist seeks beauty, while the scientist seeks truth. While S.
Chandrasekhar argued 18 that this is not the case, in this work we sought both
truth and beauty in artwork. The title of this article is (a translation by A.
Shafarenko) from poetry of Arseny Tarkovsky – a Soviet and Russian poet, who
was arrested and sentenced to execution for publishing a poem which contained
an acrostic about Lenin, but managed to escape – and encapsulates the richness
of the language one may speak with fluid phenomena.

While fluid phenomena themselves are a great source of inspiration for artists
19; 20, as we saw from the analyzed paintings, artists often used fluids as a
conceptual means to communicate certain emotions: Aivazovsky’s “The Ninth
Wave” (Late Romanticism) – destructive power of nature, Roerich’s “Merchants
from Overseas” (Russian realism) – forcefulness of incoming Vikings, Ingres’
“La Source” (Neoclassicism) – poetic inspiration, Chardin’s “Soap bubbles”
(Rococo) – fragility of human existence, Bouguerau’s “L’Aurore” (Academism) –
unbearable lightness of being, Caravaggio’s “Bacchus” (Baroque) – the
liveliness of wine. While in some cases, such as “Bacchus”, one can say that
the artist corrected nature, in most of fluid phenomena misrepresentations the
mistakes were genuine based on the lack of understanding of fluid phenomena;
otherwise it would go against aspirations to portray the world truthfully as
discussed at the beginning of the article. However, all the contortions could
not mar the harmonious beauty of the canvas. Besides debunking the
scientifically troubled renderings of fluid phenomena, we also tried to
highlight sometimes tragic background or hypocritical motivation: the
deceptive nature of aesthetically valued artwork goes well beyond picturing
fluid dynamics.

##  VIII Instead of an afterword

To reinforce this point of view, and to append our telegrammatic excursion in
the wide range of epochs and movements, we invite the reader to ponder over a
few more paintings from impressionism and post-impressionism, to contrast
with, say, academism. In fact, the Impressionists had the gall to open their
first exhibition in 1874 in opposition to ingrained academic favorites like
Bouguereau, who was then at the height of his career as a Salon master.
Turning to Van Gogh, whose expressive exaggerations are perhaps most clearly
seen in his last works, painted in 1889, during his stay in a psychiatric
facility hospital in Saint-Rémy, near Arles. Both paintings are shown in Fig.
7 and fluid depictions there aim to deliver vivacity of the sky in van Gogh’s
“The Starry Night” as well as human addiction in his “Self-portrait with
bandaged ear and pipe”.“The Starry Night”, which was a view from Vincent’s
window, distorted to express the artist’s existential turmoils. In both
paintings van Gogh tried to capture the essence of the apparently turbulent
fluid flows, which we formally call coherent structures nowadays.
Impressionism and coherent structures, discovered in science only relatively
recently 21, go hand-in-hand despite that the idea of depicting complex
phenomena via the key features goes back at least to Leonardo da Vinci, who
was engaged with fluid flows throughout his life 22; 23. It has been noticed
by fluid mechanicians 24; 25 that Leonardo “seems to be thinking about ways of
separating flow into steady and turbulent components” and presage the concept
of coherent structures 26; 27 as the words eddies and eddying motions in
turbulenza percolate throughout Leonardo’s treatise of liquid flows222“Observe
the motion of the surface of the water, which resembles that of hair, which
has two motions, of which one is caused by the weight of the hair, the other
by the direction of the curls; thus the water has eddying motions, one part of
which is due to the principal current, the other to the random and reverse
motion.” (Leonardo, c. 1500).

Figure 7: (left) “The Starry Night” (1889), Museum of Modern Art, New York
City. (right) “Self-portrait with Bandaged Ear and Pipe” (1889), Courtauld
Institute, London by Vincent van Gogh (1853-1890).

One would expect that impressionists, including van Gogh, could have learned
not only from da Vinci, but also benefited from the contemporary science 29.
The question to ask is if in “The Starry Night” the artist respected the
regular structure of van Karman street by arbitrarily terminating the von
Karman vortices or, may be, he captured the transient stage of the street
formation. Compare this to the famous photographs of the Kelvin-Helmholtz
instability, e.g. in “The Album of Fluid Motion”19. Another canvas, cf. Fig.
7(b), is his self-portrait.

Leaving aside that in reality he cut his left ear, the smoke coming out of the
pipe is able to develop fanciful large scale instabilities (vortices) a short
distance from the pipe: compare to the smoke plume photographed by Werner
Wolff – Black Star 30. In the latter case, the combustion and release of hot,
opaque compounds form a column of smoke, warmer than the surrounding air,
which develops upward, forming a laminar flow. At some point, however, free
convection weakens as the smoke column cools down and thus cannot maintain its
stability by advecting irregularly in the surrounding air forming eddies.

Figure 8: “The Scream” (1893) by Edvard Munch (1863-1944), National Gallery of
Norway.

Echoed van Gogh’s mood is the final powerful canvas “The Scream” of our
discussion shown in Fig. 8, symbolizing the anxiety of the human condition or
distress based on a personal experience of the artist described in his
autobiographical prose poem333 I was walking along the road with two friends.
The sun was setting. I felt a breath of melancholy
–Suddenly the sky turned blood-red. I stopped,
and leaned against the railing, deathly tired–
looking out across the flaming clouds that hung
like blood and a sword over the blue-black fjord
and town. My friends walked on–I stood there,
trembling with fear. And I sensed a great, infinite
scream pass through nature. However, the figure in the painting is not just
“sensing” but appears to produce the scream with the reverberation of the
sound throughout the vast landscape suggested by the long curving strokes.

While we see tall ships in the background, it is an invisible fluid dynamics –
acoustics – which is our focus now. It is believed to be inspired by the 1883
Krakatoa eruption, which exploded with such a force that it destroyed its
island and atoll, releasing 20 million tons of ash into the atmosphere that it
made sunsets reddish for years. It also produced the loudest “scream” ever
recorded: measuring instruments 160 km (100 mi) away from Krakatoa picked up
an immense 172 decibels, quite loud considering the distance. For comparison,
humans start to experience acute pain at 125 decibels.

Given the inverse square law for the sound pressure $L$ reduction with
distance $R$ from a source

$\displaystyle L_{2}-L_{1}=-20\log{\frac{R_{2}}{R_{1}}},$

we can find the range of sound level $L_{1}$ of the screamer given that the
two strangers behind the screamer do not cover their ears. Taking
$R_{1}=0.1\,\mathrm{m}$ as the distance from his mouth to ear and
$R_{2}=10\,\mathrm{m}$, we find $L_{2}-L_{1}=-40\,\mathrm{db}$, i.e. since
$\max{L_{2}}=125\,\mathrm{db}$, we get $L_{1}\leq 165\,\mathrm{db}$. If one
naturally assumes that the screamer covers his ears due to acute pain, then
$\min{L_{1}}=125\,\mathrm{db}$ and hence at the strangers location
$85\,\mathrm{db}\leq L_{2}\leq 125\,\mathrm{db}$. This is a believable range
as the world’s loudest shouting record with $121.7\,\mathrm{db}$ (equivalent
to a jet engine) was produced by Irish teacher Annalisa Flanagan.

###### Acknowledgements.

The author is grateful to Professor A. Stasenko at my Alma Mater (Moscow
Institute of Physics and Technology), whose article 12, originally published
in Russian in 1973, sparked a life-long interest in the subject of depicting
fluid phenomena in art – the illustrations of which have been used by the
author in fluid mechanics classes both at the University of California at
Santa Barbara and University of Alberta. Special thanks to Dr. A. Zelnikov for
many stimulating discussions as well as Dr. A. Herczynski for his feedback on
an earlier version of this manuscript.

*

##  Appendix A Auxiliary figures

Figure 9: On critical distance. Figure 10: Air flow in the bubble.

## References

  * Ruskin (1907) J. Ruskin,  _Elements of drawing and perpective_ (J. M. Dent & Sons, 1907).
  * Gombrich (1951) E. H. Gombrich,  _The story of art_ (Phaidon Press, 1951).
  * Ruskin (1984) J. Ruskin,  _The critical heritage_ (Routledge, 1984).
  * Cartwright and Nakamura (2009) J. H. E. Cartwright and H. Nakamura, “What kind of a wave is Hokusai’s great wave off Kanagawa?” Notes Rec. R. Soc. 63, 119–135 (2009).
  * Dudley  _et al._ (2013) J. M. Dudley, V. Sarano, and F. Dias, “On Hokusai’s great wave off Kanagawa: localization, linearity and a rogue wave in sub-antarctic waters,” Notes Rec. R. Soc. 67, 159–164 (2013).
  * Ornes (2014) S. Ornes, “Science and culture: dissecting the great wave,” PNAS 111, 13245 (2014).
  * Cortat and Miklavcic (2003) F. P. A. Cortat and S. J. Miklavcic, “How closely can a solid approach an air-water surface without becoming wet?” Phys. Rev. E 68, 052601 (2003).
  * Stiles (1976) R. N. Stiles, “Frequency and displacement amplitude relations for normal hand tremor,” J. Appl. Physiology 40, 44–54 (1976).
  * Patterson  _et al._ (2017) M. C. Patterson, T. B. Cole, E. Siegel, and P. A. Mackowiak, “A patient as art: Andrew Wyeth’s portrayal of Christina Olson’s neurologic disorder in Christina’s World,” J. Child Neurology 32, 647–649 (2017).
  * Lazzeri  _et al._ (2016) D. Lazzeri, M. F. Castello, M. Matucci-Cerinic, D. Lippi, and G. M. Weisz, “Osteoarthritis in the hands of Michelangelo Buonarroti,” Journal of the Royal Society of Medicine 109, 180–183 (2016).
  * Belikov and Knyazeva (1973) P. F. Belikov and V. P. Knyazeva,  _Nikolai Konstantonovich Roerich_ (Moskow, Molodaya Gvardiya, 1973).
  * Stasenko (1990) A. Stasenko, “It’s beautiful – but is it science?” Quantum 1, 8–12 (1990).
  * Rayleigh (1879) L. Rayleigh, “On the capillary phenomena of jets,” Proc. R. Soc. London 29, 71–97 (1879).
  * Ma and Zhu (2020) Y. Ma and D. Z. Zhu, “Axis switching of free-falling elliptical water jets,” J. Hydraul. Eng. 146, 06020009 (2020).
  * Behroozi (2008) F. Behroozi, “Soap bubbles in paintings: Art and science,” Am. J. Physics 76, 1087–1091 (2008).
  * Hosseini (2020) A. Hosseini, “Air bubbles. colors of light: Newton’s observations and Chardin’s representations,” Venti Journal 1, 75–85 (2020).
  * Note (1) Besides him, soap bubbles appear in works of many other artists in various countries including contemporaries Gerrit Dou, Karel Dujardin, Charles Van Loo, Jan Steen, Edward Wilson, John Everett Millais, Thomas Couture, Edouard Manet, Hans Heyerdahl, to name a few. Not so long time ago, the art world was shaken with the stolen painting “Child with a Soap Bubble” in 1999, which is believed to belong to Rembrandt. The painting was returned to the municipal museum of Draguignan in March, 2014\.
  * Chandrasekhar (1987) S. Chandrasekhar,  _Truth and beauty: aesthetics and motivations in science_ (The University of Chicago Press, 1987).
  * van Dyke (1982) M. van Dyke,  _An album of fluid motion_ (The Parabolic press, 1982).
  * Zabusky (2015) N. J. Zabusky, “Fluids in motion: Inspiration and realization for artists and stemists,” Leonardo 48, 138–146 (2015).
  * Brown and Roshko (1974) G. L. Brown and A. Roshko, “On density effects and large structure in turbulent mixing layers,” J. Fluid Mech. 64, 775–816 (1974).
  * MacCurdy (1938) E. MacCurdy,  _The Notebooks of Leonardo Da Vinci._ (Reynal & Hitchcock, 1938).
  * Macagno (2000) E. O. Macagno,  _Leonardian fluid mechanics: “Libro dell’Acqua” II. Monogr. 120_ (Iowa Inst. Hydraul. Res., Iowa City, IA, 2000).
  * Taylor (1974) G. I. Taylor, “The interaction between experiment and theory in fluid mechanics,” Annu. Rev. Fluid Mech. 6, 1–17 (1974).
  * Lumley (1992) J. Lumley, “Some comments on turbulence,” Phys. Fluids A 4, 203–211 (1992).
  * Gad-el-hak (1998) M. Gad-el-hak, “Fluid mechanics from the beginning to the third millennium,” Int. J. Engng. Ed. 14, 177–185 (1998).
  * Marusic and Broomhall (2021) I. Marusic and S. Broomhall, “Leonardo da Vinci and fluid mechanics,” Annu. Rev. Fluid Mech. 53, 1–25 (2021).
  * Note (2) “Observe the motion of the surface of the water, which resembles that of hair, which has two motions, of which one is caused by the weight of the hair, the other by the direction of the curls; thus the water has eddying motions, one part of which is due to the principal current, the other to the random and reverse motion.” (Leonardo, c. 1500).
  * Guébhard (1881) A. Guébhard, “Studies of vortex-rings,” Popular Science 20, 1–16 (1881).
  * Corrsin (1961) S. Corrsin, “Turbulent flow,” American Scientist 49, 300–325 (1961).
  * Note (3) I was walking along the road with two friends.
The sun was setting. I felt a breath of melancholy
–Suddenly the sky turned blood-red. I stopped,
and leaned against the railing, deathly tired–
looking out across the flaming clouds that hung
like blood and a sword over the blue-black fjord
and town. My friends walked on–I stood there,
trembling with fear. And I sensed a great, infinite
scream pass through nature.


# Differentiable Safe Controller Design through Control Barrier Functions

Shuo Yang∗, Shaoru Chen∗, Victor M. Preciado, Rahul Mangharam Shuo Yang,
Shaoru Chen, Victor M. Preciado, Rahul Mangharam are with the Department of
Electrical and Systems Engineering, University of Pennsylvania, Philadelphia,
PA 19104, USA. Email: {yangs1, srchen, preciado,
<EMAIL_ADDRESS>contributed equally.

###### Abstract

Learning-based controllers, such as neural network (NN) controllers, can show
high empirical performance but lack formal safety guarantees. To address this
issue, control barrier functions (CBFs) have been applied as a safety filter
to monitor and modify the outputs of learning-based controllers in order to
guarantee the safety of the closed-loop system. However, such modification can
be myopic with unpredictable long-term effects. In this work, we propose a
safe-by-construction NN controller which employs differentiable CBF-based
safety layers, and investigate the performance of safe-by-construction NN
controllers in learning-based control. Specifically, two formulations of
controllers are compared: one is projection-based and the other relies on our
proposed set-theoretic parameterization. Both methods demonstrate improved
closed-loop performance over using CBF as a separate safety filter in
numerical experiments.

###### Index Terms:

Safety-critical control, learning-based control, barrier functions, control
synthesis.

##  I Introduction

Learning-based control has become increasingly popular for controlling complex
dynamical systems [1] since it requires little expert knowledge and can be
carried out in an automatic, data-driven manner. However, due to the black-box
nature of most learning models, learning-based controllers such as neural
network (NN) controllers lack formal guarantees which significantly limit
their deployment in safety-critical applications.

The integration of control-theoretical approaches and machine learning has
provided a promising solution to safe learning control, where trainable
machine learning modules are embedded into a control framework that guarantees
the safety or stability of the dynamical system [2, 3]. Of wide applicability
is the control barrier function (CBF) framework [4, 5] which explicitly
specifies a safe control set and guards the system inside a safe invariant
set. This is achieved by constructing a CBF-based safety filter that projects
any reference control input (possibly generated by a NN controller) onto the
safe control set online. When a continuous-time control-affine system is
considered, such projection reduces to a convex quadratic program (QP) which
is referred to as CBF-QP. Due to its simplicity, flexibility, and formal
safety guarantees, CBFs have been applied in safe learning control with many
successful applications [6, 7].

(a) Gauge map-based safe NN controller architecture.

(b) CBF-QP-based safe NN controller architecture.

Figure 1: Safe-by-construction NN controllers that utilize CBFs to construct
differentiable safety layers (yellow blocks).

Compared with model predictive control (MPC) [8], which needs to handle a
possibly nonconvex optimization problem in the face of nonlinear dynamical
systems, CBF-QP is computationally efficient to solve online. However, unlike
MPC, the QP-based safety filter only operates in a minimally invasive manner,
i.e., it generates the safe control input closest (in the Euclidean norm) to
the reference control input, unaware of the long-term effects of its action.
This indicates that the effects of the safety filter on the performance of the
closed-loop system are hard to predict. Therefore, the application of the
safety filter may give rise to myopic controllers [9] that induce sub-par
performance in the long term.

To address the issue of myopic CBF-based safety filters, in this work we
propose to utilize CBF to construct safe-by-construction NN controllers that
allow end-to-end learning. Incorporating safety layers in the NN controller
allows the learning agent to take the effects of safety filters into account
during training in order to maximize long-term performance. Inspired by [10],
we design a differentiable safety layer using the gauge map which establishes
a bijection between the polytopic output set of a NN (e.g., an $\ell_{\infty}$
norm ball) and the safe control set characterized by CBFs. We denote the
proposed architecture as NN-gauge (Fig. 1(a)). We compare NN-gauge with an
alternative safe NN controller, NN-diff-QP, which consists of a NN followed by
a differentiable CBF-QP layer (Fig. 1(b)). In the online execution, NN-gauge
requires closed-form evaluation or solving a linear program (LP) while NN-
diff-QP solves a quadratic program. Both methods are significantly cheaper to
run than MPC.

###  I-A Related works

Safe controller design: The use of gauge map in safe learning control was
proposed in [10, 11] which only considers linear dynamics for which control
invariant sets and an interior safe control policy are achievable. In this
work, by proposing NN-gauge, we significantly extend the scope of this
framework to handle nonlinear dynamics with CBFs and overcome the arising
computational difficulties by applying an implicit interior policy
parameterization. The construction of NN-diff-QP naturally follows from the
use of differentiable optimization layers [12, 13, 14], and NN-diff-QP is
applied in [15] and [16] in the context of stochastic control and
reinforcement learning, respectively. CBFs are applied as penalty functions to
promote the safety of NN controllers in [17]. However, unlike NN-gauge or NN-
diff-QP, the resulting NN controllers do not have formal safety guarantees.

Differential CBF: Introducing learning modules in the parameterization of CBFs
can improve the feasibility and performance of the safety filter [18, 19] even
in the face of changing environments [20, 21]. These works focus on learning
or improving CBFs such that the safe control set is enlarged and the safety
filter can work better with a reference controller. Instead, in this paper, we
consider NN controller synthesis with a given CBF.

The contributions of the paper are summarized below:

  1. 1.

We propose a novel, differentiable, safe-by-construction NN controller using
CBF with efficient online evaluation.

  2. 2.

We compare the proposed controller with the projection-based safe NN
controller numerically, showing that both of them have similar performance
while the online evaluation of our proposed method is more computationally
efficient.

  3. 3.

We demonstrate that learning safe-by-construction NN controllers leads to
better long-term closed-loop performance than filtering a trained unsafe NN
controller.

In Section II, background on CBF is provided and the control problem is
introduced. The construction of NN-gauge and NN-diff-QP is shown in Section
III, followed by numerical examples in Section IV. Section V concludes the
paper.

##  II Preliminary and problem formulation

###  II-A System model

In this paper we are interested in a continuous-time nonlinear control-affine
system:

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

where $f$ and $g$ are locally Lipschitz, $x\in D\subseteq\mathbb{R}^{n}$ is
the state and $u\in\mathbb{R}^{m}$ is the control input subject to bounded
polytopic control constraints:

$u\in U=\\{u\in\mathbb{R}^{m}\mid A_{u}u\leq b_{u}\\}.$ (2)

###  II-B Control barrier functions

Safety can be framed in the context of enforcing set invariance in the state
space, i.e., the state cannot exit a safe set. The safe set $\mathcal{C}$ is
represented by the superlevel set of a continuously differentiable function
$h(x)$. The algebraic expressions of the safe set $\mathcal{C}$, boundary of
the safe set $\partial\mathcal{C}$, and interior of the safe set
$\text{Int}(\mathcal{C})$ are given by:

$\displaystyle\mathcal{C}$ $\displaystyle=\\{x\in\mathbb{R}^{n}:h(x)\geq
0\\},$ $\displaystyle\partial\mathcal{C}$ $\displaystyle=\\{x\in
D\subset\mathbb{R}^{n}:h(x)=0\\},$ (3) $\displaystyle\text{Int}(\mathcal{C})$
$\displaystyle=\\{x\in D\subset\mathbb{R}^{n}:h(x)>0\\}.$

For a locally Lipschitz continuous control law $u=k(x)$, we have that
$\dot{x}=f(x)+g(x)k(x)$ is locally Lipschitz continuous. Thus, for any initial
condition $x_{0}\in D$, there exists a maximum time interval of existence
$I(x_{0})=[0,\tau_{max})$, such that $x(t)$ is the unique solution to the
ordinary differential equation (1) on $I(x_{0})$. We frame the safety of
system (1) in terms of set invariance as shown below.

######  Definition 1.

(Forward invariance and safety) The set $\mathcal{C}$ is _forward invariant_
if for every $x_{0}\in\mathcal{C}$, $x(t)\in\mathcal{C}$ holds for
$x(0)=x_{0}$ and all $t\in I(x_{0})$. If $\mathcal{C}$ is forward invariant,
we say system (1) is safe.

To verify invariance of $\mathcal{C}$ under the control input constraints (2),
a control barrier function is constructed as a certificate which characterizes
the admissible set of control inputs that render $\mathcal{C}$ forward
invariant.

######  Definition 2.

(Extended class $\mathcal{K}_{\infty}$ function) A continuous function
$\alpha$: $\mathbb{R}\rightarrow\mathbb{R}$ is said to be an extended
$\mathcal{K}_{\infty}$ if it is strictly increasing and $\alpha(0)=0$.

######  Definition 3.

(Control barrier function) Let $\mathcal{C}\subset D\subset\mathbb{R}^{n}$ be
the superlevel set of a continuously differentiable function
$h:D\rightarrow\mathbb{R}$, then $h$ is a control barrier function if there
exists an extended class $\mathcal{K}_{\infty}$ function $\alpha(\cdot)$ such
that for the control system (1):

$\displaystyle\mathsf{sup}_{u\in U}[L_{f}h(x)+L_{g}h(x)u]\geq-\alpha(h(x)),$
(4)

for all $x\in D$ where $L_{f}h(x)$ and $L_{g}h(x)$ denote the Lie derivatives.

Given the CBF $h(x)$, the set of all control values that render $\mathcal{C}$
safe is given by:

$K_{cbf}(x)=\\{u\in U:L_{f}h(x)+L_{g}h(x)u+\alpha(h(x))\geq 0\\}$ (5)

which we denote as the safe control set. The following theorem shows that the
existence of a control barrier function implies that the control system (1) is
safe:

######  Theorem 1.

([4, Theorem 2]) Assume $h(x)$ is a CBF on $D\supset\mathcal{C}$ and
$\frac{\partial h}{\partial x}(x)\neq 0$ for all $x\in\partial\mathcal{C}$.
Then any Lipschitz continuous controller $u(x)$ such that $u(x)\in K_{cbf}(x)$
for all $x\in\mathcal{C}$ will render the set $\mathcal{C}$ forward invariant.

One important feature of $K_{cbf}(x)$ is that it is a polytope 111The
polytopic safe control set $K_{cbf}(x)$ is non-empty for all $x\in\mathcal{C}$
by definition of CBF. In this work, we assume a valid CBF for system (1) is
given and use it for NN controller design, although synthesizing CBFs can be
challenging itself and is an active area of research.  for all states. This
enables the construction of a QP-based safety filter that modifies any given
reference controller $k(x)$ in a minimally invasive fashion [5] as follows:

$\displaystyle u(x)=\underset{u}{\text{argmin}}$ $\displaystyle\quad\lVert
u-k(x)\rVert_{2}^{2}$ (6) subject to $\displaystyle\quad
L_{f}h(x)+L_{g}h(x)u+\alpha(h(x))\geq 0,$ $\displaystyle\quad u\in U.$

Although the safety filter (6) guarantees the forward invariance of the safe
set, it does not take into account the consequences of the projection on
future states and the performance of the closed-loop system. This issue is
inevitable when the reference controller $k(x)$ and the safety filter are
designed separately, and we propose to fix it by designing safe-by-
construction controllers that are amenable to any learning or optimization
framework. Particularly, in this work, we consider optimizing NN controllers
using modern machine learning solvers (such as stochastic gradient descent
(SGD) and Adam [22]) with known system dynamics.

###  II-C Problem formulation

Following the definition of safe control set (5), we define the set of safe
control policies as

$\displaystyle\Pi:=\\{\pi:\mathbb{R}^{n}\mapsto\mathbb{R}^{m}|\pi(x)\in
K_{cbf}(x),\forall x\in\mathcal{C}\\}.$ (7)

Our task is to design a controller for system (1) such that a performance
objective is optimized and the closed-loop system always stays inside the safe
set $\mathcal{C}$. In other words, we consider finding a policy $\pi(x)$ that
solves the following optimal control problem within a horizon $T<\tau_{max}$:

$\displaystyle\underset{\pi\in\Pi}{\text{minimize}}$
$\displaystyle\quad\mathbb{E}_{x(0)\in\mathcal{C}}\Big{[}\frac{1}{T}\int_{0}^{T}c(x(t),u(t))\,dt\Big{]}$
(8) subject to $\displaystyle\quad\dot{x}=f(x)+g(x)u,$ $\displaystyle\quad
u(t)=\pi(x(t)),$

where $\mathbb{E}_{x(0)\in\mathcal{C}}$ is the expectation with respect to the
initial state, $c(x(t),u(t))$ is the cost associated with occupying state
$x(t)$ and action $u(t)$. The cost function $c(x(t),u(t))$ could be of any
form and is problem-specific, e.g., it can formulate the penalty or barrier
functions of constraints on the state $x(t)$.

Despite the complex dynamics, cost functions, and safe policy constraints,
problem (8) can still be effectively approached by parameterizing an NN policy
and applying SGD/Adam which is empowered by automatic differentiation and
modern machine learning solvers [23]. This procedure is simple and has been
shown effective to synthesize high-performance NN controllers [24, 25]. Next,
we study different parameterizations of safe NN policies $\pi\in\Pi$ while
solving (8) and evaluate their performances. Notably, with safe NN policies,
the effects of the CBF-based safety layer are automatically taken into
account.

##  III Safe Controller Design

A natural way to construct a safe NN controller $\pi(x)$ is to restrict the
output of the controller into the safe control set (5) for all states. NN-
diff-QP (Fig. 1(b)) achieves this by concatenating a differentiable projection
layer which can be implemented in a NN using toolboxes such as cvxpylayers
[14] and qpth [12]. In this section, we propose a different parameterization
of a safe NN controller which achieves better online computational efficiency
using gauge maps.

###  III-A Gauge map

The authors in [10] observe that, while it is challenging to directly restrict
the output of a NN inside a general polytope, adding a hyperbolic tangent
activation layer easily constrains the NN output into the unit
$\ell_{\infty}$-norm ball $\mathbb{B}_{\infty}$. This motivates the
application of the gauge map that establishes a bijection between
$\mathbb{B}_{\infty}$ and a general polytope which in this paper we consider
as $K_{cbf}(x)$. The notion of gauge map is facilitated by the concept of
C-set.

######  Definition 4.

(C-set [26]) A C-set is a convex and compact set including the origin as an
interior point.

The gauge function (or Minkowski function) of a vector $v\in\mathbb{R}^{m}$
with respect to a C-set $Q\subset\mathbb{R}^{m}$ is given by

$\displaystyle\gamma_{Q}(v)=\mathsf{inf}\\{\lambda\geq 0|v\in\lambda Q\\}.$
(9)

When $Q$ is a polytopic C-set defined by $\\{w\in\mathbb{R}^{m}|F^{T}_{i}w\leq
g_{i},i=1,2,\cdots,r\\}$, the gauge function can be written in closed-form as

$\displaystyle\gamma_{Q}(v)=\mathsf{max}_{i}\\{F_{i}^{T}v/g_{i}\\}.$ (10)

For any $v\in\mathbb{B}_{\infty}$, the gauge map is defined as

$\displaystyle G(v|Q)=\frac{||v||_{\infty}}{\gamma_{Q}(v)}\cdot v$ (11)

which constructs a bijection (see [10, Lemma 1]) between the unit ball
$\mathbb{B}_{\infty}$ and the C-set $Q$. As shown in Fig. 2, all points
$v\in\mathbb{B}_{\infty}$ are mapped “proportionally” to the same level set of
the polytope $Q$ by the gauge map, while their projections onto $Q$ tend to
concentrate on the boundary of $Q$.

Figure 2: Illustration of the gauge map from the $\ell_{\infty}$ ball
$B_{\infty}$ to a polytopic set $Q$.

###  III-B Interior safe policy

To apply the gauge map as a safety layer that maps the output of an NN to the
safe control set $K_{cbf}(x)$, we have to first find an interior safe policy
$\pi_{int}(x)$ such that $\pi_{int}(x)\in\text{Int}(\mathcal{C})$ for all
$x\in\mathcal{C}$. This is necessary since $K_{cbf}(x)$ may not be a C-set,
and we need to shift $K_{cbf}(x)$ by $\pi_{int}(x)$ such that it is recentered
around the origin. An explicit construction of $\pi_{int}(x)$ is achievable
for linear dynamical systems through multi-parametric programming [11]. The
recent work [27] demonstrates how to extract an explicit safe policy from CBFs
for system (1). When $\pi_{int}(x)$ is available, we have that the shifted
safe control set $\hat{K}_{cbf}(x):=K_{cbf}(x)-\pi_{int}(x)$ is a C-set for
which we can apply the gauge map.

Although theoretically feasible, extracting the interior safe policy
$\pi_{int}(x)$ can be challenging when the dynamics in (1) becomes complex. To
address this issue, we propose an alternative method that implicitly
constructs an interior safe policy by choosing $\pi_{int}(x)$ as the Chebyshev
center [28] of $K_{cbf}(x)$. Specifically, we have $\pi_{int}(x)=u^{*}$ where
$u^{*}$ is the solution to the following LP:

$\displaystyle\underset{u,R\geq 0}{\text{maximize}}$ $\displaystyle\quad R$
(12) $\displaystyle\mathsf{subject\>to}$ $\displaystyle\quad-
L_{g}h(x)u+R\lVert-L_{g}h(x)\rVert_{2}$
$\displaystyle\quad\quad\quad\quad\quad\quad\quad\leq L_{f}h(x)+\alpha(h(x)),$
$\displaystyle\quad a_{u,i}^{\top}u+R\lVert a_{u,i}\rVert_{2}\leq b_{u,i},\
i=1,\cdots,d_{u},$

where $a_{u,i}$ denotes the $i$-th row of $A_{u}$, $b_{u,i}$ denotes the
$i$-th entry of $b_{u}$, and $d_{u}$ represents the number of linear
constraints defining $U$. The Chebyshev center formulation of $\pi_{int}(x)$
pushes $\pi_{int}(x)$ away from the boundary of $K_{cbf}(x)$. It also
facilitates training of the upstream NN, since it makes the target function
that the NN needs to learn smoother. By the validity of CBF, the safe control
set $K_{cbf}(x)$ is a non-empty polytope for all $x\in\mathcal{C}$. By
choosing $\pi_{int}(x)$ as the Chebyshev center from (12), we readily have
that $\pi_{int}(x)\in\text{Int}(K_{cbf}(x))$ for $R>0$ and the shifted safe
control set $\hat{K}_{cbf}(x)$ is a C-set.

###  III-C Control policy architecture

With an interior safe policy $\pi_{int}(x)$, we can now construct a safe NN
controller using the gauge map, as shown in the following theorem.

######  Theorem 2.

Let $\phi_{\theta}:\mathcal{C}\rightarrow\mathbb{B}_{\infty}$ be a neural
network parameterized by $\theta$ and $\pi_{int}(x)$ be an interior safe
policy. Then, for any system state $x$ in the set $\mathcal{C}$, the policy:

$\pi_{\theta}(x):=G(\phi_{\theta}(x)|\hat{K}_{cbf}(x))+\pi_{int}(x)$ (13)

has the following properties:

  1. 1.

$\pi_{\theta}$ is safe.

  2. 2.

The policy $\pi_{\theta}(\cdot)$ is trainable with respect to the NN
parameters $\theta$.

###### Proof.

1) By the construction of the gauge map, we have
$G(\phi_{\theta}(x)|\hat{K}_{cbf}(x))\in\hat{K}_{cbk}(x)$ which is the safe
control set $K_{cbf}(x)$ shifted by the interior safe control input
$\pi_{int}(x)$. Therefore, we have $\pi_{\theta}(x)\in K_{cbf}(x)$ for all
$x\in\mathcal{C}$, and conclude that $\pi_{\theta}(x)$ is safe.

2) As shown in [10, Theorem 1], automatic differentiation can be applied to
compute the subgradients of $\pi_{\theta}$ with respect to $\theta$, making it
possible to train $\pi_{\theta}$. ∎

The neural network $\phi_{\theta}$ can embed any architecture and learns the
residual control policy added to the interior safe policy $\pi_{int}(x)$. By
the construction of the gauge map, the applied controller $\pi_{\theta}(x)$
shown in (13) belongs to the safe control set $K_{cbf}(x)$, and the
performance of $\pi_{\theta}(x)$ is usually no worse than $\pi_{int}(x)$ after
training. The online evaluation of $\pi_{\theta}$ can be done in closed-form
if an explicit interior safe policy $\pi_{int}(x)$ is given, or by solving an
LP if the implicit construction (12) is used.

######  Remark 1.

Our analysis of safe NN controllers can be readily applied to incorporate high
order CBFs [29] in which case the safe control set $K_{cbf}(x)$ is still a
polytope in the control input. In addition to the CBF-based safe control sets,
NN-gauge and NN-diff-QP can easily encode other forms of polytopic safe
control sets such as $u\in\\{u\mid Gu\leq b\\}$ and $u\in\\{u\mid Fx+Gu\leq
b\\}$.

##  IV Numerical Experiments

In this section, we demonstrate the application of safe NN controllers in
adaptive cruise control and aircraft collision avoidance. The following
control methods are considered:

  1. 1.

MPC: Model predictive control could guarantee the safety of the closed-loop
system with good performance, but it is computationally costly to run online
when the horizon is large or the dynamics is nonlinear.

  2. 2.

NN: A feedforward NN controller is trained to optimizes (8) with regularizers
penalizing violations of safety constraints. No safety filter is applied in
the online evaluation of this NN controller. This method enjoys fast online
evaluation but suffers the risk of safety violations.

  3. 3.

NN-QP: The above NN controller is equipped with the CBF-QP safety filter
during online evaluation.

  4. 4.

NN-diff-QP and NN-gauge: The safe-by-construction NN controllers introduced in
Section III that are trained directly to optimize the control performance (8).

For NN-diff-QP, we use cvxpylayers [14] to construct the differentiable QP
layer. For NN-gauge, the implicit interior policy (12) is applied. All
training is performed using PyTorch on Google Colab with Adam [22] as the
optimizer.

###  IV-A Adaptive cruise control

Adaptive cruise control (ACC) is a common example to validate safe control
strategies [4, 30]. The control goal of ACC is to let the ego car achieve the
desired cruising speed while maintaining a safe distance from the leading car.
We consider the scenario where the ego car tries to follow the leading car on
a straight road. The dynamics of the problem is given by (model adapted from
[30]):

$\displaystyle\dot{x}(t)=\begin{bmatrix}0&1&0\\\ 0&-0.1&0\\\
0&-1&0\end{bmatrix}x(t)+\begin{bmatrix}0\\\ 2.5\\\ 0\end{bmatrix}u,$ (14)

where
$x=\begin{bmatrix}x_{1}&x_{2}&x_{3}\end{bmatrix}^{\top}=\begin{bmatrix}p_{f}&v_{f}&d\end{bmatrix}^{\top}$,
$p_{f}$ is the position, $v_{f}$ is the velocity of the ego car, $d$ is the
distance between the ego and leading cars, and $u\in[-1,1]$ is the control
input denoting the acceleration. To prevent collision between the two cars,
the following CBF is applied:

$\displaystyle h(x)=x_{3}-1.8x_{2}$ (15)

and the trajectory cost for the ego car is given by

$\displaystyle\texttt{Cost}=\int_{0}^{T}(0.01(x_{2}(t)-v_{f}^{*})^{2}+0.05u(t)^{2})dt$
(16)

with a desired cruising speed $v_{f}^{*}=30m/s$. The leading car travels at a
constant speed of $16m/s$, so we expect the ego car’s speed to converge to
$16m/s$ at the steady state since a speed greater than $16m/s$ will lead to a
violation of safe distance.

All NN controllers, namely NN, NN-gauge and NN-diff-QP, are trained to
optimize (8) with horizon $T=1$s by Adam with randomly sampled initial states.
The system dynamics (14) is discretized with sampling rate $\Delta t=0.1$s
during training. The trainable NN modules $\phi_{\theta}$ in NN, NN-gauge, and
NN-diff-QP have the same architecture.

We test all NN controllers on 5 randomly sampled initial states over a horizon
$T=20$s in order to evaluate their long-term performance. The results are
reported in Table I. For one of the testing initial condition
$x(0)=\begin{bmatrix}p_{f0}&v_{f0}&d_{0}\end{bmatrix}^{\top}=\begin{bmatrix}0&30&100\end{bmatrix}^{\top}$,
we plot the values of the CBF function and the velocity along the trajectory
of the ego car in Fig. 3(b) and Fig. 3(a), respectively.

We observe that both NN-gauge and NN-diff-QP achieve similar closed-loop
performance, comparable to MPC 222While MPC solves a finite horizon optimal
control problem to optimality, its closed-loop performance is not guaranteed
to be optimal. In Table I, NN-diff-QP achieved better performance than MPC in
the long term.. While NN achieves reasonable performance, it violates safety
constraints (Fig. 3(a)). Directly applying the CBF-QP safety filter on it
enforces safe control, but deteriorates the long-term closed-loop performance
of the NN controller as shown in Fig. 3(b) where the optimal behavior of the
ego car is supposed to have a steady-state velocity of $16m/s$. NN-gauge has
an edge over NN-diff-QP in training and online evaluation time due to its use
of LP-based safety layers.

TABLE I: Closed-loop performance comparison of controllers in the adaptive cruise control example. Average values from $5$ randomly sampled initial states are reported. | Safety | Trajectory | Training time | Solve
---|---|---|---|---
|  | cost | per epoch | time
MPC | Safe | 269.3 | N/A | 3.11s
NN | Unsafe | 640.3 | 0.03s | 0.04s
NN-QP | Safe | 818.2 | 0.03s | 0.36s
NN-diff-QP | Safe | 258.3 | 11.0s | 0.35s
NN-gauge | Safe | 270.9 | 8.3s | 0.28s

(a) CBF values.

(b) Velocity of the ego car.

Figure 3: CBF values (left) and velocity of the ego car (right) under
different controllers are evaluated in closed-loop for $20$s. A CBF value
below zero indicates unsafety, and the optimal behavior of the ego car is
expected to have a steady state velocity of $16m/s$ same as the leading car.

###  IV-B Aircraft collision avoidance

We apply our framework to the aircraft collision avoidance problem which is
adapted from [31]. Specifically, we consider a dynamical system with states
$x=\begin{bmatrix}x_{a}^{\top}&x_{b}^{\top}\end{bmatrix}^{\top}$, where
$x_{a}=\begin{bmatrix}p_{a,x}&p_{a,y}&\theta_{a}\end{bmatrix}^{\top}\in\mathbb{R}^{3}$
is the state of aircraft $a$ with $p_{a,x}$ and $p_{a,y}$ denoting the
position and $\theta_{a}$ denoting the orientation. The state $x_{b}$ of
aircraft $b$ is defined similarly. The control inputs are
$u=\begin{bmatrix}v_{a}&w_{a}&v_{b}&w_{b}\end{bmatrix}^{\top}\in\mathbb{R}^{4}$
where $v_{a}$ and $w_{a}$ are speed and turning rate of aircraft $a$,
respectively, and $v_{b},w_{b}$ are defined similarly. The dynamics of the
aircraft $a$ (and similarly for aircraft $b$) is given by:

$\displaystyle\dot{x}_{a}(t)=\begin{bmatrix}v_{a}(t)\text{cos}(\theta_{a}(t))\\\
v_{a}(t)\text{sin}(\theta_{a}(t))\\\ w_{a}(t)\end{bmatrix}.$ (17)

As shown in Fig. 4, our goal is to drive aircraft $a$ to the left and aircraft
$b$ to the right while maintaining a minimum safe distance of $0.5$ between
them. Aircrafts $a$ and $b$ try to stay close to $p_{a,x}=-5$ and $p_{b,x}=5$,
respectively. A quadratic cost function is defined accordingly over a horizon
$T$, and we adopt the constructive CBF developed in [31] to encode the safe
set in the state space which also considers the input constraint

$\displaystyle(v_{a},w_{a},v_{b},w_{b})\in[0.1,1]\times[-1,1]\times[0.1,1]\times[-1,1].$

Note that the minimum admissible velocity of aircraft $a$ and $b$ is $0.1$
which means that they cannot stop exactly at $p_{a,x}=-5$ or $p_{b,x}=5$.

TABLE II: Closed-loop performance comparison of controllers in the aircraft collision avoidance example. Average values from $5$ randomly sampled initial states are reported. | Safety | Trajectory | Training time | Solve
---|---|---|---|---
|  | cost | per epoch | time
MPC | Safe | 1346.8 | N/A | 44.1s
NN | Unsafe | 712.0 | 0.24s | 0.06s
NN-QP | Safe | 2557.6 | 0.24s | 1.28s
NN-diff-QP | Safe | 1900.7 | 133.0s | 1.30s
NN-gauge | Safe | 2157.8 | 101.5s | 0.74s

All NN controllers are trained similarly as in the ACC example with sampling
rate $\Delta_{t}=0.1$s and horizon $T=2$s, and are tested on 5 randomly
sampled initial states with horizon $T=20$s. The results are shown in Table
II. One set of closed-loop system trajectories starting from the initial state
$x(0)=\begin{bmatrix}0.5&0&\pi&-0.5&0&0\end{bmatrix}^{\top}$ under different
controllers are plotted in Fig. 4 together with the induced costs. With this
initial condition, aircrafts $a$ and $b$ start close to each other with
orientations leading to a head-on collision.

Figure 4: Trajectories of the aircraft $a$ and $b$ under different
controllers. The induced costs of each trajectory are labeled accordingly.

From Table II and Fig. 4, we observe that NN achieves the best performance but
is unsafe. Adding the CBF-QP as a safety filter (i.e., NN-QP) drastically
deteriorates the performance of the NN controller. Among the NN controllers
with safety guarantee, NN-diff-QP performs the best and NN-gauge achieves a
similar level of performance. The MPC controller has the best performance with
safety guarantee, but it has a significantly higher online solve time than NN
controllers.

##  V Conclusion

In this paper, we showed that CBF-based safety filters can deteriorate closed-
loop performance if their long-term effects are not considered during
learning. To address this issue, we proposed a novel safe-by-construction NN
controller which utilizes CBF and gauge map to construct a differentiable
safety layer. The proposed gauge map-based NN controller achieves comparable
performances as the projection-based NN controller while being computationally
more efficient to evaluate online. Both the gauge map-based and projection-
based safe NN controllers demonstrate improved performance with respect to
filtered NN controllers in numerical examples.

## References

  * [1] H. A. Pierson and M. S. Gashler.  Deep learning in robotics: a review of recent research.  Advanced Robotics, 31(16):821–835, 2017.
  * [2] P. L. Donti, M. Roderick, M. Fazlyab, and J. Z. Kolter.  Enforcing robust control guarantees within neural network policies.  In International Conference on Learning Representations, 2020.
  * [3] L. Furieri, C. L. Galimberti, M. Zakwan, and G. Ferrari-Trecate.  Distributed neural network control with dependability guarantees: a compositional port-hamiltonian approach.  In Learning for Dynamics and Control Conference, pages 571–583. PMLR, 2022.
  * [4] 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, 62(8):3861–3876, 2016.
  * [5] A. D. Ames, S. Coogan, M. Egerstedt, G. Notomista, K. Sreenath, and P. Tabuada.  Control barrier functions: Theory and applications.  In 2019 18th European control conference (ECC), pages 3420–3431. IEEE, 2019.
  * [6] A. Anand, K. Seel, V. Gjærum, A. Håkansson, H. Robinson, and A. Saad.  Safe learning for control using control lyapunov functions and control barrier functions: A review.  Procedia Computer Science, 192:3987–3997, 2021.
  * [7] C. Dawson, Z. Qin, S. Gao, and C. Fan.  Safe nonlinear control using robust neural lyapunov-barrier functions.  In Conference on Robot Learning, pages 1724–1735. PMLR, 2022.
  * [8] E. F. Camacho and C. B. Alba.  Model predictive control.  Springer science & business media, 2013.
  * [9] M. H. Cohen and C. Belta.  Approximate optimal control for safety-critical systems with control barrier functions.  In 2020 59th IEEE Conference on Decision and Control (CDC), pages 2062–2067. IEEE, 2020.
  * [10] D. Tabas and B. Zhang.  Computationally efficient safe reinforcement learning for power systems.  arXiv preprint arXiv:2110.10333, 2021.
  * [11] D. Tabas and B. Zhang.  Safe and efficient model predictive control using neural networks: An interior point approach.  arXiv preprint arXiv:2203.12196, 2022.
  * [12] B. Amos and J. Z. Kolter.  Optnet: Differentiable optimization as a layer in neural networks.  In International Conference on Machine Learning, pages 136–145. PMLR, 2017.
  * [13] B. Amos, I. Jimenez, J. Sacks, B. Boots, and J. Z. Kolter.  Differentiable mpc for end-to-end planning and control.  Advances in neural information processing systems, 31, 2018.
  * [14] A. Agrawal, B. Amos, S. Barratt, S. Boyd, S. Diamond, and J. Z. Kolter.  Differentiable convex optimization layers.  Advances in neural information processing systems, 32, 2019.
  * [15] M. A. Pereira, Z. Wang, I. Exarchos, and E. A. Theodorou.  Safe optimal control using stochastic barrier functions and deep forward-backward sdes.  arXiv preprint arXiv:2009.01196, 2020.
  * [16] Y. Emam, P. Glotfelter, Z. Kira, and M. Egerstedt.  Safe model-based reinforcement learning using robust control barrier functions.  arXiv preprint arXiv:2110.05415, 2021.
  * [17] W. S. Cortez, J. Drgona, A. Tuor, M. Halappanavar, and D. Vrabie.  Differentiable predictive control with safety guarantees: A control barrier function approach.  arXiv preprint arXiv:2208.02319, 2022.
  * [18] H. Parwana and D. Panagou.  Recursive feasibility guided optimal parameter adaptation of differential convex optimization policies for safety-critical systems.  In 2022 International Conference on Robotics and Automation (ICRA), pages 6807–6813. IEEE, 2022.
  * [19] K. Xu, W. Xiao, and C. G. Cassandras.  Feasibility guaranteed traffic merging control using control barrier functions.  arXiv preprint arXiv:2203.04348, 2022.
  * [20] H. Ma, B. Zhang, M. Tomizuka, and K. Sreenath.  Learning differentiable safety-critical control using control barrier functions for generalization to novel environments.  arXiv preprint arXiv:2201.01347, 2022.
  * [21] W. Xiao, R. Hasani, X. Li, and D. Rus.  Barriernet: A safety-guaranteed layer for neural networks.  arXiv preprint arXiv:2111.11277, 2021.
  * [22] D. P. Kingma and J. Ba.  Adam: A method for stochastic optimization.  In 3rd International Conference on Learning Representations, ICLR 2015, 2015.
  * [23] A. G. Baydin, B. A. Pearlmutter, A. A. Radul, and J. M. Siskind.  Automatic differentiation in machine learning: a survey.  Journal of Marchine Learning Research, 18:1–43, 2018.
  * [24] J. Drgona, A. Tuor, and D. Vrabie.  Learning constrained adaptive differentiable predictive control policies with guarantees.  arXiv preprint arXiv:2004.11184, 2020.
  * [25] S. Mukherjee, J. Drgoňa, A. Tuor, M. Halappanavar, and D. Vrabie.  Neural lyapunov differentiable predictive control.  arXiv preprint arXiv:2205.10728, 2022.
  * [26] F. Blanchini and S. Miani.  Set-theoretic methods in control, volume 78.  Springer, 2008.
  * [27] H. Wang, K. Margellos, and A. Papachristodoulou.  Explicit solutions for safety problems using control barrier functions.  arXiv preprint arXiv:2204.09380, 2022.
  * [28] S. Boyd and L. Vandenberghe.  Convex optimization.  Cambridge university press, 2004.
  * [29] W. Xiao and C. Belta.  Control barrier functions for systems with high relative degree.  In 2019 IEEE 58th conference on decision and control (CDC), pages 474–479. IEEE, 2019.
  * [30] J. Zeng, B. Zhang, Z. Li, and K. Sreenath.  Safety-critical control using optimal-decay control barrier function with guaranteed point-wise feasibility.  In 2021 American Control Conference (ACC), pages 3856–3863. IEEE, 2021.
  * [31] E. Squires, P. Pierpaoli, and M. Egerstedt.  Constructive barrier certificates with applications to fixed-wing aircraft collision avoidance.  In 2018 IEEE Conference on Control Technology and Applications (CCTA), pages 1656–1661. IEEE, 2018.