# Configurational entropy of polydisperse systems can never reach zero

Vasili Baranau<EMAIL_ADDRESS>Ulrich Tallarek  Department of
Chemistry, Philipps-Universität Marburg, Hans-Meerwein-Strasse 4, 35032
Marburg, Germany

###### Abstract

We present examples of systems whose configurational entropy $S_{\text{conf}}$
can never reach zero and is instead limited from below by the entropy of
mixing $S_{\text{mix}}$ of the corresponding ideal gas. We use
$S_{\text{conf}}$ defined through the local minima of the potential energy
landscape, $S_{\text{conf}}^{\text{PEL}}$. We show that this happens in mean-
field models, in collections of hard spheres with infinitesimal
polydispersity, and for one-dimensional hard rods. We demonstrate that these
results match recent advances in understanding the configurational entropy
defined in the free energy landscape, $S_{\text{conf}}^{\text{FEL}}$. We
demonstrate that if $\min(S_{\text{conf}}^{\text{FEL}})=0$, then for an
arbitrary system $\min(S_{\text{conf}}^{\text{PEL}})=AN+S_{\text{mix}}$, where
$N$ is the number of particles and $A$ is some constant determined by the
interaction potential. We discuss which implications these results have on the
Adam–Gibbs (AG) and RFOT relations and show that the latter retain a
physically meaningful shape for both configurational entropies,
$S_{\text{conf}}^{\text{FEL}}$ and $S_{\text{conf}}^{\text{PEL}}$.

††preprint: AIP/123-QED

##  I Introduction

Recently, Ozawa and Berthier studied a problemOzawa and Berthier (2017) of
diverging and discontinuous configurational entropy per particle
$S_{\text{conf}}/N$ in systems with continuous polydispersity. Divergence
happens due to the diverging and discontinuous (with respect to small
compositional changes) mixing entropy per particle $S_{\text{mix}}/N$.
Specifically, the total entropy $S_{\text{tot}}$ contains the mixing entropy
$S_{\text{mix}}$, while the vibrational entropy $S_{\text{vib}}$ does not.
Thus, the configurational entropy $S_{\text{conf}}$ from the relationSpeedy
(1998a); Stillinger, Debenedetti, and Truskett (2001); Sastry (2001); Angelani
and Foffi (2007); Donev, Stillinger, and Torquato (2007a); Foffi and Angelani
(2008); Starr, Douglas, and Sastry (2013)
$S_{\text{tot}}=S_{\text{conf}}+S_{\text{vib}}$ shall contain
$S_{\text{mix}}$. It is problematic because some observables are believed to
depend on $S_{\text{conf}}$Adam and Gibbs (1965); Cavagna (2009); Starr,
Douglas, and Sastry (2013); Kirkpatrick and Thirumalai (2015); Cavagna (2009);
Starr, Douglas, and Sastry (2013); Bouchaud and Biroli (2004) but from
physical considerations shall be continuous with slight updates of the system
composition, as well as finite and non-zero for continuous polydispersity.
Ozawa and BerthierOzawa and Berthier (2017) suggested a redefinition of
$S_{\text{vib}}$ (a “scheme” of merging the basins of attraction) so that
$S_{\text{mix}}$ is moved into $S_{\text{vib}}$ and $S_{\text{conf}}$ thus
remains finite. The authors refined their approach and provided a more robust
and unambiguous procedure in a follow-up paper.Ozawa, Parisi, and Berthier
(2018) Their approach essentially defines $S_{\text{conf}}$ through the free
energy landscape (FEL).Du and Weeks (2016); Ikeda, Miyazaki, and Ikeda (2016)

In an attempt to resolve this paradox in the potential energy landscape (PEL),
we studiedBaranau and Tallarek (2017) $S_{\text{vib}}^{\text{PEL}}$ using a
variant of thermodynamic integrationFrenkel and Ladd (1984); Frenkel and Smit
(2002); Stillinger, Debenedetti, and Truskett (2001); Donev, Stillinger, and
Torquato (2007b, a); Asenjo, Paillusson, and Frenkel (2014) and suggested that
it is in fact normal for $S_{\text{conf}}^{\text{PEL}}$ to contain
$S_{\text{mix}}$. We suggested that all the observables that depend on
$S_{\text{conf}}^{\text{PEL}}$ and are believed to be well-behaving shall
instead depend on $S_{\text{conf}}^{\text{PEL}}-S_{\text{mix}}$ or in general
$S_{\text{conf}}^{\text{PEL}}-\min(S_{\text{conf}}^{\text{PEL}})$. For this,
we implicitly assumed that $\min(S_{\text{conf}}^{\text{PEL}})$ can be non-
zero, in particular $S_{\text{mix}}$, i.e., that
$S_{\text{conf}}^{\text{PEL}}$ is bounded from below solely by the system
composition.

In this paper, we investigate the lower bound on
$S_{\text{conf}}^{\text{PEL}}$ in more detail and explicitly demonstrate that
for some systems $S_{\text{conf}}^{\text{PEL}}$ is indeed bounded from below
by $S_{\text{mix}}$. Additionally, we investigate connections between
$S_{\text{conf}}^{\text{PEL}}$ and $S_{\text{conf}}$ defined through the free
energy landscape, $S_{\text{conf}}^{\text{FEL}}$. We show that for our model
systems $\min(S_{\text{conf}}^{\text{FEL}})=0$ and that the condition
$\min(S_{\text{conf}}^{\text{FEL}})=0$ implies for an arbitrary system
$\min(S_{\text{conf}}^{\text{PEL}})=NA+S_{\text{mix}}$, where $A$ is some
constant. We discuss which implications it has on the Adam–Gibbs (AG) and RFOT
relations and show that they retain a physically meaningful shape for both
configurational entropies, $S_{\text{conf}}^{\text{FEL}}$ and
$S_{\text{conf}}^{\text{PEL}}$. We show that the two versions of the AG
relation will differ by a well-behaving discrepancy term.

##  II Preliminaries

###  II.1 Different entropies

We study a canonical $NVT$ ensemble of classical systems and use a standard
notation with $N$, $V$, $T$, $U$, $d$, and $\Lambda$ standing for the number
of particles, volume, temperature, internal energy, dimensionality, and de
Broglie wavelength, respectively. $S$, $F$, $Z$, and $U_{N}(\vec{r})$ stand
for the entropy, Helmholtz free energy, partition function, and total
potential energy in the configuration $\vec{r}$. We study polydisperse systems
and denote with $N_{m}$ the number of particles of type $m$.

One can routinely derive from standard equations $S=(U-F)/T$, $F=-k_{B}T\ln
Z$, and
$Z=\frac{1}{\Pi_{i=1}^{M}N_{i}!}\frac{1}{\Lambda^{dN}}\int_{V^{N}}e^{-U_{N}(\vec{r})/k_{B}T}\text{d}\vec{r}$
for the total entropy of the system.Carnahan and Starling (1970); Adams
(1974); Speedy (1998a); Donev, Stillinger, and Torquato (2007b) It can be
shownOzawa and Berthier (2017); Baranau and Tallarek (2017) that such an
entropy will contain the mixing entropy
$S_{\text{mix}}=-k_{B}N\sum_{i}^{M}\frac{N_{i}}{N}\ln\frac{N_{i}}{N}$. It can
also be shown that $S_{\text{mix}}/N$ diverges for continuous particle type
distributions (in the simplest case, $N_{i}=1~{}\forall i$,
$S_{\text{mix}}=k_{B}N\ln N$).Ozawa and Berthier (2017); Baranau and Tallarek
(2017) $S_{\text{mix}}/N$ is also obviously discontinuous if the system
composition is abruptly changed.

If we define the configurational entropy through the PEL,
$S_{\text{conf}}^{\text{PEL}}$ counts (up to particle permutations) the
density of potential energy minimaDonev, Stillinger, and Torquato (2007a);
Ashwin  _et al._ (2012); Asenjo, Paillusson, and Frenkel (2014); Martiniani
_et al._ (2016); Stillinger, DiMarzio, and Kornegay (1964); Stillinger (1995);
Torquato and Jiao (2010) $N_{\text{min}}(N,u_{\text{min}})$ at a given energy
per particle $u_{\text{min}}=U_{\text{min}}/N$. Specifically,
$N_{\text{min}}(N,u_{\text{min}})=e^{S_{\text{conf}}^{\text{PEL}}(N,u_{\text{min}})}/w$,
where $w$ is the characteristic interval that comes from the saddle point
approximation.Speedy (1998a); Parisi and Zamponi (2005, 2010); Berthier,
Jacquin, and Zamponi (2011) It is believed to be subexponential, i.e.,
$\ln(w)/N\to 0$ with $N\to\infty$.

The vibrational partition function in the PEL is expressed as
$Z_{\text{vib}}^{\text{PEL}}=\frac{\Pi_{i=1}^{M}N_{i}!}{\Pi_{i=1}^{M}N_{i}!}\frac{1}{\Lambda^{dN}}\int_{\text{basin}}e^{-U_{N}(\vec{r})/k_{B}T}\text{d}\vec{r}$,
where integration is performed over a certain basin of attractionStillinger,
DiMarzio, and Kornegay (1964); Stillinger (1995); Debenedetti and Stillinger
(2001); Torquato and Jiao (2010) in the PEL. The $\Pi_{i=1}^{M}N_{i}!$ term in
the denominator comes from particle indistinguishability as for the total
partition function.Frenkel (2014); Meng  _et al._ (2010); Asenjo, Paillusson,
and Frenkel (2014); Martiniani  _et al._ (2016) The same term in the nominator
comes from basin multiplicity, i.e., indistinguishability of basins due to
particle indistinguishability.Stillinger and Salsburg (1969); Speedy (1993,
1998b); Ozawa and Berthier (2017) Due to cancellation of these two terms,
$S_{\text{vib}}^{\text{PEL}}$ does not contain $S_{\text{mix}}$. We
showedBaranau and Tallarek (2017) on the basis of thermodynamic integration
that $S_{\text{mix}}$ is also not hidden in the integral over a basin and that
$S_{\text{vib}}^{\text{PEL}}$ still contains the term
$\ln\left(\frac{V}{\Lambda^{d}N}\right)$, which also enters $S_{\text{tot}}$.

The three entropies are connected by the relation
$S_{\text{tot}}=S_{\text{conf}}^{\text{PEL}}+S_{\text{vib}}^{\text{PEL}}$,Speedy
(1998a); Stillinger, Debenedetti, and Truskett (2001); Sastry (2001); Angelani
and Foffi (2007); Donev, Stillinger, and Torquato (2007a); Foffi and Angelani
(2008); Starr, Douglas, and Sastry (2013) which comes from the saddle point
approximation and is exact in the thermodynamic limit if the definitions from
above are used.Speedy (1998a); Baranau and Tallarek (2017) In this relation,
$S_{\text{conf}}^{\text{PEL}}$ is the equilibrium configurational entropy
(equilibrium complexity), i.e., computed for the minima of those basins that
are equilibrium (dominant) for the given temperature. At a given temperature
$T$ basins with a certain average energy of minima $u_{\text{eq,min}}(T)$
dominate the phase space. Then, the relation for entropies shall fully be read
as

$\displaystyle S_{\text{tot}}(N,V,T)=$ $\displaystyle
S_{\text{conf}}^{\text{PEL}}(N,V,u_{\text{eq,min}}(T))$ (1)
$\displaystyle+S_{\text{vib}}^{\text{PEL}}(N,V,u_{\text{eq,min}}(T),T).$

Because $S_{\text{tot}}$ contains $S_{\text{mix}}$ and
$S_{\text{vib}}^{\text{PEL}}$ does not, $S_{\text{mix}}$ shall be contained in
$S_{\text{conf}}^{\text{PEL}}$. If the configurational and vibrational
entropies are defined through the FEL, the equation also holds.

The configurational entropy per particle $S_{\text{conf}}/N$ (however it is
defined) enters several relations for the relaxation time $\tau$ of colloids
and glassy systems at high volume fraction or low temperature. Both versions
of $S_{\text{conf}}$ are used by different authors: either
$S_{\text{conf}}^{\text{PEL}}$ or $S_{\text{conf}}^{\text{FEL}}$. The well-
known Adam–Gibbs relation reads as
$\tau=\tau_{0}\exp\left(\frac{C}{S_{\text{conf}}/N}\right)$, where $\tau_{0}$
and $C$ are presumed to be approximately constant.Adam and Gibbs (1965);
Cavagna (2009); Starr, Douglas, and Sastry (2013) The other relation stems
from the Random First Order Theory (RFOT).Kirkpatrick and Thirumalai (2015);
Cavagna (2009); Starr, Douglas, and Sastry (2013); Bouchaud and Biroli (2004)
If $S_{\text{conf}}/N$ is divergent for continuous particle type distributions
and discontinuous with a small system composition change, $\tau\equiv 0$ for
systems with continuous particle type distributions and $\tau$ will be
discontinuous when the system composition is slightly changed. It is
paradoxical because $\tau$ shall be finite and continuous.

Even worse, particle type distribution can be induced “artificially”, by
assigning arbitrary properties to particles without changing their interaction
potential.Frenkel (2014) For example, colloidal particles can be colored
(dyed) arbitrarily. This will increase $S_{\text{conf}}^{\text{PEL}}$
arbitrarily without changing $\tau$.

#### Hard spheres

A valid collection of $N$ frictionless hard spheres with predefined radii
$R_{i}$ always has a zero potential energy, while an invalid collection has an
infinite energy. Still, a pseudo-PEL can be introduced for hard
spheres.Torquato and Jiao (2010); Baranau and Tallarek (2014); Zinchenko
(1994) At any system configuration $\vec{x}\in\mathbb{R}^{3N}$ particle radii
can be scaled proportionally to ensure that at least one pair of particles is
in contact. The resulting solid volume fraction $\varphi$ defines the pseudo-
PEL $U(\vec{x})=-\varphi$ (or rather, $U(\vec{x})=-N\varphi$ for $U(\vec{x})$
to be extensive). Local minima in this PEL correspond to mechanically stable
(jammed) configurations with some jamming densities $\varphi_{J}$, while
steepest descents correspondTorquato and Jiao (2010); Baranau and Tallarek
(2014); Zinchenko (1994) to proportionally scaling sphere radii and moving the
spheres as little as possible in the configuration space to ensure the absence
of intersections. Eq. (1) will then look like
$S_{\text{tot}}(N,\varphi,T)=S_{\text{conf}}^{\text{PEL}}(N,\varphi_{\text{eq,J}}(\varphi))+S_{\text{vib}}^{\text{PEL}}(N,\varphi_{\text{eq,J}}(\varphi),T)$.

###  II.2 Motivation

We arguedBaranau and Tallarek (2017) that it is normal for
$S_{\text{conf}}^{\text{PEL}}$ to be divergent and discontinuous. Still, the
relaxation time shall be finite for continuous particle type distributions and
shall be continuous when slightly changing particle types or coloring the
particles. We suggested that the equations for $\tau$ and for any other “well-
behaving” observables shall be changed instead. If such an observable depends
on $S_{\text{conf}}^{\text{PEL}}$, it shall always depend on
$S_{\text{conf}}^{\text{PEL}}-S_{\text{mix}}$ and can not depend on
$S_{\text{conf}}^{\text{PEL}}$ only. Thus, arbitrary changes in
$S_{\text{conf}}^{\text{PEL}}$ will be compensated.

Indeed, if $S_{\text{mix}}$ is the lower bound for
$S_{\text{conf}}^{\text{PEL}}$, then only changes with respect to
$S_{\text{mix}}$ can matter for physical observables. A more detailed
derivation for the updated Adam–Gibbs equation can be found in our previous
paper.Baranau and Tallarek (2017) For this derivation, we implicitly assumed
that $S_{\text{conf}}^{\text{PEL}}$ can not decrease below $S_{\text{mix}}$.
In this paper, we provide a more detailed explanation and physical insight
into this lower bound for $S_{\text{conf}}^{\text{PEL}}$.

We do not discuss in this paper whether $S_{\text{conf}}^{\text{PEL}}$ is a
better candidate for such equations than $S_{\text{conf}}^{\text{FEL}}$. But
we point out that there indeed is a body of work showing that
$S_{\text{conf}}^{\text{PEL}}$ can be a possible candidate for
$S_{\text{conf}}$ in the Adam–Gibbs and RFOT relations.Speedy (1998a);
Angelani and Foffi (2007); Starr, Douglas, and Sastry (2013)

##  III Main part

In the worst case, the lowest equilibrium complexity
$S_{\text{conf}}^{\text{PEL}}(N,V,u_{\text{eq,min}}(T))$ at some low
temperature $T_{0}$ will correspond to only one configuration of particles (up
to particle permutations). This configuration will be the local minimum of a
basin of attraction that is equilibrium at this $T_{0}$ and dominates the
phase space. If particles are considered distinguishable and there are $N_{m}$
particles of each type $m$, there are $\Pi_{m}N_{m}!$ such local minima.
Still, there are $N!$ particle permutations in total. Each of $N!$
configurations (if particles are considered distinguishable) lies in a certain
basin of attraction in the potential energy landscape and each of these basins
has a corresponding energy minimum (inherent structure). The question is
whether these basins will also be equilibrium at $T_{0}$ or not.

###  III.1 Mean-field limit

To determine if a certain configuration among these $N!$ is in equilibrium at
$T_{0}$, we compare its chemical potential $\mu$ with the chemical potential
from the original equilibrium configuration. To compute the chemical
potential, we rely on the Widom particle insertion method.Widom (1963); Adams
(1974); Speedy (1981); Speedy and Reiss (1991); Frenkel and Smit (2002);
Baranau and Tallarek (2016) This method relates the chemical potential $\mu$
to the potential $\Psi_{N+1,i}(\vec{r})$ that is experienced by a $(N+1)$th
test particle of type $i$ if we try to insert this particle into the system at
a position $\vec{r}$, sampled uniformly. Specifically,

$-\frac{\Delta\mu_{i}}{kT}=\ln\left[\left\langle\exp\left(-\frac{\Psi_{N+1,i}}{kT}\right)\right\rangle\right],$
(2)

where $\Delta\mu_{i}$ is the excess chemical potential for particles of type
$i$ and averaging is performed with respect to $\vec{r}$ of particle insertion
positions.Baranau and Tallarek (2016) If the original configuration is in
equilibrium, then we estimate the chemical potential (the equilibrium one).

Thus, our main focus is how the field $\Psi_{N+1,i}$ changes when we perform
any of the $N!$ particle permutations in the initial equilibrium
configuration. We do not have an answer to this question in the general case,
but this question can be resolved in the mean-field limit.

We will denote the potential $\Psi_{N+1,i}(\vec{r})$ perturbed by a
permutation $\Sigma_{p}$ as $\Psi_{\Sigma_{p},N+1,i}$. In the mean-field
limit, the potential $\Psi_{\Sigma_{p},N+1,i}$ does not depend on
$\Sigma_{p}$,

$\Psi_{\Sigma_{p},N+1,i}\equiv\Psi_{N+1,i}.$ (3)

Thus, all of the $N!$ configurations are equilibrium. The number of such
configurations up to particle permutations is $\frac{N!}{\Pi_{i}N_{i}!}$ and
the lowest possible equilibrium complexity is thus

$\min(S_{\text{conf}}^{\text{PEL}})=k_{B}\ln\left(\frac{N!}{\Pi_{i}N_{i}!}\right).$
(4)

It can be easily shown with the help of the Stirling approximation $N!\sim
N\ln N-N$ that $\min(S_{\text{conf}}^{\text{PEL}})$ is

$\min(S_{\text{conf}}^{\text{PEL}})=-k_{B}N\sum_{i}\frac{N_{i}}{N}\ln\frac{N_{i}}{N}=S_{\text{mix}}.$
(5)

###  III.2 Infinitesimally polydisperse systems

Let us assume that a system is monodisperse and the worst case is realized for
it, i.e., there is only one equilibrium basin (up to particle permutations) at
the lowest temperature $T_{0}$ where the system can still be equilibrated.
Thus, there are $N!$ equilibrium basins and respective minima if particles are
considered distinguishable. Let us now assume that particle types are slightly
changed, e.g., particle charges that govern pair potentials or radii of
colloidal particles are infinitesimally changed.

For sufficiently small particle type changes, the structure of the phase space
and its basins will remain unchanged. Thus, all of the $N!$ basins will still
be equilibrium at this $T_{0}$— but there are now $\Pi_{m}N_{m}!$ particle
permutations that keep the system in the same state, instead of $N!$. Hence,
the minimum configurational entropy becomes
$\min(S_{\text{conf}}^{\text{PEL}})=k_{B}\ln\left(\frac{N!}{\Pi_{i}N_{i}!}\right)=S_{\text{mix}}$.

###  III.3 $1D$ hard rods

Our last example is very brief and almost obvious. A collection of $1D$ hard
rods has $N!$ jammed configurations (pseudo-PEL minima) if particles are
considered distinguishable. Because there are $\Pi_{m}N_{m}!$ particle
permutations if particles are considered indistinguishable,
$S_{\text{conf}}^{\text{PEL}}=k_{B}\ln\left(\frac{N!}{\Pi_{i}N_{i}!}\right)=S_{\text{mix}}$.
Particle types can stem either from actual rod lengths or from artificial
labeling assigned by the observer arbitrarily. This is actually the only value
of $S_{\text{conf}}^{\text{PEL}}$ available for $1D$ hard rods.

###  III.4 Remark: colloids and hard spheres

The procedure outlined for the mean-field case looks as follows for hard
spheres. Assume that the largest density at which the system can be
equilibrated is denoted as $\varphi_{K}$. The corresponding jamming density
can be denoted as $\varphi_{\text{GCP}}$, the Glass Close Packing
density.Parisi and Zamponi (2010) Then, one has to take an equilibrium
configuration at $\varphi=\varphi_{K}$ (maybe even the one with
$\varphi=\varphi_{\text{GCP}}$ but with particle radii proportionally
decreased to achieve $\varphi_{K}$) and performs particle permutations. Some
pairs of particles may intersect after permutations. To avoid permutations,
one decreases particle radii proportionally until the closest pair of
particles just touches each other; the density is decreased in this process
($\varphi<\varphi_{K}$). Then, one increases particle radii proportionally and
simultaneously moves particles to ensure minimal possible movement in the
configuration space (moves particles along contact normals). This movement
corresponds to the steepest descent in the pseudo-PEL. One performs this
particle scaling until $\varphi=\varphi_{K}$. Then, one has to check if the
resulting configuration is equilibrium. Numerically, this can be done by
employing the Widom particle insertion method. If the system becomes jammed
before reaching $\varphi_{K}$, then the new basin (after permutation) was
definitely not in equilibrium.

##  IV Relationship with $S_{\text{conf}}$ defined through the free energy
landscape

We would like to discuss how $S_{\text{conf}}^{\text{PEL}}$ relates to the
configurational entropy defined in the free energy landscape,
$S_{\text{conf}}^{\text{FEL}}$. The free energy landscape to compute
$S_{\text{conf}}^{\text{FEL}}$ can stem, e.g., from the density functional
theory. We will show that for our model systems the result
$\min(S_{\text{conf}}^{\text{PEL}})=S_{\text{mix}}$ complies with the
requirement $\min(S_{\text{conf}}^{\text{FEL}})=0$ for a certain definition of
$S_{\text{conf}}^{\text{FEL}}$.Ozawa, Parisi, and Berthier (2018)

###  IV.1 $S_{\text{conf}}^{\text{FEL}}$ of Ozawa et al.Ozawa, Parisi, and
Berthier (2018)

A recent advancement in understanding $S_{\text{conf}}^{\text{FEL}}$ can be
found in a paper by Ozawa, Parisi, and Berthier.Ozawa, Parisi, and Berthier
(2018) Their approach is explicit and mathematically precise and can be used
in computer simulations. Essentially, while $S_{\text{vib}}^{\text{PEL}}$ is
computed over a certain PEL basin, the authors compute
$S_{\text{vib}}^{\text{FEL}}$ by including in the vibrational partition
function all possible particle permutations around some reference
configurations and weighting the permutations according to their statistical
weights. Reference configurations are sampled from the entire configuration
space, but are also taken with corresponding statistical weights. We will rely
on this approach for the discussion below.

We briefly repeat the essential equations from that paper (note that the paper
uses $k_{B}=1$, so we follow this convention in this section as well). The
authors use the modified Frenkel–Ladd method of computing $S_{\text{vib}}$;
i.e., they attach particles with springs to the reference positions and
integrate over the spring constant $\alpha$. The authors utilize the following
types of statistical averages:

$\displaystyle\left\langle(\cdots)\right\rangle_{\alpha}^{\rm T,S}$
$\displaystyle=$
$\displaystyle\frac{\frac{1}{N!}\sum_{\pi}\int_{V}\mathrm{d}{\bf
r}^{N}(\cdots)e^{-\beta U_{\alpha}(\Sigma_{\pi}^{N},{\bf r}^{N},{\bf
r}_{0}^{N})}}{\frac{1}{N!}\sum_{\pi}\int_{V}\mathrm{d}{\bf r}^{N}e^{-\beta
U_{\alpha}(\Sigma_{\pi}^{N},{\bf r}^{N},{\bf r}_{0}^{N})}},$ (6)
$\displaystyle\left\langle(\cdots)\right\rangle_{\alpha}^{\rm T}$
$\displaystyle=$ $\displaystyle\frac{\int_{V}\mathrm{d}{\bf
r}^{N}(\cdots)e^{-\beta U_{\alpha}({\bf r}^{N},{\bf
r}_{0}^{N})}}{\int_{V}\mathrm{d}{\bf r}^{N}e^{-\beta U_{\alpha}({\bf
r}^{N},{\bf r}_{0}^{N})}},$ (7) $\displaystyle\overline{(\cdots)}$
$\displaystyle=$ $\displaystyle\frac{\int_{V}\mathrm{d}{\bf
r}_{0}^{N}(\cdots)e^{-\beta U_{0}({\bf r}_{0}^{N})}}{\int_{V}\mathrm{d}{\bf
r}_{0}^{N}e^{-\beta U_{0}({\bf r}_{0}^{N})}}.$ (8)

The superscripts T and S represent the statistical average over positions (T)
and permutations (S), respectively. Here, $U_{\alpha}(\Sigma_{\pi}^{N},{\bf
r}^{N},{\bf r}_{0}^{N})$ is the potential energy at a configuration ${\bf
r}^{N}$ around a reference configuration ${\bf r}_{0}^{N}$ when particles are
additionally permuted by a permutation $\pi$ (note that $\Sigma_{\pi}$, not
$\Sigma_{\pi}^{N}$, represents a usual sum over permutations).

The glassy entropy (in our terms $S_{\text{vib}}^{\text{FEL}}$) is expressed
according to their approach as

$\displaystyle S_{\text{vib}}^{\text{FEL}}$ $\displaystyle=$ $\displaystyle
S_{\rm
glass}=\frac{Nd}{2}-N\ln\Lambda^{d}-\frac{Nd}{2}\ln\left(\frac{\alpha_{\rm
max}}{\pi}\right)$ $\displaystyle+N\lim_{\alpha_{\rm min}\to
0}\int_{\alpha_{\rm min}}^{\alpha_{\rm
max}}\mathrm{d}\alpha\Delta_{\alpha}^{\rm T,S}+S_{\rm
mix}-\overline{\mathcal{S}_{\rm mix}({\bf r}_{0}^{N},\beta)},$

where $\Delta_{\alpha}^{\rm T,S}$ is a mean-squared displacement defined by

$\Delta_{\alpha}^{\rm
T,S}=\frac{1}{N}\overline{\left\langle\sum_{i=1}^{N}|{\bf r}_{i}-{\bf
r}_{0i}|^{2}\right\rangle_{\alpha}^{\rm T,S}},$ (10)

and $\mathcal{S}_{\rm mix}({\bf r}_{0}^{N},\beta)$ is a mixing entropy
contribution defined by

$\mathcal{S}_{\rm mix}({\bf
r}_{0}^{N},\beta)=-\ln\left(\frac{1}{N!}\sum_{\pi}e^{-\beta\left(U_{0}(\Sigma_{\pi}^{N},{\bf
r}_{0}^{N})-U_{0}({\bf r}_{0}^{N})\right)}\right).$ (11)

$S_{\text{conf}}^{\text{FEL}}$ is then given as usual by
$S_{\text{tot}}-S_{\text{vib}}^{\text{FEL}}$.

Compare Eq. (LABEL:eq:S_glass_final) to a “usual” expression for
$S_{\text{vib}}^{\text{PEL}}$, when particle permutations are not included in
the definition:

$\displaystyle S_{\text{vib}}^{\text{PEL}}$ $\displaystyle=$
$\displaystyle\frac{Nd}{2}-N\ln\Lambda^{d}-\frac{Nd}{2}\ln\left(\frac{\alpha_{\rm
max}}{\pi}\right)$ (12) $\displaystyle+N\lim_{\alpha_{\rm min}\to
0}\int_{\alpha_{\rm min}}^{\alpha_{\rm
max}}\mathrm{d}\alpha\Delta_{\alpha}^{\rm T}.$

Here,

$\Delta_{\alpha}^{\rm T}=\frac{1}{N}\overline{\left\langle\sum_{i=1}^{N}|{\bf
r}_{i}-{\bf r}_{0i}|^{2}\right\rangle_{\alpha}^{\rm T}}.$ (13)

We note that a typical definition for $S_{\text{vib}}^{\text{PEL}}$ would use
$\Delta_{\alpha}^{\rm T}=\frac{1}{N}\left\langle\sum_{i=1}^{N}|{\bf
r}_{i}-{\bf r}_{0i}|^{2}\right\rangle_{\alpha}^{\rm T}$, where ${\bf r}_{0}$
represents a PEL minimum that is equilibrium at a given temperature. It is an
approximation for averaging over all equilibrium PEL minima, which is itself
an approximation for averaging over all PEL minima, given that they are taken
with corresponding statistical weights. We note that integrals in Eqs. (8) and
(13) will be dominated by (equilibrium) PEL minima, so using Eq. (13) for
$S_{\text{vib}}^{\text{PEL}}$ is legitimate.

###  IV.2 Model systems from our paper

For our examples (a mean-field system and an infinitesimally polydisperse
system), particle permutations do not change the potential energy. Hence,
$U_{\alpha}(\Sigma_{\pi}^{N},{\bf r}^{N},{\bf r}_{0}^{N})=U_{\alpha}({\bf
r}^{N},{\bf r}_{0}^{N})$. It then follows that $\mathcal{S}_{\rm mix}({\bf
r}_{0}^{N},\beta)=0$ (cf. Eq. (11)) and $\Delta_{\alpha}^{\rm
T,S}=\Delta_{\alpha}^{\rm T}$ (cf. Eqs. (6) and (7)). Hence, for these systems
$S_{\text{vib}}^{\text{FEL}}=S_{\text{vib}}^{\text{PEL}}+S_{\text{mix}}$ (cf.
Eqs. (LABEL:eq:S_glass_final) and (12)) and finally
$S_{\text{conf}}^{\text{FEL}}=S_{\text{conf}}^{\text{PEL}}-S_{\text{mix}}$.
Thus, our bound $\min(S_{\text{conf}}^{\text{PEL}})=S_{\text{mix}}$ is
actually equivalent to the usual assumption
$\min(S_{\text{conf}}^{\text{FEL}})=0$:

$\min(S_{\text{conf}}^{\text{FEL}})=0\iff\min(S_{\text{conf}}^{\text{PEL}})=S_{\text{mix}}.$
(14)

###  IV.3 $\min(S_{\text{conf}}^{\text{PEL}})$ in the general case

By comparing Eqs. (LABEL:eq:S_glass_final) and (12) and utilizing
$S_{\text{tot}}=S_{\text{conf}}+S_{\text{vib}}$, we get in the general case
$S_{\text{conf}}^{\text{PEL}}=S_{\text{conf}}^{\text{FEL}}+S_{\text{vib}}^{\text{FEL}}-S_{\text{vib}}^{\text{PEL}}$
and finally

$\displaystyle S_{\text{conf}}^{\text{PEL}}=$ $\displaystyle
S_{\text{conf}}^{\text{FEL}}+N\lim_{\alpha_{\rm min}\to 0}\int_{\alpha_{\rm
min}}^{\alpha_{\rm max}}\mathrm{d}\alpha\left(\Delta_{\alpha}^{\rm
T,S}-\Delta_{\alpha}^{\rm T}\right)$ (15) $\displaystyle+S_{\rm
mix}-\overline{\mathcal{S}_{\rm mix}({\bf r}_{0}^{N},\beta)}.$

One can presume that $\lim_{\alpha_{\rm min}\to 0}\int_{\alpha_{\rm
min}}^{\alpha_{\rm max}}\mathrm{d}\alpha\left(\Delta_{\alpha}^{\rm
T,S}-\Delta_{\alpha}^{\rm T}\right)$ and $\overline{\mathcal{S}_{\rm mix}({\bf
r}_{0}^{N},\beta)}$ are decreasing functions of temperature. Indeed, the lower
the temperature, the less statistical weight permuted configurations have and
hence the closer $\Delta_{\alpha}^{\rm T,S}$ is to $\Delta_{\alpha}^{\rm T}$.
The same applies to $\min(\overline{\mathcal{S}_{\rm mix}({\bf
r}_{0}^{N},\beta)})$. It is usually believed that
$S_{\text{conf}}^{\text{FEL}}$ reaches the minimum value when the phase space
becomes non-ergodic (at $T=T_{K}$). It means that
$S_{\text{conf}}^{\text{PEL}}$ is defined in the same range of temperatures
and is not defined for $T<T_{K}$. Because Eq. (15) can be represented as a sum
of terms that are either constant or decreasing with $T$ and because
$S_{\text{conf}}^{\text{PEL}}$ is defined in the same range of temperatures as
$S_{\text{conf}}^{\text{FEL}}$, we conclude that
$\min(S_{\text{conf}}^{\text{PEL}})$ happens at the same temperature as
$\min(S_{\text{conf}}^{\text{FEL}})$. Hence, we write

$\displaystyle\min$
$\displaystyle(S_{\text{conf}}^{\text{PEL}})=\min(S_{\text{conf}}^{\text{FEL}})$
(16) $\displaystyle+N\min(\lim_{\alpha_{\rm min}\to 0}\int_{\alpha_{\rm
min}}^{\alpha_{\rm max}}\mathrm{d}\alpha\left(\Delta_{\alpha}^{\rm
T,S}-\Delta_{\alpha}^{\rm T}\right))$ $\displaystyle+S_{\rm
mix}-\min(\overline{\mathcal{S}_{\rm mix}({\bf r}_{0}^{N},\beta)})$
$\displaystyle=NA+S_{\rm mix},$

where $A$ is some constant determined purely by particle interactions.

A crucial feature in Eq. (16) is that $S_{\rm mix}$ can be increased
arbitrarily by infinitesimal changes in particle types or by artificial
“coloring” of particles (assigning types that do not influence interaction
potentials). If the distribution of “colors” is continuous, then $S_{\rm
mix}/N\to\infty$. Thus, if we presume that
$\min(S_{\text{conf}}^{\text{FEL}})=0$, then
$\min(S_{\text{conf}}^{\text{PEL}})$ can be increased arbitrarily for a
general system just by additional labeling of particles:

$\min(S_{\text{conf}}^{\text{FEL}})=0\implies\min(S_{\text{conf}}^{\text{PEL}})\text{
is arbitrary}.$ (17)

Additionally, $\min(S_{\text{conf}}^{\text{PEL}})/N$ in this case is
discontinuous with introducing a slight polydispersity to a monodisperse
system and divergent for a continuous polydispersity.

###  IV.4 Adam–Gibbs and RFOT relations

As we argued in our previous paper,Baranau and Tallarek (2017) arbitrary value
of $\min(S_{\text{conf}}^{\text{PEL}})$ and other problems with it are
actually not problematic for observables that may depend on
$S_{\text{conf}}^{\text{PEL}}$, like the relaxation time $\tau_{R}$ of glassy
systems. If one re-derives the Adam–Gibbs and RFOT relations between
$\tau_{R}$ and $S_{\text{conf}}$ assuming that $\min(S_{\text{conf}})\neq 0$,
one obtains

$\tau_{R}=f\left(\frac{S_{\text{conf}}-\min(S_{\text{conf}})}{N}\right)$ (18)

instead of the usual version $\tau_{R}=f(S_{\text{conf}}/N)$. For example, the
Adam–Gibbs relation shall be transformed from a standard version
$\tau_{R}^{\text{AG}}=\tau_{0}\exp\left(\frac{C}{TS_{\text{conf}}/N}\right)$
into

$\tau_{R}^{\text{AG}}=\tau_{0}\exp\left(\frac{C}{T\left[S_{\text{conf}}-\min(S_{\text{conf}})\right]/N}\right),$
(19)

where $\tau_{0}$ and $C$ are (approximately) constants.

It is interesting to compare what will be the difference between
$\tau_{R}^{\text{PEL}}$ and $\tau_{R}^{\text{FEL}}$, given that
$\tau_{R}^{\text{FEL}}=f(S_{\text{conf}}^{\text{FEL}}/N)$ and
$\tau_{R}^{\text{PEL}}=f([S_{\text{conf}}^{\text{PEL}}-\min(S_{\text{conf}}^{\text{PEL}})]/N)$.
By using Eqs. (15) and (16) and again presuming
$\min(S_{\text{conf}}^{\text{FEL}})=0$, we get for the Adam–Gibbs equation

$\displaystyle\tau_{R}^{\text{FEL}}=$
$\displaystyle\tau_{0}\exp\left(\frac{C}{TS_{\text{conf}}^{\text{FEL}}/N}\right),$
(20) $\displaystyle\tau_{R}^{\text{PEL}}=$
$\displaystyle\tau_{0}\exp\left(\frac{C}{T[S_{\text{conf}}^{\text{FEL}}/N+\delta]}\right),$

where $\delta$ is expressed as

$\displaystyle\delta=$ $\displaystyle\lim_{\alpha_{\rm min}\to
0}\int_{\alpha_{\rm min}}^{\alpha_{\rm
max}}\mathrm{d}\alpha\left(\Delta_{\alpha}^{\rm T,S}-\Delta_{\alpha}^{\rm
T}\right)$ (21) $\displaystyle-\min(\lim_{\alpha_{\rm min}\to
0}\int_{\alpha_{\rm min}}^{\alpha_{\rm
max}}\mathrm{d}\alpha\left(\Delta_{\alpha}^{\rm T,S}-\Delta_{\alpha}^{\rm
T}\right))$ $\displaystyle-\left(\overline{\mathcal{S}_{\rm mix}({\bf
r}_{0}^{N},\beta)}-\min(\overline{\mathcal{S}_{\rm mix}({\bf
r}_{0}^{N},\beta)})\right)/N.$

As can be seen from Eq. (21), $\delta$ possesses nice properties: (i) it is a
quantity per particle and (ii) the mixing entropy per particle
$S_{\text{mix}}/N$ is not present in $\delta$ (it was canceled out). Thus,
$\delta$ will remain continuous with introduction of a small polydispersity,
finite for a continuous polydispersity, and independent of artificial particle
labeling (defined only by interaction potentials). Thus, both forms of
$\tau_{R}^{\text{AG}}$ from Eq. (20) seem feasible. Which form better fits
actual data can be tested through simulations. We plan to present in a follow-
up paper results that demonstrate that $\tau_{R}^{\text{PEL}}$ can fit
relaxation times for polydisperse hard spheres well.

##  V Conclusion

We presented several examples of systems where the lower bound of the
configurational entropy $S_{\text{conf}}^{\text{PEL}}$ is the mixing entropy
$S_{\text{mix}}$. Our examples are mean-field systems, infinitesimally
polydisperse systems, and $1D$ hard rods. We believe these examples provide a
deep physical insight into why $S_{\text{conf}}^{\text{PEL}}$ shall in
principle have a non-zero lower bound (which still does not preclude the
existence of the ideal glass transition). Additionally, we demonstrated that
this result is equivalent for our model systems to the condition
$\min(S_{\text{conf}}^{\text{FEL}})=0$. Finally, we showed that the condition
$\min(S_{\text{conf}}^{\text{FEL}})=0$ implies for an arbitrary system that
$\min(S_{\text{conf}}^{\text{PEL}})=NA+S_{\rm mix}$, where $A$ is a constant
determined by particle interaction potentials, while $S_{\rm mix}$ can be
increased arbitrarily by introducing “artificial” particle types (e.g.,
coloring colloidal particles). Thus, $\min(S_{\text{conf}}^{\text{PEL}})$ can
be changed arbitrarily. If the distribution of “colors” is continuous,
$\min(S_{\text{conf}}^{\text{PEL}})/N\to\infty$. Also,
$\min(S_{\text{conf}}^{\text{PEL}})/N$ will be diverging by introducing slight
polydispersity to particle types.

Nevertheless, we argued that the AG and RFOT relations will still retain a
physically meaningful shape. We compared the AG relations derived through both
versions of $S_{\text{conf}}$, $S_{\text{conf}}^{\text{FEL}}$ and
$S_{\text{conf}}^{\text{PEL}}$ and demonstrated that they will be difference,
but the discrepancy term will be well-behaving. We suggested that one can
check which version better describes relaxation times through simulations. We
plan to present in a follow-up paper results that demonstrate that
$\tau_{R}^{\text{PEL}}$ can fit relaxation times for polydisperse hard spheres
well.

## Acknowledgments

We thank Misaki Ozawa and Ludovic Berthier for helpful correspondence about
the present manuscript (though they do not necessarily agree with all the
ideas presented here).

## References

  * Ozawa and Berthier (2017) M. Ozawa and L. Berthier, J. Chem. Phys. 146, 014502 (2017).
  * Speedy (1998a) R. J. Speedy, Mol. Phys. 95, 169 (1998a).
  * Stillinger, Debenedetti, and Truskett (2001) F. H. Stillinger, P. G. Debenedetti, and T. M. Truskett, J. Phys. Chem. B 105, 11809 (2001).
  * Sastry (2001) S. Sastry, Nature 409, 164 (2001).
  * Angelani and Foffi (2007) L. Angelani and G. Foffi, J. Phys.: Condens. Matter 19, 256207 (2007).
  * Donev, Stillinger, and Torquato (2007a) A. Donev, F. H. Stillinger, and S. Torquato, J. Chem. Phys. 127, 124509 (2007a).
  * Foffi and Angelani (2008) G. Foffi and L. Angelani, J. Phys.: Condens. Matter 20, 075108 (2008).
  * Starr, Douglas, and Sastry (2013) F. W. Starr, J. F. Douglas, and S. Sastry, J. Chem. Phys. 138, 12A541 (2013).
  * Adam and Gibbs (1965) G. Adam and J. H. Gibbs, J. Chem. Phys. 43, 139 (1965).
  * Cavagna (2009) A. Cavagna, Phys. Rep. 476, 51 (2009).
  * Kirkpatrick and Thirumalai (2015) T. Kirkpatrick and D. Thirumalai, Rev. Mod. Phys. 87, 183 (2015).
  * Bouchaud and Biroli (2004) J.-P. Bouchaud and G. Biroli, J. Chem. Phys. 121, 7347 (2004).
  * Ozawa, Parisi, and Berthier (2018) M. Ozawa, G. Parisi, and L. Berthier, arXiv:1805.06017 [cond-mat] (2018), arXiv: 1805.06017.
  * Du and Weeks (2016) X. Du and E. R. Weeks, Phys. Rev. E 93, 062613 (2016).
  * Ikeda, Miyazaki, and Ikeda (2016) H. Ikeda, K. Miyazaki, and A. Ikeda, J. Chem. Phys. 145, 216101 (2016).
  * Baranau and Tallarek (2017) V. Baranau and U. Tallarek, J. Chem. Phys. 147, 224503 (2017).
  * Frenkel and Ladd (1984) D. Frenkel and A. J. C. Ladd, J. Chem. Phys. 81, 3188 (1984).
  * Frenkel and Smit (2002) D. Frenkel and B. Smit,  _Understanding Molecular Simulation: From Algorithms to Applications_ , 2nd ed. (Academic Press, San Diego, 2002).
  * Donev, Stillinger, and Torquato (2007b) A. Donev, F. H. Stillinger, and S. Torquato, J. Comput. Phys. 225, 509 (2007b).
  * Asenjo, Paillusson, and Frenkel (2014) D. Asenjo, F. Paillusson, and D. Frenkel, Phys. Rev. Lett. 112, 098002 (2014).
  * Carnahan and Starling (1970) N. F. Carnahan and K. E. Starling, J. Chem. Phys. 53, 600 (1970).
  * Adams (1974) D. J. Adams, Mol. Phys. 28, 1241 (1974).
  * Ashwin  _et al._ (2012) S. S. Ashwin, J. Blawzdziewicz, C. S. O’Hern, and M. D. Shattuck, Phys. Rev. E 85, 061307 (2012).
  * Martiniani  _et al._ (2016) S. Martiniani, K. J. Schrenk, J. D. Stevenson, D. J. Wales, and D. Frenkel, Phys. Rev. E 93, 012906 (2016).
  * Stillinger, DiMarzio, and Kornegay (1964) F. H. Stillinger, E. A. DiMarzio, and R. L. Kornegay, J. Chem. Phys. 40, 1564 (1964).
  * Stillinger (1995) F. H. Stillinger, Science 267, 1935 (1995).
  * Torquato and Jiao (2010) S. Torquato and Y. Jiao, Phys. Rev. E 82, 061302 (2010).
  * Parisi and Zamponi (2005) G. Parisi and F. Zamponi, J. Chem. Phys. 123, 144501 (2005).
  * Parisi and Zamponi (2010) G. Parisi and F. Zamponi, Rev. Mod. Phys. 82, 789 (2010).
  * Berthier, Jacquin, and Zamponi (2011) L. Berthier, H. Jacquin, and F. Zamponi, Phys. Rev. E 84, 051103 (2011).
  * Debenedetti and Stillinger (2001) P. G. Debenedetti and F. H. Stillinger, Nature 410, 259 (2001).
  * Frenkel (2014) D. Frenkel, Mol. Phys. 112, 2325 (2014).
  * Meng  _et al._ (2010) G. Meng, N. Arkus, M. P. Brenner, and V. N. Manoharan, Science 327, 560 (2010).
  * Stillinger and Salsburg (1969) F. H. Stillinger and Z. W. Salsburg, J. Stat. Phys. 1, 179 (1969).
  * Speedy (1993) R. J. Speedy, Mol. Phys. 80, 1105 (1993).
  * Speedy (1998b) R. J. Speedy, J. Phys.: Condens. Matter 10, 4387 (1998b).
  * Baranau and Tallarek (2014) V. Baranau and U. Tallarek, Soft Matter 10, 3826 (2014).
  * Zinchenko (1994) A. Zinchenko, J. Comput. Phys. 114, 298 (1994).
  * Widom (1963) B. Widom, J. Chem. Phys. 39, 2808 (1963).
  * Speedy (1981) R. J. Speedy, J. Chem. Soc., Faraday Trans. 2 77, 329 (1981).
  * Speedy and Reiss (1991) R. J. Speedy and H. Reiss, Mol. Phys. 72, 999 (1991).
  * Baranau and Tallarek (2016) V. Baranau and U. Tallarek, J. Chem. Phys. 144, 214503 (2016).


# Peakless Motzkin paths of bounded height

Helmut Prodinger  Department of Mathematics, University of Stellenbosch 7602,
Stellenbosch, South Africa and NITheCS (National Institute for Theoretical and
Computational Sciences), South Africa<EMAIL_ADDRESS>
###### Abstract.

There was recent interest in Motzkin paths without peaks (peak: up-step
followed immediately by down-step); additional results about this interesting
family is worked out. The new results are the enumeration of such paths that
live in a strip $[0..\ell]$, and as consequence the asymptotics of the average
height, which is given by $2\cdot 5^{-1/4}\sqrt{\pi n}$. Methods include the
kernel method and singularity analysis of generating functions.

###### Key words and phrases:

Motzkin paths, peakless, height, generating functions, asymptotics

###### 2010 Mathematics Subject Classification:

05A15

##  1\. Introduction

Motzkin paths are cousins of the more famous Dyck paths. They appear first in
[7]. In the encyclopedia [12] they are sequence A001006, with many references
given. They consist of up-steps $U=(1,1)$, down-steps $D=(1,-1)$ and
horizontal (flat) steps $F=(1,0)$. Slightly different notations are also in
use. They start at the origin and must never go below the $x$-axis. Usually
one requires the path to end on the $x$-axis as well, but occasionally one
uses the term _Motzkin path_ also for paths that end on a different level.
Figure 1 shows all Motzkin paths of 4 steps (=length 4).

Figure 1. All 9 Motzkin of 4 steps (length 4).

An important concept is the _height_ of a path. It is the maximal
$y$-coordinate when scanning the path (from left to right, say). For the paths
in Figure 1, the heights are (in this order) $0,1,1,1,1,1,1,1,2$. The average
height of all Motzkin paths of length $n$ was computed in an early paper of
the present writer [9].

Recently, I learned from the paper [2] that there is interest in _peakless
Motzkin paths._ A peak in a Motzkin path is a sequence of an up-step followed
immediately by a down-step. Figure 2 indicates all peaks in the list of
Motzkin paths of length 4. The enumerating sequence is A004148 in [12], where
one can find several references; one recent paper about the subject is [3]. A
general discussion about forbidden patterns, concentrating on analytic
aspects, is in [1].

Figure 2. Motzkin paths with peaks indicated.

The paths without peaks are called peakless, and there are four of them, as
shown in Figure 3.

Figure 3. Peakless Motzkin paths of length 4.

##  2\. A warmup: Enumeration of peakless Motzkin paths via the kernel method

$0$$1$$2$$3$$4$$5$$6$$7$$8$ Figure 4. Graph (automaton) to recognize peakless
Motzkin paths. Starting at the origin and ending at nodes labelled 0
corresponds to Motzkin paths, and ending at a node labelled $k$ to a path that
ends at level $k$.

We will use the following generating functions: $[z^{n}]f_{i}(z)$ is the
number of peakless paths ending at state $i$ in the top layer (Figure 4),
$[z^{n}]g_{i}(z)$ is the number of peakless paths ending at state $i$ in the
bottom layer; for convenience, we mostly write $f_{i}$ and $g_{i}$. The
following recursions can be read off the automaton, by considering the last
step separately.

$\displaystyle f_{0}$ $\displaystyle=1+zf_{0}+zf_{1}+zg_{0},$ $\displaystyle
f_{i}$ $\displaystyle=zf_{i}+zf_{i+1}+zg_{i},\quad i\geq 1,$ $\displaystyle
g_{0}$ $\displaystyle=0,$ $\displaystyle g_{i+1}$
$\displaystyle=zf_{i}+zg_{i},\quad i\geq 0.$

To solve this system, one introduces double generating functions

$F(u,z)=\sum_{i\geq 0}u^{i}f_{i}(z)\quad\text{and}\quad G(u,z)=\sum_{i\geq
0}u^{i}g_{i}(z).$

Again, for convenience, we mostly write $F(u)$ and $G(u)$. By summing the
recursions, we find

$\displaystyle F(u)$ $\displaystyle=1+zF(u)+\frac{z}{u}(F(u)-F(0))+zG(u),$
$\displaystyle G(u)$ $\displaystyle=zuF(u)+zuG(u)=\frac{zuF(u)}{1-zu}.$

Eliminating one function, we are left to solve (note that $F(0)=f_{0}$)

$F(u)=1+zF(u)+\frac{z}{u}(F(u)-F(0))+\frac{z^{2}uF(u)}{1-zu}.$

Rewriting the functional equation, we find

$F(u)=\frac{(-u+zF(0))(1-zu)}{zu^{2}+(z-z^{2}-1)u+z}=\frac{(-u+zF(0))(1-zu)}{z(u-s_{1})(u-s_{2})}$

with

$s_{1}=\frac{1-z+z^{2}+\sqrt{(1+z+z^{2})(1-3z+z^{2})}}{2z}$ (1)

and

$s_{2}=\frac{1-z+z^{2}-\sqrt{(1+z+z^{2})(1-3z+z^{2})}}{2z}.$ (2)

Note that $s_{1}s_{2}=1$. Plugging $u=0$ into that equation does not help, but
one of the factors from the denominator can be cancelled. This is (a simple
instance of) the kernel method. Some twenty years ago, I collected various
related examples [10], and many more just recently [11].

Since $s_{2}=z+z^{2}+\cdots$, the factor $(u-s_{2})$ must cancel, since
$1/(u-s_{2})$ would not have a power series expansion around $u$, $z$ close to
zero. Performing the cancellation, we get

$F(u)=\frac{zs_{2}-1+zu-z^{2}F(0)}{z(u-s_{1})}\quad\text{and}\quad
F(0)=\frac{zs_{2}-1-z^{2}F(0)}{-zs_{1}}.$

Solving,

$\displaystyle F(0)=f_{0}$
$\displaystyle=\frac{s_{2}}{z}=\frac{1-z+z^{2}-\sqrt{(1+z+z^{2})(1-3z+z^{2})}}{2z^{2}}$
$\displaystyle=1+z+z^{2}+2z^{3}+4z^{4}+8z^{5}+17z^{6}+\cdots$

We are mostly interested in

$H(u):=F(u)+G(u)=\frac{-u+zF(0)}{zu+z-{z}^{2}u-u+z{u}^{2}},$

although $F(u)$ and $G(u)$ could be computed separately as well. The
coefficients of $[u^{k}]H(u)$ are just enumerating all peakless Motzkin paths
ending on level $k$. Hence

$H(u)=\frac{-u+zF(0)}{zu+z-{z}^{2}u-u+z{u}^{2}}=\frac{-1}{z(u-s_{1})}.$

Expanding and noting that $s_{1}=1/s_{2}$, we find

$[u^{k}]H(u)=\frac{s_{2}^{k+1}}{z}.$

This could be seen as well by a canonical decomposition of a peakless Motzkin
path according to the last return to the $x$-axis.

The denominator of $H(u)$ has a special significance; if we write
$h_{k}=[u^{k}]H(u)$, the recursion for these quantities can be read off from
the denominator:

$zh_{k}+(z-z^{2}-1)h_{k-1}+zh_{k-2}=0,$ (3)

this can be checked directly as well by inserting $h_{k}=s_{2}^{k+1}/z$ and
simplifying.

The sequence enumerating peakless Motzkin paths is A004148 in [12]. If we call
them $m(n)=[z^{n}]s_{2}/z$, then the software Gfun, implemented in Maple,
produces the recursion

$nm(n)-(2n+3)m(n+1)-(n+3)m(n+2)-(2n+9)m(n+3)+(n+6)m(n+4)=0,$

with initial values $m(0)=1$, $m(1)=1$, $m(2)=1$, $m(3)=2$.

A second order linear recursion with constant coefficients is driven by the
characteristic equation and its two roots. Not surprisingly, they are $s_{1}$
and $s_{2}$.

For completeness, we mention the asymptotics of the coefficients $m(n)$ of

$\frac{s_{2}}{z}=\frac{1-z+z^{2}-\sqrt{(1+z+z^{2})(1-3z+z^{2})}}{2z^{2}}.$

This is a standard application of _singularity analysis of generating
functions_ , as described in [5]. First, we must consider the closest
singularity to the origin. The candidates are the solutions of
$(1+z+z^{2})(1-3z+z^{2})=0$. There are two complex solutions of absolute value
1, which are irrelevant, and then $\frac{3\pm\sqrt{5}}{2}$. The relevant value
is $\varrho=\frac{3-\sqrt{5}}{2}$; note that $1/\varrho=\frac{3+\sqrt{5}}{2}$,
which is the square of the golden ratio $\phi=\frac{1+\sqrt{5}}{2}$ from the
Fibonacci fame.

The local expansion around $z\sim\varrho$ looks like

$\frac{s_{2}}{z}\sim\frac{1}{\varrho}-\frac{5^{1/4}}{\varrho}\sqrt{1-\frac{z}{\varrho}},$

and following the principles of singularity analysis we might translate this
to the coefficients:

$[z^{n}]\frac{s_{2}}{z}\sim\frac{5^{1/4}\varrho^{-n-1}}{2\sqrt{\pi}n^{3/2}}.$

##  3\. Peakless Motzkin paths of bounded height

We fix a parameter $\ell\geq 0$ and postulate that $[u^{j}]H(u)=0$ for
$j>\ell$. This means that states $\ell+1,\ell+2,\dots$ (on both layers) can
never be reached. The recursion, compare (3)

$zh_{k+1}+(z-z^{2}-1)h_{k}+zh_{k-1}=0$

is then best written as a matrix equation:

$\begin{pmatrix}z-z^{2}-1&z&0&0&\dots\\\ z&z-z^{2}-1&z&0&0&\dots\\\
0&z&z-z^{2}-1&z&0&\dots\\\ \vdots&\vdots&\vdots&\ddots&\ddots&\\\
&&&&z&z-z^{2}-1\end{pmatrix}\begin{pmatrix}h_{0}\\\ h_{1}\\\ h_{2}\\\
\vdots\\\ h_{\ell}\end{pmatrix}=\begin{pmatrix}-1\\\ 0\\\ 0\\\ \vdots\\\
0\end{pmatrix}$

Let $\mathscr{D}_{\ell}$ be the determinant of the $(\ell+1)\times(\ell+1)$
matrix. $\mathscr{D}_{0}=z-z^{2}-1$,
$\mathscr{D}_{1}=(z-z^{2}-1)^{2}-z^{2}=(1-z)^{2}(1+z^{2})$. It is a bit easier
to work with $\mathscr{D}_{-1}=1$ instead of $\mathscr{D}_{1}$. The solution
is, by standard methods,

$\mathscr{D}_{\ell}=-\frac{1}{W}(-zs_{1})^{\ell+2}+\frac{1}{W}(-zs_{2})^{\ell+2},$

where we use the abbreviation $W=\sqrt{(1+z+z^{2})(1-3z+z^{2})}$.

The quantity $h_{0}$, describing _all_ paths, restricted as described, can, by
Cramer’s rule, be written as

$\displaystyle h_{0}=\frac{-\mathscr{D}_{\ell-1}}{\mathscr{D}_{\ell}}$
$\displaystyle=-\frac{-\frac{1}{W}(-zs_{1})^{\ell+1}+\frac{1}{W}(-zs_{2})^{\ell+1}}{-\frac{1}{W}(-zs_{1})^{\ell+2}+\frac{1}{W}(-zs_{2})^{\ell+2}}$
$\displaystyle=\frac{1}{z}\frac{-(-s_{1})^{\ell+1}+(-s_{2})^{\ell+1}}{(-s_{1})^{\ell+2}-(-s_{2})^{\ell+2}}=\frac{1}{z}\frac{s_{1}^{\ell+1}-s_{2}^{\ell+1}}{s_{1}^{\ell+2}-s_{2}^{\ell+2}}.$

In the limit $\ell\to\infty$, $h_{0}=\dfrac{1}{zs_{1}}=\dfrac{s_{2}}{z}$. This
is, as we have seen already, the generating function of peakless Motzkin paths
without boundary. The other functions $h_{i}$ could be computed by Cramer’s
rule as well, but we concentrate only on the paths that return to the origin
and are bounded by $\ell$. They live in the strip $[0..\ell]$. At this stage,
we drop the ‘$0$’ from the notation and make the ‘$\ell$’ explicit by writing
$A_{n,\ell}$. We summarize the results:

######  Theorem 1.

The number of peakless Motzkin paths of length $n$ (returning to the
$x$-axis), bounded by $\ell$, is given by

$A_{n,\ell}=[z^{n}]\frac{1}{z}\frac{s_{1}^{\ell+1}-s_{2}^{\ell+1}}{s_{1}^{\ell+2}-s_{2}^{\ell+2}},$

with the functions $s_{1}$ and $s_{2}$ given in (1) and (2). The formula is
only correct for $\ell\geq 1$; if $\ell=0$, there is only one such path of
length of $n$, namely consisting of flat steps only; the generating function
must then be replaced by $\dfrac{1}{1-z}$.

We will use the notations ($\ell\geq 1$)

$A_{\ell}(z)=\frac{1}{z}\frac{s_{1}^{\ell+1}-s_{2}^{\ell+1}}{s_{1}^{\ell+2}-s_{2}^{\ell+2}}\quad\text{and}\quad
A_{\infty}(z)=\dfrac{s_{2}}{z}.$

As can be checked directly (best with a computer), there is the recursion of
the continued fraction type

$A_{\ell}=\frac{1}{1-z+z^{2}-z^{2}A_{\ell-1}},$

and thus

$A_{1}=\cfrac{1}{1-z+z^{2}-\cfrac{z^{2}}{1-z+z^{2}}},\quad
A_{2}=\cfrac{1}{1-z+z^{2}-\cfrac{z^{2}}{1-z+z^{2}-\cfrac{z^{2}}{1-z+z^{2}}}}$

and so on. As discussed, the formula for $A_{0}=\frac{1}{1-z}$ is slightly
different.

Continued fraction expansions are very common in the context of generating
functions of lattice paths of bounded height; the first example is (perhaps)
[4].

The total generating function $A_{\infty}$ has a pretty expansion as well,
viz.

$A_{\infty}=1+\cfrac{z}{1-\cfrac{z}{1-\cfrac{z}{1-\cfrac{z^{3}}{1-\cfrac{z}{1-\cfrac{z}{1-\cfrac{z^{3}}{1-\ddots}}}}}}}$

this can be checked directly. I am wondering whether this might have an easy
combinatorial interpretation.

Now we consider the _average height_ of peakless Motzkin paths of length $n$,
assuming that all of them are equally likely. As always, when enumerating the
average height, the relevant formula is

$\frac{[z^{n}]\sum_{\ell\geq
0}\big{(}A_{\infty}(z)-A_{\ell}(z)\big{)}}{[z^{n}]A_{\infty}(z)}.$

The difference $A_{\infty}(z)-A_{\ell}(z)$ enumerates paths of height $>\ell$,
or

$\dfrac{s_{2}}{z}-\frac{1}{z}\frac{s_{1}^{\ell+1}-s_{2}^{\ell+1}}{s_{1}^{\ell+2}-s_{2}^{\ell+2}}=\frac{1}{z}\frac{(1-s_{2}^{2})s_{2}^{\ell+1}}{s_{1}^{\ell+2}-s_{2}^{\ell+2}}=\frac{W}{z^{2}}\frac{s_{2}^{\ell+2}}{s_{1}^{\ell+2}-s_{2}^{\ell+2}}$

This must be expanded around $z=\varrho$; we decrease $\ell$ by one, since
then we enumerate paths of height $\geq\ell$. Since $s_{2}\sim 1$, we find the
approximation

$\frac{W}{\varrho^{2}}\frac{1}{s_{1}^{\ell+1}-1}=\frac{W}{\varrho^{2}}\frac{s_{2}^{\ell+1}}{1-s_{2}^{\ell+1}}=\frac{W}{\varrho^{2}}\sum_{k\geq
1}s_{2}^{(\ell+1)k}$

A local expansion yields

$\frac{W}{\varrho^{2}}\sim\frac{2\cdot
5^{1/4}}{\varrho}\Big{(}1-\frac{z}{\varrho}\Big{)}^{1/2}.$

From our earlier computation,

$s_{2}\sim 1-5^{1/4}\Big{(}1-\frac{z}{\varrho}\Big{)}^{1/2}.$

From [6] we conclude that

$\sum_{k\geq
1}\frac{s_{2}^{k}}{1-s_{2}^{k}}\sim-\frac{\log(1-s_{2})}{1-s_{2}},$

as $s_{2}\to 1$. Our paper [6] has many more technical details about a similar
scenario. Putting both expansions together (we don’t care about a missing term
in the sum as we only want to work out the leading term),

$\displaystyle\frac{W}{\varrho^{2}}\sum_{k,\ell\geq 1}s_{2}^{(\ell+1)k}$
$\displaystyle\sim\frac{2\cdot
5^{1/4}}{\varrho}\Big{(}1-\frac{z}{\varrho}\Big{)}^{1/2}\cdot\frac{-\log(1-s_{2})}{1-s_{2}}$
$\displaystyle\sim\frac{2\cdot
5^{1/4}}{\varrho}\Big{(}1-\frac{z}{\varrho}\Big{)}^{1/2}\cdot\frac{-\log\Big{(}5^{1/4}\big{(}1-\dfrac{z}{\varrho}\big{)}^{1/2}\Big{)}}{5^{1/4}\big{(}1-\dfrac{z}{\varrho}\big{)}^{1/2}}$
$\displaystyle\sim-\frac{2}{\varrho}\log\Big{(}5^{1/4}\big{(}1-\dfrac{z}{\varrho}\big{)}^{1/2}\Big{)}\sim-\frac{2}{\varrho}\log\Big{(}1-\dfrac{z}{\varrho}\Big{)}^{1/2}\sim-\frac{1}{\varrho}\log\Big{(}1-\dfrac{z}{\varrho}\Big{)}.$

By singularity analysis (transfer theorem) we find that

$[z^{n}]\frac{W}{\varrho^{2}}\sum_{k,\ell\geq
1}s_{2}^{(\ell+1)k}\sim-[z^{n}]\frac{1}{\varrho}\log\Big{(}1-\dfrac{z}{\varrho}\Big{)}\sim\frac{\varrho^{-n-1}}{n}.$

As discussed before, the total number of peakless Motzkin paths of length $n$
is asymptotic to

$[z^{n}]\frac{s_{2}}{z}\sim\frac{5^{1/4}\varrho^{-n-1}}{2\sqrt{\pi}n^{3/2}}.$

and for the average height we have to consider the quotient of the last two
expressions, which is

$\frac{\varrho^{-n-1}}{n}\frac{2\sqrt{\pi}n^{3/2}}{5^{1/4}\varrho^{-n-1}}=\frac{2\sqrt{\pi
n}}{5^{1/4}}.$

The numerical constant $2/5^{1/4}=1.337480610$. This can be compared with the
average height of _all_ Motzkin paths of length $n$, which is asymptotic to
$\sqrt{\frac{\pi n}{3}}$, see [8]; $3^{-1/2}=0.5773502693$.

## References

  * [1] Andrei Asinowski, Axel Bacher, Cyril Banderier, and Bernhard Gittenberger.  Analytic combinatorics of lattice paths with forbidden patterns, the vectorial kernel method, and generating functions for pushdown automata.  Algorithmica, 82(3):386–428, 2020.
  * [2] Jean-Luc Baril and José Luis Ramírez.  Fibonacci and Catalan paths in a wall.  2023\.
  * [3] Naiomi Cameron and Everett Sullivan.  Peakless Motzkin paths with marked level steps at fixed height.  Discrete Mathematics, 344.1:112154, 2021.
  * [4] Nicolaas G. de Bruijn, Donald E. Knuth, and Stephen O. Rice.  The average height of planted plane trees.  In Graph Theory and Computing, pages 15–22. Academic Press, New York, 1972.
  * [5] Philippe Flajolet and Robert Sedgewick.  Analytic Combinatorics.  Cambridge University Press, Cambridge, 2009.
  * [6] Clemens Heuberger, Helmut Prodinger, and Stephan Wagner.  The height of multiple edge plane trees.  Aequationes Math., 90(3):625–645, 2016.
  * [7] Theodore S. Motzkin.  Relations between hypersurface cross ratios, and a combinatorial formula for partitions of a polygon, for permanent preponderance, and for non-associative products.  Bulletin of the American Mathematical Society, 54 (4):352–360, 1948\.
  * [8] Helmut Prodinger.  The average height of a stack where three operations are allowed and some related problems.  J. Combin. Inform. System Sci., 5(4):287–304, 1980.
  * [9] Helmut Prodinger.  The height of planted plane trees revisited.  Ars Combin., 16(B):51–55, 1983.
  * [10] Helmut Prodinger.  The kernel method: a collection of examples.  Sém. Lothar. Combin., 50:Art. B50f, 19, 2003/04.
  * [11] Helmut Prodinger.  A walk through my lattice path garden.  Sém. Lothar. Combin., 87b:Art. B87b,49, 2023.
  * [12] Neil J. A. Sloane and The OEIS Foundation Inc.  The On-Line Encyclopedia of Integer Sequences, 2023.


11institutetext: Radboud University, The Netherlands 22institutetext: Leiden
University, The Netherlands
22email<EMAIL_ADDRESS>
# Learning to Learn from Mistakes: Robust Optimization for Adversarial Noise

Alex Serban 11 Erik Poll 11 Joost Visser 22

###### Abstract

Sensitivity to adversarial noise hinders deployment of machine learning
algorithms in security-critical applications. Although many adversarial
defenses have been proposed, robustness to adversarial noise remains an open
problem. The most compelling defense, adversarial training, requires a
substantial increase in processing time and it has been shown to overfit on
the training data. In this paper, we aim to overcome these limitations by
training robust models in low data regimes and transfer adversarial knowledge
between different models. We train a meta-optimizer which learns to robustly
optimize a model using adversarial examples and is able to transfer the
knowledge learned to new models, without the need to generate new adversarial
examples. Experimental results show the meta-optimizer is consistent across
different architectures and data sets, suggesting it is possible to
automatically patch adversarial vulnerabilities.

###### Keywords:

Adversarial Examples Meta Learning Machine Learning.

##  1 Introduction

Machine learning (ML) algorithms exhibit low robustness against intentionally
crafted perturbations [8] or high invariance for distinct inputs [12]. From a
security standpoint this means an attacker can craft inputs that look similar
but cause ML algorithms to misbehave or find distinct inputs that give the
same result. From a safety standpoint, this means that ML algorithms are not
robust against perturbations close to an input or are inflexible to changes in
the operational environment.

Creating defenses (or finding robust counterparts for ML models) has received
increased attention, in particular in the field of adversarial examples.
Although many defenses have been proposed, with the most notable results using
formal methods to guarantee robustness [10] most solutions overfit on the
training data and behave poorly against data outside this distribution [30,
20]. Theoretical investigations suggest these results are expected because
training robust models requires more data [6], more computational resources
[3] or accepting a trade off between accuracy and robustness [24]. Moreover,
solutions to one vulnerability have a negative impact on others [11].

On a different path, designing ML algorithms capable to rapidly adapt to
changes in the operational environment, to adapt to distribution shifts or
capable to learn from few samples is an active research field. Particularly,
the field of meta-learning investigate optimization algorithms learned from
scratch for faster training [1], with less resources [19] and capable of rapid
adaptation [7].

In this paper, we show that meta-learning algorithms can be used to extract
knowledge from a model’s vulnerability to adversarial examples and transfer it
to new models. We train a meta-optimizer to learn how to robustly optimize
other models using adversarial training. Later, when asked to optimize new
models without seeing adversarial examples, the trained meta-optimizer can do
it robustly. This process is analogous to learning a regularization term for
adversarial examples, instead of manually designing one. The experimental
results suggest a broader horizon, in which algorithms learn how to
automatically repair or treat vulnerabilities without explicit human design.

The rest of the paper is organized as follows. In Section 2 we introduce
prerequisites and related work. Section 3 formalizes meta-learning and the
adversarial training problems and gives implementation details. Section 4
presents experimental results on two distinct data sets, followed by a
discussion in Section 5 and conclusions in Section 6.

##  2 Background and Related Work

We focus on the task of supervised classification, i.e. given a set of inputs
from the input space ${\mathcal{X}}$ and a set of labels from the space
${\mathcal{Y}}$ a ML algorithm attempts to find a mapping from
$f:{\mathcal{X}}\rightarrow{\mathcal{Y}}$ which minimizes the number of
misclassified inputs. We assume that ${\mathcal{X}}$ is a metric space such
that a distance function $d(\cdot)$ between two points of the space exists.
The error made by a prediction $f({\bm{x}})=\hat{y}$ when the true label is
$y$ is measured by a loss function
$l:{\mathcal{Y}}\times{\mathcal{Y}}\rightarrow{\mathbb{R}}$ with non-negative
values when the labels are different and zero otherwise. $f(\cdot)$ is defined
over a hypothesis space ${\mathcal{F}}$ which encompasses any mapping from
${\mathcal{X}}$ to ${\mathcal{Y}}$ and can take any form – e.g. a linear
function or a neural network. Through learning, a ML algorithm selects
$f^{*}(\cdot)$ from ${\mathcal{F}}$ such that the expected empirical loss on a
training set ${\mathcal{S}}$ consisting of pairs of samples
$({\bm{x}}_{i},y_{i})\sim{\mathcal{Z}}={\mathcal{X}}\times{\mathcal{Y}}$ is
minimal.

A robust solution to the minimization problem above involves immunizing it
against uncertainties in the input space. In the adversarial examples setting,
uncertainties are modeled in the space around an input ${\bm{x}}_{u}$:
${\mathcal{U}}_{{\bm{x}}}=\\{{\bm{x}}|d({\bm{x}},{\bm{x}}_{c})\leq\epsilon\\}$,
where $d(\cdot)$ is a distance function defined on the metric space
${\mathcal{X}}$. The robust counterpart of the learning problem is defined as:

$f=\operatorname*{arg\,min}_{f\in{\mathcal{F}}}\mathbb{E}_{({\bm{x}},y)\sim{\mathcal{S}}}[\max_{{\bm{x}}_{u}\in{\mathcal{U}}_{{\bm{x}}}}l(f({\bm{x}}_{u}),y)],$
(1)

where ${\bm{x}}_{u}$ is a realization of ${\bm{x}}$ in the uncertainty set
described by $d(\cdot)$. A common distance function used in the field of
adversarial examples is the $l_{p}$-norm. Arguably, this metric can not
capture task-specific semantic information, but it is suitable for comparative
benchmarks. Moreover, the lack of robust p-norm solutions for complex tasks
makes it hard to believe that other notions of robustness can lead to better
results [4]. In this paper we use the $l_{\infty}$-norm distance.

Finding models robust to adversarial examples is an open question, spanning
two research directions: (1) finding solutions robust to an upper bound on the
perturbation size $\epsilon$ (i.e. no perturbation higher than $\epsilon$ can
cause a misclassification) [8, 16] and (2) finding lower bounds on robustness
(i.e. no adversarial example can be found in a space around an input) [26,
17].

Notable results are obtained by solving the inner maximization problem from
Eq. (1) with projected gradient descent (PGD) [16], by training with over-
approximations of the input [17] or extensively searching for spaces where no
adversarial examples can be found using exact solvers [13]. In all cases, the
benchmarks show overfitting on the training distribution and can be bypassed
by samples outside it [30].

Until recently, solving the outer minimization objective from Eq. (1) relied
on static, hand designed algorithms, such as stochastic gradient descent or
ADAM. This line of research is driven by the no free lunch theorem of
optimization which states that, on average, in combinatorial optimization no
algorithm can do better than a random strategy; suggesting that designing
specific algorithms for a class of problems is the only way to improve
performance.

Recent advancements in the field of meta-learning have taken a different
approach, posing the problem of algorithm design dynamically and modeling it
with neural networks [1] or as a reinforcement learning (RL) problem [14]. In
both cases, the algorithms show empirically faster convergence and the ability
to adapt to different tasks.

For example, in [1] the hand designed update rules are replaced with a
parameterized function modeled with a recurrent neural network (RNN). During
training, both the _optimizer_ and the _optimizee_ parameters are adjusted.
The algorithm design now becomes a learning problem, allowing to specify
constraints through data examples. The method has roots in previous research
[22], where initial RNNs were designed to update their own weights. The
results of backprobagation from one network were only later fed to a second
network which learns to update the first [28]. Andrychowicz et al. build on
previous work using a different learning architecture. Similarly, [19] used
this paradigm to train neural networks in a few shot regime and [21] augmented
it with memory networks. Meta-learning has also shown promising results in
designing algorithms for fast adaptation using gradient information [7].

Instead of using RNNs, [14] formulate the optimization problem as a RL problem
where the reward signal represents successful minimization and train an agent
using guided policy search. Later, the authors refine their method for
training neural networks [15]. In both cases the agent learns to optimize
faster than hand designed algorithms and exhibits better stability.

Recent research in adversarial examples has also tackled the need to decrease
training resources by either accumulating perturbations [23, 27] or by
restricting back-propagation to some layers of the model [29]. While the
latter method requires forward and backward passes, the former reduces the
need to do backward passes in order to generate adversarial examples. and not
to abstract or transfer information from adversarial training.

##  3 Learning to Optimize Robustly

Meta-learning frames the learning problem at two levels: acquiring higher
level (meta) information about the optimization landscape and applying it to
optimize one task. In this paper, we are interested to learn robust update
rules and transfer this knowledge to new tasks without additional constraints.
We focus on training robust ML models through adversarial training, which is
one of the most effective defenses against upper bounded adversarial examples
[16]. Because generating adversarial examples during training is time
consuming, especially for iterative procedures and can only provide robustness
for the inputs used during training [30], through meta-learning we learn to
optimize robustly without explicit regularization terms and transfer the
knowledge to new tasks, without the need to generate new adversarial examples.

At a high level, adversarial training discourages _local_ sensitive behavior –
in the vicinity of each input in the training set – by guiding the model to
behave constantly in regions surrounding the training data. The regions are
defined by a chosen uncertainty set (as in Eq. (1)). This procedure is
equivalent to adding a prior that defines local constancy, for each model we
want to train. In most cases, specifying this prior is not trivial and
requires the design of new regularization methods [18], new loss functions
[26] or new training procedures [17]. In this paper we take a different
approach and try to learn a regularization term automatically, using meta-
learning. During the meta-knowledge acquisition phase, the meta-optimizer
learns to perform the updates robustly using adversarial training. Later, the
knowledge acquired is transferred to new models using the meta-optimizer to
train new models, without generating adversarial examples. In the next
paragraphs we describe the meta-optimizer, some implementation details and the
adversarial training procedure.

#### Learning to Optimize.

The classification function defined in Section 2, $f(\cdot)$, is parametrized
by a set of parameters $\theta$. Upon seeing new data, we update the
parameters in order to minimize the prediction errors. The update consists in
moving one step in the opposite direction of the gradient:

$\theta_{t+1}=\theta_{t}-\eta\nabla_{\theta_{t}}l(\cdot),$ (2)

where $\eta$ (the learning rate) determines the size of the step. Different
choices of $\eta$ or ways to automatically adapt it results in different
optimization algorithms such as stochastic gradient descent or ADAM.

In order to avoid overfitting or impose additional constraints such as
constancy around inputs, it is common to add a regularization term to the loss
function which will be back-propagated and reflected on all parameter updates.
Instead of looking for regularization terms manually, we use a method to
automatically learn robust update steps with regularization included.

As discussed in Section 2, a parameterized update rule has been previously
represented with a RNN [1, 19] or as a RL problem [14]. In this paper, we
follow an approach similar to [1] and model the update rule with a RNN with
long short-term memory (LSTM) cells:

$\theta_{t}=\theta_{t-1}+c_{t},$

where :

$c_{t}=f_{t}c_{t-1}+i_{t-1}\tilde{c}_{t},$ (3)

is the output of an LSTM network $m$ with input
$\nabla_{\theta_{t}}(l(\cdot))$:

$\displaystyle\begin{bmatrix}c_{t}\\\
h_{t+1}\end{bmatrix}=m(\nabla_{t},h_{t},\phi),$ (4)

and $\phi$ are the LSTM’s parameters.

In [19], the authors consider each term in Eq. (3) equivalent to each term in
Eq. (2) – e.g. $f_{t}=1$, $c_{t-1}=\theta_{t}$ – and disentangle the internal
state of the LSTM, $h_{t}$, with special terms for individual updates of
$f_{t}$ and $i_{t}$. This type of inductive bias brings benefits in some cases
and will be further discussed in Section 5. However, in this paper we try to
avoid such biases whenever possible.

#### Parameters Sharing and Gradient Preprocessing.

In order to limit the number of parameters of the optimizer, we follow a
procedure similar to [1, 19] in which for each parameter of the function we
want to optimize, $\theta_{i:n}$, we keep an equivalent internal state of the
optimizer, $h_{i:n}$, but share the weights $\phi$ between all states. This
procedure allows a more compact optimizer to be used and makes the update rule
dependent only on its respective past representation, thus being able to
simulate hand designed optimization features such as momentum.

Moreover, since gradient coordinates can take very distinct values, we apply a
common normalization step in which the gradients are scaled and the
information about their magnitude and their direction is separated:

$\displaystyle\nabla$
$\displaystyle\rightarrow\begin{cases}\left(\frac{\log(|\nabla|)}{p}\,,\operatorname{sign}(\nabla)\right)&\text{if
}|\nabla|\geq e^{-p}\\\ (-1,e^{p}\nabla)&\text{otherwise.}\end{cases}$ (5)

We experiment with different values for $p$ by grid search and observe that
increasing the size of $p$ yields better results when the perturbations are
larger. However, for consistency, we use $p=10$ for all experiments.

The meta-optimizer’s parameters are updated using an equivalent to Eq. (2).
Since its inputs are based on the gradient information of the function to be
optimized (also called the optimizee) the updates will require second order
information about it (taking the gradient of the gradient). This information
is commonly used for meta-learning – e.g. in [7, 25] – and will be further
discussed in Section 5. However, in this paper only first order information is
used, corresponding to limiting the propagation of the optimizer’s gradient on
the optimizee parameters (or stopping the gradient flow in the computational
graph).

#### Adversarial Training.

Adversarial training is still one of the most effective defenses against
adversarial examples. This procedure is equivalent to adding a regularization
term to the loss function of the optimizee, corresponding to:

$\tilde{l}(\cdot)=\alpha
l({\bm{\theta}},{\bm{x}},y)+(1-\alpha)l({\bm{\theta}},{\bm{x}}^{\prime},y),$
(6)

where ${\bm{x}}^{\prime}$ is an adversarial example generated from input
${\bm{x}}$ and $\alpha$ is a parameter which controls the contribution of the
adversarial loss. In some cases adversarial training is performed using only
the adversarial loss, corresponding to $\alpha=0$ in Eq. (6), e.g. when
training with the worst adversarial loss as in Eq. (1).

Several methods to generate adversarial examples have been proposed, ranging
from fast, less precise, attacks [8] to strong, adaptive attacks [5]. However,
due to processing constraints, only some algorithms are used in adversarial
training: the fast gradient sign method (FGSM) attack – which moves one step
in the direction of the gradient – and the PGD attack – which approximates the
uncertainty set ${\mathcal{U}}_{{\bm{x}}}$ running several steps of the FGSM
and projecting the outcomes on the norm ball surrounding an input, defined by
$\epsilon$.

In this paper, we use both FGSM and PGD during the meta-learning phase to
generate adversarial examples and incorporate them in the training procedure.
Formally, the FGSM attack is defined as:

${\bm{x}}^{\prime}={\bm{x}}+\epsilon\operatorname{sign}(\nabla_{x}l({\bm{\theta}},{\bm{x}},y)).$
(7)

The term $\epsilon$ controls the size of the perturbation added in order to
cause a misclassification. We are aware that defenses against FGSM attacks
sometimes lead to a false sense of protection and to gradient masking [2],
however, we remind that hereby the goal is to reduce the data needed to train
robust models and increase the generalization outside the training
distribution.

The PGD attack is defined as follows:

${\bm{x}}^{\prime}=\prod_{n}({\bm{x}}^{\prime}+\epsilon\operatorname{sign}(\nabla_{x}l({\bm{\theta}},{\bm{x}},y))).$
(8)

If the number of iteration suffices to approximate the space we want to
provide robustness to, training with adversarial examples generated with this
method should protect against any perturbation in the space defined by
$\epsilon$. In order to obtain better approximations, we can increase the
iteration number, $n$. However, this comes at increased computational costs.

##  4 Results

In all experiments the optimizer consists of a two-layer LSTM network with a
hidden state size of 20. We compare the results on training two types of
neural networks on two distinct data sets with the adaptive optimizer ADAM.

We focus on two experiments related to training neural networks, as in prior
work on meta-learning [1, 19, 14, 15]. More experiments with minimizing other
functions – e.g. logistic regression – and an integration with the Cleverhans
framework are available in the project’s repository, which can be found at
https://github.com/NullConvergence/Learning2LearnAdv. In all cases, an
optimizer is trained using normal and adversarial examples on a data set and
tested by training a robust optimizee without generating adversarial examples.
Several perturbation sizes are analyzed, as introduced below.

###  4.1 MNIST

We begin by training a small, fully connected, neural network with 20 units
and ReLU activation on the MNIST data set. The perturbations take different
values in the set $\epsilon\in\\{0.05,0.1,0.2,0.3\\}$ for both attacks
introduced earlier (Eq. (7) and Eq. (8)). We experiment with different
learning rates by grid search and find the best to be $0.001$ for the meta-
optimizer. Training is performed using the common cross entropy loss function,
with a batch size of 128. We shuffle the training data set (consisting of
60.000 examples) and divide in two parts equally: the first is used to train a
meta-optimizer using both normal data and adversarial examples and the second
is used to test its performance while training with normal data and testing
with perturbed data. Each experiment ran for 100 steps and the average results
are illustrated in Figures 1 and 2. All experiments are done using
$\alpha=0.5$ in Eq. (6) during training and $\alpha=0.0$ during the meta-
optimizer transfer phase, as first introduced in [8]. In addition to ADAM’s
performance compared to the meta-optimizer, we evaluate the performance of the
meta-optimizer during training and the performance of training a meta-
optimizer using $\alpha=1$ and testing with $\alpha=0$ (L2L and Transfer-NOT
labels in Figures 1 and 2). Figure 1 illustrates the results from generating
adversarial examples using the FGSM method (Eq. (7)). The loss functions start
from approximately the same value because the networks are always initialized
with the same values.

(a) $\epsilon$ = 0.05

(b) $\epsilon$ = 0.1

(c) $\epsilon$ = 0.2

(d) $\epsilon$ = 0.3

Figure 1: The loss landscape when training a neural network on the MNIST data
set, perturbed with FGSM and different perturbation sizes ($\epsilon$). The
meta-optimizer is trained with adversarial examples – label L2L – transferred
to a scenario where training is performed with normal and adversarial data,
but tested with adversarial examples – label L2L-Transfer – and compared with
a meta-optimizer trained with normal data and transferred to adversarial
settings – label Transfer-NOT – and with ADAM. Best seen in color.

(a) $\epsilon$ = 0.05

(b) $\epsilon$ = 0.1

(c) $\epsilon$ = 0.2

(d) $\epsilon$ = 0.3

Figure 2: The loss landscape when training a neural network on the MNIST data
set perturbed with the PGD method and different perturbation sizes
($\epsilon$). The legend is detailed in the caption of Figure 1.

In all cases, the meta-optimizer is able to transfer the information learned
during training and has comparable performance to ADAM (in some cases
performing better). We remind that during testing the optimizer uses normal
data, but the plots are generated by feeding adversarial perturbed data to the
optimizee. This implies that the meta-optimizer proposes update rules which
lead to smooth surfaces around the tested inputs. Moreover, it is able to
learn a robust regularization term during training and transfer it to new
tasks without the need to generate new data. Also, the trained meta-optimizer
exhibits more stable behavior. This brings evidence that adversarial training
leads to more interpretable gradients [24].

When the optimizer is trained only with normal examples, but used to optimize
the model using adversarial examples – Transfer-NOT label in Figure 1 – its
performance decrease significantly. This implies that a meta-optimizer is
domain specific and does not have the general behavior of ADAM, an
observations which will be further discussed in Section 5.

In Figure 2 we illustrate the results from running similar experiments, but
generate adversarial examples using the PGD method from Eq. (8). Training with
PGD is generally performed only using the perturbed examples (corresponding to
$\alpha=0$ in Eq. (6)), as in the original paper [16]. We take a similar
approach in this paper.

The results are consistent with the FGSM method, although the gap between ADAM
and the transferred meta-optimizer is smaller. A constant decrease in
performance is also observed, possibly corresponding to the decrease in
performance specific to adversarial training [24]. Nevertheless, the results
are consistent and bring evidence that the meta-optimizer is able to learn
robust update rules.

###  4.2 CIFAR-10

(a) $\epsilon$ = 0.05

(b) $\epsilon$ = 0.1

(c) $\epsilon$ = 0.2

(d) $\epsilon$ = 0.3

Figure 3: The loss landscape when training a neural network on the CIFAR-10
data set, perturbed with the FGSM method and different perturbation sizes
($\epsilon$). The legend is detailed in the caption of Figure 1.

(a) $\epsilon$ = 0.05

(b) $\epsilon$ = 0.1

(c) $\epsilon$ = 0.2

(d) $\epsilon$ = 0.3

Figure 4: The loss landscape when optimizing a neural network on the CIFAR-10
data set, perturbed with PGD and different perturbation sizes ($\epsilon$).
The legend is detailed in the caption of Figure 1.

We present the results from training a model using both convolution and fully
connected layers on the CIFAR-10 data set. The network consists of three
convolutional layers with kernel size 3, a fully connected layer with 32
hidden units and a logits layer of size 10. All activation functions are ReLU,
the loss is cross-entropy and batch normalization is used in the convolutional
layers. The meta-optimizer is trained using a learning rate of $0.001$
determined through grid search.

Since there are striking differences between the convolution and the linear
layer, we use two sets of parameters for the meta-optimizer – one for
optimizing all convolutional layers and one for the linear layers. Moreover,
since the CIFAR-10 problem is more difficult, we train the optimizee for 1000
steps and present the results in Figure 3 for perturbations generated with
FGSM and in Figure 4 for perturbations generated with PGD. The evaluation was
performed as earlier, using 2-fold cross validation for training a meta-
optimizer in adversarial settings and transfer it to training a model in
normal settings which is tested with adversarial examples.

We observe that in the case of FGSM, the transferred meta-optimizer (label
L2L-Transfer, Figure 3) exhibits similar behavior as in the MNIST experiments:
it has similar and sometimes better performance than ADAM. We remind that, in
this case, no adversarial examples are used during training. The meta-
optimizer trained normally, but tested with adversarial examples (Transfer-NOT
label, Figure 3) performs visibly worse, which strengthens the observation
that meta-learning optimization is domain specific.

Figure 4 shows results from running the same experiment using perturbations
generated with the PGD method, with a number of 7 steps, as in the original
paper [16]. In all cases, the loss improvements are small, although the meta-
optimizer exhibits better performance than ADAM both during training and
testing. However, the improvements in training time are significant since
after training an adversarial meta-optimizer, it can be applied to different
models without the need to execute 7 forward and back propagation steps for
each batch of data.

##  5 Discussion

Typically used to rapidly adapt to new tasks or generalize outside the i.i.d
assumption, meta-learning algorithms show promising results to reduce the
training samples needed for adversarial training. The results presented in
this paper suggest these algorithms can be used in the future to build
adversarial defenses with less computational resources and capable to adapt to
new data.

We hereby note some weaknesses discovered during the process. Although capable
of achieving better performance than hand crafted optimizers [1, 15] and, as
discussed in Section 4, showing promising results in transferring information
about adversarial examples, meta learning algorithms still suffer from broader
generalization. In particular, trained optimizers can not generalize to
different activation functions, or between architectures with noticeable
differences [25]. This means that an optimizer trained for ReLU can not be
used for sigmoid (or other) activation functions. Moreover, if the meta-
optimizer is not trained with specific data that will be later used, it does
not exhibit general behavior. For example, if the meta-optimizer does not use
any adversarial examples during training, but it encounters such examples
during testing, it faces difficulties. This behavior is illustrated in the
figures above with the label Transfer-NOT.

Second order information (taking the gradient of the gradient, as introduced
in Section 3) was not used in this paper. As shown in [25], this information
can help the meta-optimizer better generalize and induce more stable behavior.
However, it also introduces more complexity. Analyzing the trade-off between
the optimizer’s complexity and its ability to learn and transfer knowledge
related to adversarial vulnerabilities is left for future work.

##  6 Conclusions and Future Research

We introduce a method to learn how to optimize ML models robust to adversarial
examples, in low data regimes. Instead of specifying custom regularization
terms, they are learned automatically by an adaptive optimizer. Acquiring meta
information about the optimization landscape under adversarial constraints
allows the optimizer to reuse it for new tasks.

For future research, we propose to train the meta-optimizer concomitantly with
different perturbation types – e.g. $l_{1},l_{2}$-norm – and test if the
optimizer can learn to robustly optimize under all constraints. Other
perturbations, such as naturally occurring perturbations [9] can also be
included. Another research direction is to use the meta-optimizer to refine a
trained model and evaluate if it is possible to robustly regularize it with
less data.

## References

  * [1] Andrychowicz, M., Denil, M., Gomez, S., Hoffman, M.W., Pfau, D., Schaul, T., Shillingford, B., De Freitas, N.: Learning to learn by gradient descent by gradient descent. In: NeurIPS. pp. 3981–3989 (2016)
  * [2] Athalye, A., Carlini, N., Wagner, D.: Obfuscated gradients give a false sense of security: Circumventing defenses to adversarial examples. ICML (2018)
  * [3] Bubeck, S., Price, E., Razenshteyn, I.: Adversarial examples from computational constraints. arXiv:1805.10204 (2018)
  * [4] Carlini, N., Athalye, A., Papernot, N., Brendel, W., Rauber, J., Tsipras, D., Goodfellow, I.J., Madry, A., Kurakin, A.: On evaluating adversarial robustness. arXiv:1902.06705 (2019), http://arxiv.org/abs/1902.06705
  * [5] Carlini, N., Wagner, D.: Towards evaluating the robustness of neural networks. In: 2017 IEEE Symposium on Security and Privacy (S&P). pp. 39–57. IEEE (2017)
  * [6] Cullina, D., Bhagoji, A.N., Mittal, P.: Pac-learning in the presence of evasion adversaries. arXiv:1806.01471 (2018)
  * [7] Finn, C., Abbeel, P., Levine, S.: Model-agnostic meta-learning for fast adaptation of deep networks. In: ICML. pp. 1126–1135 (2017)
  * [8] Goodfellow, I.J., Shlens, J., Szegedy, C.: Explaining and harnessing adversarial examples. arXiv:1412.6572 (2014)
  * [9] Hendrycks, D., Dietterich, T.: Benchmarking neural network robustness to common corruptions and perturbations. ICLR (2019)
  * [10] Huang, X., Kroening, D., Kwiatkowska, M., Ruan, W., Sun, Y., Thamo, E., Wu, M., Yi, X.: Safety and trustworthiness of deep neural networks: A survey. arXiv:1812.08342 (2018)
  * [11] Jacobsen, J.H., Behrmann, J., Carlini, N., Florian Tramer, N.: Exploiting excessive invariance caused by norm-bounded adversarial robustness. ICLR (2019)
  * [12] Jacobsen, J.H., Behrmann, J., Zemel, R., Bethge, M.: Excessive invariance causes adversarial vulnerability. ICLR (2019)
  * [13] Katz, G., Barrett, C., Dill, D.L., Julian, K., Kochenderfer, M.J.: Reluplex: An efficient smt solver for verifying deep neural networks. In: CAV. pp. 97–117. Springer (2017)
  * [14] Li, K., Malik, J.: Learning to optimize. arXiv:1606.01885 (2016)
  * [15] Li, K., Malik, J.: Learning to optimize neural nets. arXiv:1703.00441 (2017)
  * [16] Madry, A., Makelov, A., Schmidt, L., Tsipras, D., Vladu, A.: Towards deep learning models resistant to adversarial attacks. ICLR (2018)
  * [17] Mirman, M., Gehr, T., Vechev, M.: Differentiable abstract interpretation for provably robust neural networks. ICML (2018)
  * [18] Miyato, T., Maeda, S.i., Koyama, M., Nakae, K., Ishii, S.: Distributional smoothing with virtual adversarial training. arXiv:1507.00677 (2015)
  * [19] Ravi, S., Larochelle, H.: Optimization as a model for few-shot learning. ICLR (2017)
  * [20] Rice, L., Wong, E., Kolter, J.Z.: Overfitting in adversarially robust deep learning. arXiv preprint arXiv:2002.11569 (2020)
  * [21] Santoro, A., Bartunov, S., Botvinick, M., Wierstra, D., Lillicrap, T.: Meta-learning with memory-augmented neural networks. In: ICML. pp. 1842–1850 (2016)
  * [22] Schmidhuber, J.: A neural network that embeds its own meta-levels. In: ICNN. pp. 407–412. IEEE (1993)
  * [23] Shafahi, A., Najibi, M., Ghiasi, A., Xu, Z., Dickerson, J., Studer, C., Davis, L.S., Taylor, G., Goldstein, T.: Adversarial training for free! arXiv:1904.12843 (2019)
  * [24] Tsipras, D., Santurkar, S., Engstrom, L., Turner, A., Madry, A.: There is no free lunch in adversarial robustness (but there are unexpected benefits). arXiv:1805.12152 (2018)
  * [25] Wichrowska, O., Maheswaranathan, N., Hoffman, M.W., Colmenarejo, S.G., Denil, M., de Freitas, N., Sohl-Dickstein, J.: Learned optimizers that scale and generalize. In: ICML. pp. 3751–3760 (2017)
  * [26] Wong, E., Kolter, J.Z.: Provable defenses against adversarial examples via the convex outer adversarial polytope. ICML (2018)
  * [27] Wong, E., Rice, L., Kolter, J.Z.: Fast is better than free: Revisiting adversarial training. arXiv preprint arXiv:2001.03994 (2020)
  * [28] Younger, A.S., Hochreiter, S., Conwell, P.R.: Meta-learning with backpropagation. In: IJCNN. vol. 3. IEEE (2001)
  * [29] Zhang, D., Zhang, T., Lu, Y., Zhu, Z., Dong, B.: You only propagate once: Accelerating adversarial training via maximal principle. arXiv:1905.00877 (2019)
  * [30] Zhang, H., Chen, H., Song, Z., Boning, D., Dhillon, I.S., Hsieh, C.J.: The limitations of adversarial training and the blind-spot attack. ICLR (2019)

  *[ML]: Machine learning
  *[PGD]: projected gradient descent
  *[RL]: reinforcement learning
  *[RNN]: recurrent neural network
  *[RNNs]: recurrent neural network
  *[LSTM]: long short-term memory
  *[FGSM]: fast gradient sign method