# Optical properties of radially-polarised twisted light

K. Koksal1,2, M. Babiker1 1Department of Physics, University of York, YO10
5DD, UK  2Physics Department, Bitlis Eren University, Bitlis 13000, Turkey  V.
E. Lembessis3 3Quantum Technology Group, Department of Physics and Astronomy,
College of Science, King Saud University, Riyadh 11451, Saudi Arabia

###### Abstract

We show that, in general, any type of radially-polarised paraxial twisted
optical mode carries only an axial total optical angular momentum (AM)
${\bar{\cal J}_{z}}=\ell{\cal L}_{0}$ where $\ell$ is the winding number and
${\cal L}_{0}$ is a constant. This mode, however, is shown to have zero spin
angular momentum (SAM), so it is endowed only with orbital angular momentum
(OAM) and no SAM. The helicity is found to be proportional to $\ell$, hence
radially-polarised modes display chirality. When applied to a Laguerre-
Gaussian (LG) mode our treatment leads to a total helicity equal to
$(\ell/|\ell|){\cal Q}$, where ${\cal Q}$ is the action constant. The factor
$(\ell/|\ell|)=\pm 1$, depends on the sign, not the magnitude of $\ell$ and so
the result holds for any radially-polarised LG mode however large the
magnitude of its winding number $\ell$ is. The magnitude of the action
constant ${\cal Q}$ and hence the helicity are diminished for all such LG
modes of large beam waist $w_{0}$.

##  I Introduction

A great deal of work has already been carried out on twisted light and on
exploring its interaction with matter Allen  _et al._ (1999, 2003); Zhan
(2004); Franke-Arnold  _et al._ (2008); Andrews and Babiker (2012); Torres and
L. Torner (2011), including atoms and molecules Babiker  _et al._ (2018).
Cylindrical twisted light modes come in a variety of forms, depending on how
they are generated, but all the familiar forms are commonly characterised by
the ubiquitous azimuthal phase factor $e^{i\ell\phi}$ in the generic amplitude
function, which we denote by ${\cal F}$. Here $\ell$ is the winding number and
$\phi$ is the azimuthal angle in cylindrical polar coordinates. Wave
polarisation is an additional, equally significant, ingredient adding another
layer of complexity to twisted light. Familiar wave polarisation types include
linear, circular and elliptical, together with superpositions of those when
dealing with interfering modes. Other forms of polarisation are generally
referred to as vector modes Holmes and Galvez (2019) and include, in
particular, the radial and azimuthal polarisations and superpositions. Here we
are concerned with radially-polarised modes which, since their generation in
laser oscillations Mushiake  _et al._ (1972); Kozawa and Sato (2005), have
been of considerable interest and, especially so, more recently in the context
of twisted light Zhan (2009). A prominent characteristic of radially-polarised
light is that its beams focus into very small waists compared with uniformly-
polarised modes Kozawa and Sato (2007); Dorn  _et al._ (2003); Bashkansky  _et
al._ (2010); Levy  _et al._ (ts).

However, as far as we know, the general properties of radially-polarised
twisted light have not been explored. In particular, there is need to
determine the energy-momentum, spin angular momentum, total angular momentum
and the helicity and chirality of such modes. The purpose of this article is
to set out the formalism needed to determine these properties. We focus on the
paraxial regime and aim to evaluate the properties to leading order. The main
ingredients of the formalism emphasise the need to ensure that the
electromagnetic fields include the longitudinal components for both the
electric and magnetic fields and that Maxwell’s equations are satisfied. This
means that once the magnetic field format is determined, the electric field
follows by application of Maxwell’s curl equation and vice versa the magnetic
field follows from the electric field by application of the second curl
equation.

We aim to deal with monochromatic paraxial radially-polarised vortex modes,
mainly without specifying the kind of vortex mode and consider the cycle-
averaged properties. These are the spin angular momentum, the orbital angular
momentum and the helicity and chirality. We aim to show how the general
results for an arbitrary optical vortex lead to well-defined properties when
the type of vortex mode is specified and we evaluate the properties for the
special case of a Laguerre-Gaussian optical vortex mode.

We find that, typically any radially-polarised paraxial twisted optical mode
is endowed with a total angular momentum (AM) which is only $\ell$ quantised
where $\ell$ is the winding number and it always has zero spin angular
momentum (SAM), so it just carries orbital angular momentum (OAM) and no SAM.
The optical helicity is shown to be proportional to $\ell$, confirming that
radially-polarised modes display chirality as it changes sign with the sign of
the winding number $\ell$. This is a manifestation of geometrical chirality
arising from the geometrical structure of the beam Nechayev  _et al._ (2021).
Our general results are then applied to the case of a Laguerre-Gaussian (LG)
mode and we find that the total helicity is equal to $(\ell/|\ell|){\cal Q}$,
where $Q$ is the action constant, while the factor $(\ell/|\ell|)=\pm 1$ is
the Hopf index of the radially-polarise LG mode. We show that this result is
applicable to any radially-polarised LG mode however large the magnitude of
its winding number $\ell$ is. We also find that the helicity is diminished for
all such LG modes of large beam waist $w_{0}$.

##  II Radially-polarised vortex fields

In cylindrical coordinates the electric and magnetic fields of a radially-
polarised paraxial twisted light mode are derivable from a vector potential in
the form

${\bf A}(\rho,\phi,z)={\mbox{\boldmath$\hat{\rho}$}}{\cal
F}_{\ell,m}(\rho,\phi)e^{ik_{z}z}$ (1)

where carets denote unit vectors, $k_{z}$ is the axial wavevector with the
light travelling along the $+z$ axis and ${\cal F}_{\ell,m}$ is the mode
function which includes both the amplitude and phase functions in terms of the
plane-polar coordinates $(\rho,\phi)$. The mode is labelled by the indices
$\ell$ and $m$, with $\ell$ the winding number and $m$ could be an azimuthal
number, as in the case of Laguerre-Gaussian (LG) optical vortex modes, or
could be redundant, as in the case of optical Bessel modes, but the treatment
is not restricted to LG or Bessel modes and is applicable in general for other
radially-polarised paraxial optical vortex modes.

Using the standard form of the curl of a vector in cylindrical coordinates we
obtain for the magnetic field ${\bf
B}(\rho,\phi,z)={\mbox{\boldmath$\nabla$}}\times{\bf A}(\rho,\phi,z)$. With
the vector potential as defined for the radially-polarised mode in Eq.(1) we
obtain the paraxial magnetic field to leading order

${\bf B}(\rho,\phi,z)=ik_{z}{\mbox{\boldmath$\hat{\phi}$}}{\cal
F}e^{ik_{z}z}-{\mbox{\boldmath$\hat{z}$}}\frac{1}{\rho}\frac{\partial{\cal
F}}{\partial\phi}e^{ik_{z}z}$ (2)

The electric field follows from the magnetic field using Maxwell’s equation

${\mbox{\boldmath$\nabla$}}\times{\bf B}=\frac{1}{c^{2}}\frac{\partial{\bf
E}}{\partial t}$ (3)

We obtain, to the same leading order as for the magnetic field

${\bf E}(\rho,\phi,z)=ick_{z}{\mbox{\boldmath$\hat{\rho}$}}{\cal
F}e^{ik_{z}z}-{\mbox{\boldmath$\hat{z}$}}c\frac{1}{\rho}\frac{\partial(\rho{\cal
F})}{\partial\rho}e^{ik_{z}z}$ (4)

The general paraxial radially-polarised vector modes we are concerned with
here have the following mode function

${\cal F}_{\ell,m}(\rho,\phi)={\tilde{\cal
F}}_{\ell,m}(\rho)e^{i\ell\phi}e^{ik_{z}z}$ (5)

where we have separated the phase dependence, so now ${\tilde{\cal
F}}_{\ell,m}$ is a radial amplitude function (as it depends on only the radial
coordinate $\rho$). Since we are dealing with optical vortex modes the index
$\ell$ is the familiar winding number, while $m$ could be a radial number as
in Laguerre-Gaussian modes, or it could be redundant as for Bessel modes. We
shall develop the analysis for a general ${\tilde{\cal F}}$, which could be
appropriate for any optical vortex.

##  III cycle-averaged optical properties

The optical vortex mode properties of significance which we shall be concerned
with here are the three cycle-averaged properties, namely the optical spin
angular momentum (SAM) density ${\bar{\bf s}}$, the angular momentum (AM)
density $\bar{j}$ and the helicity density ${\bar{\eta}}$. These cycle-
averaged densities are generally defined as follows.

$\bar{s}$ $\displaystyle=$
$\displaystyle\frac{1}{4\omega}\Im\left\\{[\epsilon_{0}{\bf E}^{*}\times{\bf
E}]+\frac{1}{\mu_{0}}[{\bf B}^{*}\times{\bf B}]\right\\};$ (6)
$\displaystyle=$
$\displaystyle{\mbox{\boldmath$\bar{s}$}}_{E}+{\mbox{\boldmath$\bar{s}$}}_{B}\;\;\;\;\;\;{\rm{(SAM\;density)}}$
${\mbox{\boldmath$\bar{j}$}}={\bf
r}\times{\bar{\mbox{\boldmath$\pi$}}};\;\;\;\;({\rm{AM\;density}})$ (7)
${\bar{\eta}}({\bf r})=-\frac{\epsilon_{0}c}{2\omega}\Im{[{\bf E}^{*}\cdot{\bf
B}]},\;\;\;\;({\rm{Helicity\;density}});$ (8)

where in the above ${\bar{\mbox{\boldmath$\pi$}}}=\frac{1}{c^{2}}{\bar{\bf
w}}$ is the linear momentum density with ${\bar{\bf
w}}=\frac{1}{2\mu_{0}}\Re[{\bf E}^{*}\times{\bf B}]$ the energy density. The
symbols ${\Re}[...]$ and $\Im[...]$ stand for real and imaginary parts of […]
and the superscript * in ${\bf E}^{*}$ stands for the complex conjugate of
${\bf E}$. As stated above we deal in turn with the evaluations of above
densities specifically in relation to the radially-polarised optical vortex
modes. The final tasks involve evaluating the total (integrated) properties,
namely total helicity, total SAM and total angular momentum with each
evaluated as the space integral of the density variations over the x-y plane.

###  III.1 Evaluation of SAM

We begin with the evaluation of the spin angular momentum density. We have for
the electric field part ${\mbox{\boldmath$\bar{s}$}}_{E}$

$\displaystyle{\mbox{\boldmath$\bar{s}$}}_{E}$ $\displaystyle=$
$\displaystyle\frac{\epsilon_{0}}{4\omega}\Im[{\bf E}^{*}\times{\bf E}]$ (9)

Substituting for the electric fields , we have

$\displaystyle{\mbox{\boldmath$\bar{s}$}}_{E}$ $\displaystyle=$
$\displaystyle\frac{\epsilon_{0}}{4\omega}\Im\left(-ick_{z}{\mbox{\boldmath$\hat{\rho}$}}{\cal
F}^{*}e^{-ik_{z}z}-{\mbox{\boldmath$\hat{z}$}}c\left\\{\frac{{\cal
F^{*}}}{\rho}+\left[\frac{\partial{\cal
F}}{\partial\rho}\right]^{*}\right\\}e^{-ik_{z}z}\right)\times\left(ick_{z}{\mbox{\boldmath$\hat{\rho}$}}{\cal
F}e^{ik_{z}z}-{\mbox{\boldmath$\hat{z}$}}c\left\\{\frac{{\cal
F}}{\rho}+\frac{\partial{\cal F}}{\partial\rho}\right\\}e^{ik_{z}z}\right)$
$\displaystyle=$
$\displaystyle-\frac{c^{2}k\epsilon_{0}}{2\omega}\left(\frac{|{\tilde{\cal
F}}|^{2}}{\rho}+{\tilde{\cal F}}{\tilde{\cal
F}}^{\prime}\right){\mbox{\boldmath$\hat{\phi}$}}$ $\displaystyle=$
$\displaystyle-\frac{c^{2}k\epsilon_{0}}{2\omega}\left(\frac{|{\tilde{\cal
F}}|^{2}}{\rho}+{\tilde{\cal F}}{\tilde{\cal
F}}^{\prime}\right)(-\sin\phi{\mbox{\boldmath$\hat{x}$}}+\cos\phi{\mbox{\boldmath$\hat{y}$}})$
(11)

Thus we have found that the electric field part of the SAM density of a
radially-polarised optical vortex mode is oriented azimuthally, so has both
$x-$ and $y-$ components which vary with $\phi$.

Consider next the magnetic field contribution
${\mbox{\boldmath$\bar{s}$}}_{B}=\frac{1}{4\omega\mu_{0}}\Im[{\bf
B}^{*}\times{\bf B}]$. We have

$\displaystyle{\mbox{\boldmath$\bar{s}$}}_{B}$ $\displaystyle=$
$\displaystyle\frac{1}{4\mu_{0}\omega}\Im\left\\{\left(-ik_{z}{\mbox{\boldmath$\hat{\phi}$}}{\cal
F}e^{-ik_{z}z}-{\mbox{\boldmath$\hat{z}$}}\frac{1}{\rho}\left[\frac{\partial{\cal
F}}{\partial\phi}\right]^{*}e^{-ik_{z}z}\right)\times\left(ik_{z}{\mbox{\boldmath$\hat{\phi}$}}{\cal
F}e^{ik_{z}z}-{\mbox{\boldmath$\hat{z}$}}\frac{1}{\rho}\left[\frac{\partial{\cal
F}}{\partial\phi}\right]e^{ik_{z}z}\right)\right\\}$ (12) $\displaystyle=$
$\displaystyle\Im\left(\frac{ik_{z}}{4\mu_{0}\omega\rho}\right)\left[{\cal
F}^{*}\left(\frac{\partial{\cal F}}{\partial\phi}\right)+{\cal
F}\left(\frac{\partial{\cal
F}}{\partial\phi}\right)^{*}\right]{\mbox{\boldmath$\hat{\phi}$}}\times{\mbox{\boldmath$\hat{z}$}}={\bf
0}$

The last equality follows once we evaluate the angular derivatives. Thus the
magnetic field contribution to the SAM density is zero and we have only the
electric field contribution ${\mbox{\boldmath$\bar{s}$}}_{E}$, given by
Eq.(11).

The total (space-integrated) SAM vanishes identically

${\mbox{\boldmath$\bar{S}$}}=\int_{0}^{2\pi}d\phi\int_{0}^{\infty}\rho\;d\rho\;[{\mbox{\boldmath$\bar{s}$}}_{E}+{\mbox{\boldmath$\bar{s}$}}_{B}]={\bf
0};$ (13)

where each of the $x-$ and $y-$ components of the SAM density yields a zero
result by virtue of the angular integration.

###  III.2 Evaluation of angular momentum

Next we evaluate the cycle-averaged angular momentum density of the radially-
polarised vortex which is given by

${\mbox{\boldmath$\bar{j}$}}={\bf
r}\times{\bar{\mbox{\boldmath$\pi$}}}=\frac{1}{2c^{2}\mu_{0}}\left[{\bf
r}\times\left(\Re[{\bf E}^{*}\times{\bf B}\right)\right]$ (14)

Consider the cross product in Eq.(14) which we evaluate as follows

$\displaystyle{\bf E}^{*}\times{\bf B}$ $\displaystyle=$
$\displaystyle\left(-ick_{z}{\mbox{\boldmath$\hat{\rho}$}}{\cal
F}^{*}e^{-ik_{z}z}-{\mbox{\boldmath$\hat{z}$}}c\left\\{\frac{{\cal
F^{*}}}{\rho}+\left[\frac{\partial{\cal
F}}{\partial\rho}\right]^{*}\right\\}e^{-ik_{z}z}\right)\times\left(ik_{z}{\mbox{\boldmath$\hat{\phi}$}}{\cal
F}e^{ik_{z}z}-e^{ik_{z}z}{\mbox{\boldmath$\hat{z}$}}\frac{\partial{\cal
F}}{\rho\partial\phi}\right)$ (15) $\displaystyle=$ $\displaystyle
ick_{z}\left(\frac{|{\tilde{\cal F}}|^{2}}{\rho}+{\tilde{\cal F}}{\tilde{\cal
F}}^{\prime}\right){\mbox{\boldmath$\hat{\rho}$}}+\left\\{\ell
ck_{z}\frac{|{\tilde{\cal
F}}|^{2}}{\rho}\right\\}{\mbox{\boldmath$\hat{\phi}$}}+k_{z}^{2}c|{\tilde{\cal
F}}|^{2}{\mbox{\boldmath$\hat{z}$}}$

Therefore on taking the real part, we have

$\Re[{\bf E}^{*}\times{\bf B}]=k_{z}c\ell\frac{|{\tilde{\cal
F}}|^{2}}{\rho}{\mbox{\boldmath$\hat{\phi}$}}+k_{z}^{2}c|{\tilde{\cal
F}}|^{2}{\mbox{\boldmath$\hat{z}$}}$ (16)

The angular momentum density follows

$\bar{j}$ $\displaystyle=$ $\displaystyle\frac{1}{2c^{2}\mu_{0}}{\bf
r}\times\Re[{\bf E}^{*}\times{\bf B}]$ $\displaystyle=$
$\displaystyle\left(\frac{k_{z}}{2c\mu_{0}}\right)\\{\rho\mbox{\boldmath${\hat{\rho}}$}\\}\times\left\\{\ell\frac{|{\tilde{\cal
F}}|^{2}}{\rho}{\mbox{\boldmath$\hat{\phi}$}}+k_{z}{|{\tilde{\cal
F}}|^{2}}{\mbox{\boldmath$\hat{z}$}}\right\\}$ $\displaystyle=$
$\displaystyle\left(\frac{k_{z}}{2c\mu_{0}}\right)\left\\{\ell{|{\tilde{\cal
F}}|^{2}}{\mbox{\boldmath$\hat{z}$}}-k_{z}|{\tilde{\cal
F}}|^{2}\rho{\mbox{\boldmath$\hat{\phi}$}}\right\\}$

Since we have
${\mbox{\boldmath$\hat{\phi}$}}=-{\mbox{\boldmath$\hat{x}$}}\sin\phi+{\mbox{\boldmath$\hat{y}$}}\cos\phi$
the angular momentum density vector has all three Cartesian components.
However, the transverse ($x-$ and $y-$) components are $\phi-$ dependent and,
as we point out shortly, will result in zero on angular integration.

Once again we consider the total angular momentum as the space integral of the
angular momentum density.

${\mbox{\boldmath$\bar{\cal
J}$}}=\int_{0}^{2\pi}d\phi\int_{0}^{\infty}\rho\;d\rho\;{\mbox{\boldmath$\bar{j}$}}$
(18)

The $x-$ and $y-$ components give zero each due to vanishing angular
integration. We are left only with the z-component, so we have

${\mbox{\boldmath$\bar{\cal
J}$}}={\mbox{\boldmath$\hat{z}$}}\ell\left(\frac{k_{z}\pi}{c\mu_{0}}\right)I_{P}$
(19)

where the integral $I_{P}$ is related to the applied power $\cal P$ of the
mode, evaluated as the space integral over the beam cross-section of the
z-component of the Poynting vector. We have

${\cal P}=\frac{1}{2\mu_{0}}\int_{0}^{2\pi}d\phi\int_{0}^{\infty}|({\bf
E}^{*}\times{\bf B})_{z}|\rho d\rho$ (20)

with

$({\bf E}^{*}\times{\bf B})_{z}=ck_{z}^{2}|{\tilde{\cal F}}|^{2}$ (21)

We can then write for $I_{P}$

$I_{P}=\int_{0}^{\infty}|{\tilde{\cal F}}|^{2}\rho\;d\rho$ (22)

Thus we obtain for the power ${\cal P}$

${\cal P}=\left(\frac{\pi ck_{z}^{2}}{\mu_{0}}\right)I_{P}$ (23)

Substituting for $I_{P}$, we have for the total angular momentum per unit
length

${\mbox{\boldmath$\bar{\cal J}$}}=\ell{\cal L}_{0}{\mbox{\boldmath$\hat{z}$}}$
(24)

where ${\cal L}_{0}$ has the dimensions of angular momentum per unit length
and is given by

${\cal L}_{0}=\left(\frac{{\cal P}}{k_{z}c^{2}}\right)$ (25)

Thus we find that $\bar{\cal J}$ is axial and proportional to $\ell$ which
confirms that the angular momentum carried by the radially-polarised LG mode
is purely an orbital angular momentum. Note that we have determined the
angular momentum without specifying the type of mode. The result is therefore
general and it agrees with the result for linearly-polarised Laguerre-Gaussian
light Koksal  _et al._ (2022a). There is no contribution to be associated with
spin angular momentum, which we have already confirmed to be zero.

###  III.3 Evaluation of helicity and chirality

In optical physics, the properties of helicity and chirality are such that two
optical modes which differ only in the sign of the winding number $\ell$ are
distinguishable. In that case, one beam is a phase-inverted mirror image of
the other beam. Such kind of optical mode is then said to exhibit chirality
which has been a subject of considerable investigation Rañada (1989); Rañada
and Trueba (1995); Ranada and Trueba (1996); Tang and Cohen (2010); Cameron
_et al._ (2012); Fernandez-Corbaton  _et al._ (2012); Afanasiev and
Stepanovsky (1996). More recent investigations have spurred interest in
chirality which include the references Bliokh  _et al._ (2017); Crimin  _et
al._ (2019a, b); Woźniak  _et al._ (2019); Nechayev  _et al._ (2021);
Lembessis  _et al._ (2021); Forbes and Jones (2021).

In their recent work, Nechayev et al Nechayev  _et al._ (2021) suggested that
there exist two kinds of chirality. The first is non-geometrical which
includes optical chirality due to elliptical polarisation and the second,
which is termed the Kelvin chirality, depends on the geometrical structure of
the mode. This work followed the experiment by Wozniak et al Woźniak  _et al._
(2019) in which linearly-polarised Laguerre-Gaussian modes were shown to
display chirality. Here we seek to evaluate the chirality density and its
space integral and aim to find out whether the chirality of radially-polarised
optical modes can be classed as Kelvin chirality, or is it the usual optical
non-geometrical chirality.

We begin by considering helicity, rather than chirality, as the two differ by
a proportionality factor Koksal  _et al._ (2022b). The cycle-averaged helicity
density of the radially-polarised mode is as defined generally in Eq.(8).
Substituting for the fields using Eqs.(2) and (4), we have for the dot product

$\displaystyle[{\bf E}^{*}\cdot{\bf B}]$ $\displaystyle=$
$\displaystyle\left\\{-ick_{z}{\mbox{\boldmath$\hat{\rho}$}}{\cal
F}e^{-ik_{z}z}-{\mbox{\boldmath$\hat{z}$}}c\frac{1}{\rho}\left(\frac{\partial(\rho{\cal
F})}{\partial\rho}\right)^{*}e^{-ik_{z}z}\right\\}\cdot\left\\{ik_{z}{\mbox{\boldmath$\hat{\phi}$}}{\cal
F}e^{ik_{z}z}-{\mbox{\boldmath$\hat{z}$}}\frac{1}{\rho}\frac{\partial{\cal
F}}{\partial\phi}e^{ik_{z}z}\right\\}$ (26) $\displaystyle=$ $\displaystyle
c\frac{1}{\rho^{2}}\left(\frac{\partial(\rho{\cal
F})}{\partial\rho}\right)^{*}\left(\frac{\partial{\cal
F}}{\partial\phi}\right)$ $\displaystyle=$
$\displaystyle\frac{ic\ell}{\rho^{2}}\left(\frac{\partial(\rho{\tilde{\cal
F})}}{\partial\rho}\right){\tilde{\cal F}}=i\ell
c\left[\frac{1}{\rho}{\tilde{\cal F}}^{\prime}{\tilde{\cal
F}}+\frac{1}{\rho^{2}}|{\tilde{\cal F}}|^{2}\right]$

where ${\tilde{\cal F}^{\prime}}=d{\tilde{\cal F}}/d\rho$. We therefore have
for the helicity density of the radially-polarised general vortex mode

${\bar{\eta}}({\bf
r})=-\ell\frac{\epsilon_{0}c^{2}}{2\omega}\left[\frac{1}{\rho}{\tilde{\cal
F}^{\prime}}{\tilde{\cal F}}+\frac{1}{\rho^{2}}|{\tilde{\cal F}}|^{2}\right]$
(27)

This result is applicable to any paraxial radially-polarised optical vortex
and, if required for a particular case, all we need then is to specify the
amplitude function ${\tilde{\cal F}}$. Note that the helicity density is
proportional to $\ell$ and so in addition to changing with the magnitude of
$\ell$, it also changes with the sign of $\ell$, hence exhibiting the
chirality feature that is common to all optical vortex modes.

###  III.4 Integrated Helicity

The total integral of the helicity density over the $x-y$ plane is

$\displaystyle{\cal{\bar{C}}}_{\ell,p}$ $\displaystyle=$
$\displaystyle-\ell\frac{\pi\epsilon_{0}c^{2}}{\omega}\int_{0}^{\infty}\rho
d\rho\left[\frac{1}{\rho}{\tilde{\cal F}}_{p\ell}^{\prime}{\tilde{\cal
F}}_{p\ell}+\frac{1}{\rho^{2}}|{\tilde{\cal F}}_{p\ell}|^{2}\right]$ (28)
$\displaystyle=$ $\displaystyle{\cal I}_{1}+{\cal I}_{2}$

We now show that ${\cal I}_{1}$ is identically zero for all ${\tilde{\cal
F}}$. We have

$\displaystyle{\cal I}_{1}$ $\displaystyle=$
$\displaystyle-\ell\frac{\epsilon_{0}\pi
c^{2}}{2\omega}\int_{0}^{\infty}\left(\frac{d}{d\rho}{\tilde{\cal
F}}^{2}\right)\;d\rho$ (29) $\displaystyle=$
$\displaystyle-\ell\frac{\epsilon_{0}\pi c^{2}}{2\omega}\left[{\tilde{\cal
F}}^{2}(\rho)\right]_{0}^{\infty}$ $\displaystyle=$ $\displaystyle 0$

since the ${\tilde{\cal F}}(0)=0={\tilde{\cal F}}(\infty)$. We are thus left
with the second term so that for any ${\tilde{\cal F}}$ the helicity per unit
length is given by

${\cal{\bar{C}}}_{\ell,p}={\cal
I}_{2}=-\ell\frac{\pi\epsilon_{0}c^{2}}{\omega}\int_{0}^{\infty}\rho
d\rho\left[\frac{1}{\rho^{2}}|{\tilde{\cal F}}_{\ell,p}|^{2}\right]$ (30)

This is the general expression for the total helicity per unit length of a
mode of any paraxial twisted light mode that is radially-polarised. Although
we are able to make definitive statements about the total SAM and total
angular momentum for all forms of ${\cal F}$ characterising radially-polarised
twisted light, we are not able to proceed further to evaluate the total
helicity in Eq.(30) without specific knowledge about the form of ${\tilde{\cal
F}}$. In the next section we focus on the form of ${\tilde{\cal F}}$
appropriate for a Laguerre-Gaussian optical vortex and aim to evaluate the
helicity density given by Eq.(27) and the total helicity given by Eq.(30).

##  IV Applications to Laguerre-Gaussian modes

A paraxial Laguerre-Gaussian mode of winding number $\ell$, radial number $p$
and waist $w_{0}$ has an amplitude function given by

${\tilde{\cal F}}_{\ell,p}(\rho)={\cal
E}_{0}\sqrt{\frac{p!}{(p+|\ell|)!}}e^{-\frac{\rho^{2}}{w_{0}^{2}}}\left(\frac{\sqrt{2}\rho}{w_{0}}\right)^{|\ell|}L^{|\ell|}_{p}\left(\frac{2\rho^{2}}{w_{0}^{2}}\right)$
(31)

where ${\cal E}_{0}$ is a normalisation factor and we have identified $m$ as
the radial number $p$ in LG modes. The factor ${\cal E}_{0}$ is determined in
terms of the applied power $\cal P$ of the mode which we have already
evaluated above in terms of the integral $I_{P}$. For a Laguerre-Gaussian mode
the integral in (22) is standard and gives

$I_{P}=\int_{0}^{\infty}|{\tilde{\cal F}}|^{2}\rho\;d\rho=\frac{1}{4}{\cal
E}_{0}^{2}w_{0}^{2}$ (32)

and so we can now determine the overall factor ${\cal E}_{0}$ for the
Laguerre-Gaussian beam. We have

${\cal E}_{0}^{2}=\frac{4{\cal P}}{\pi\epsilon_{0}c^{3}k_{z}^{2}w_{0}^{2}}$
(33)

Consider first the variations of the helicity density for representative
modes, namely $\ell=1,2$. The general expression for the helicity density
Eq.(27) is circularly-symmetric as it is a function only of the radial
coordinate $\rho$. For a LG mode, we simply substitute for ${\tilde{\cal F}}$
given by Eq.(31) and we obtain the variation of the helicity density, as shown
in Fig. 1 for $\ell=1$ and Fig. 2 for $\ell=2$. The plots in each case show
the contributions from the first term, the second term and their sum.
Confirmation of the vanishing integral of the first term are shown in Figs. 3
for $\ell=1$ and Fig.4 for $\ell=2$ where the areas under the curves
corresponding to the areas enclosed by the integrands due to the first term
are zero for both $\ell=1$ and $\ell=2$. We have confirmed that the variations
shown for the case $\ell=2$ define the trend for $\ell>2$

Note, in particular, that the variations of the helicity density for the case
$\ell=1$, shown in Fig. 1 differ significantly from those of $\ell\geq 2$ in
Fig.2, primarily in that the helicity density does not vanish at the core
where $\rho=0$ for $\ell=1$, while it does vanish at $\rho=0$ for $\ell\geq
2$. This behaviour can be explained by inspecting the general form of the
helicity density Eq.(27). When applied to the Laguerre-Gaussian ${\cal F}$ for
$\ell=1$, we have from Eq.(31)

${\tilde{\cal F}}_{\ell=1}\propto\rho
e^{-\rho^{2}/w_{0}^{2}}L^{1}_{p}\left(\frac{2\rho^{2}}{w_{0}^{2}}\right)$ (34)

Hence

$|{\tilde{\cal
F}}_{\ell=1}|^{2}\propto\rho^{2}e^{-2\rho^{2}/w_{0}^{2}}\left[L^{1}_{p}\left(\frac{2\rho^{2}}{w_{0}^{2}}\right)\right]^{2}$
(35)

Also, since ${\cal F}^{\prime}$ does not vanish at $\rho=0$ we can write

$[{\tilde{\cal F}}^{\prime}{\tilde{\cal F}}]_{\ell=1}\propto\rho
e^{-\rho^{2}/w_{0}^{2}}L^{1}_{p}\left(\frac{2\rho^{2}}{w_{0}^{2}}\right){\tilde{\cal
F}}^{\prime}_{\ell=1}$ (36)

When substituted in the helicity density expression Eq.(27) we see that the
$1/\rho$ in the first term cancels with the factor $\rho$ in the numerator.
Similarly the factor $1/\rho^{2}$ in the second term cancels the factor
$\rho^{2}$ in the numerator of the second term. The overall variation amounts
to a non zero value of the helicity at $\rho=0$ only in the case $\ell=1$.
This variation contrasts with the case $\ell\geq 2$ in which the numerators in
the two terms have higher powers of $\rho$, always guaranteeing that the
helicity density vanishes at $\rho=0$.

Figure 1: Variation with the radial coordinate $\rho$ (in units of the beam
waist $w_{0}$) of the helicity density due to LG modes for which (a) $\ell=1$.
The dashed red curve shows the contribution of the first term in Eq.(27)
involving the derivative ${\tilde{\cal F}^{\prime}}$. The blue dotted curve
represents the contribution of the second term and the solid black curve is
the sum. Note in particular that for $\ell=1$ the helicity density does not
vanish at the core $\rho=0$

.

Figure 2: The case of LG mode with $\ell=2$. Variation with the radial
coordinate $\rho$ (in units of the beam waist $w_{0}$) of the helicity density
due to LG modes for which (a) $\ell=2$. The dashed red curve shows the
contribution of the first term in Eq.(27) involving the derivative
${\tilde{\cal F}^{\prime}}$. The blue dotted curve represents the contribution
of the second term and the solid black curve is the sum. Note in particular
that for $\ell=2$ the helicity density vanishes at the core $\rho=0$

.

Figure 3: The case of LG mode with $\ell=1$. Variations of the helicity
density integrand terms in ${\cal I}_{1}$ and ${\cal I}_{2}$, (defined in
Eq.(28)) with the radial coordinate (in units of $w_{0}$). The area enclosed
by each curve corresponds to the contribution of the term to the total
(integrated) helicity. The area enclosed by the dashed red curve is verified
to be zero, consistent with Eq.(29). Figure 4: The case of LG mode with
$\ell=2$. Variations of the helicity density integrand terms in ${\cal I}_{1}$
and ${\cal I}_{2}$, (defined in Eq.(28)) with the radial coordinate (in units
of $w_{0}$). The area enclosed by each curve corresponds to the contribution
of the term to the total (integrated) helicity. The area enclosed by the
dashed red curve is verified to be zero, consistent with Eq.(29).

It is straightforward to proceed to evaluate the helicity per unit length for
${\tilde{\cal F}}$ corresponding to a Laguerre-Gaussian mode. Using the
integration variable $x=2\rho^{2}/w_{0}^{2}$ we have

$\displaystyle{\cal{\bar{C}}}_{\ell,p}$ $\displaystyle=$
$\displaystyle-\ell\frac{\pi\epsilon_{0}c^{2}}{\omega}{\cal
E}_{0}^{2}\frac{p!}{2(p+|\ell|)!}\int_{0}^{\infty}x^{|\ell|-1}e^{-x}[L^{|\ell|}_{p}(x)]^{2}dx$
(37) $\displaystyle=$
$\displaystyle-\ell\frac{\pi\epsilon_{0}c^{2}}{\omega}{\cal
E}_{0}^{2}\frac{1}{2|\ell|}$

The details of the evaluation of the integral in Eq.(37) are shown in Appendix
A. We now have

$\displaystyle{\cal{\bar{C}}}_{\ell,p}$ $\displaystyle=$
$\displaystyle-\frac{{\cal
E}_{0}^{2}\pi\epsilon_{0}c^{2}}{2\omega}\frac{\ell}{|\ell|}$ (38)
$\displaystyle=$
$\displaystyle-\left(\frac{\ell}{|\ell|}\right)\left(\frac{{\cal
P}}{k_{z}c^{2}}\right)\left\\{\frac{2}{k_{z}^{2}w_{0}^{2}}\right\\}$
$\displaystyle=$ $\displaystyle\pm{\cal
L}_{0}\left\\{\frac{2}{k_{z}^{2}w_{0}^{2}}\right\\}$

where ${\cal L}_{0}$ is the constant angular momentum per unit length as
defined in Eq.(25) for a fixed power ${\cal P}$ and we have substituted for
${\cal E}_{0}$ using Eq.(33). Note that this result is independent of $p$ and
clearly depends only on the sign (not the magnitude) of $\ell$. We know that
the mode is not circularly polarised, but we have found that the total
helicity is similar to, but not the same as, that of circular polarisation,
characterised by the pre-factor $\sigma=\pm 1$. It is easy to check that the
helicity has the dimensions of angular momentum per unit length, but the
factor $1/k_{z}^{2}w_{0}^{2}$ is small for $w_{0}^{2}>>1/k_{z}^{2}$ (which
amounts to $w_{0}>>{\bar{\lambda}}$ where ${\bar{\lambda}}=\lambda/2\pi$ is a
reduced wavelength). Thus the helicity is significant in the case of LG beams
only for small beam waits $w_{0}$ and diminishes for progressively larger
$w_{0}$.

##  V Comments and Conclusions

Our primary aim in this paper involved the derivation of the optical
properties of paraxial radially-polarised twisted light modes, namely the spin
orbital angular momentum (SAM), the total angular momentum (which is for mally
the sum of spin and orbital angular momentum for paraxial light), and their
helicity and chirality. We set out to keep the type of mode unspecified and
arrived at the results for the angular momentum and SAM, but for the helicity
we are able to arrive at general results and applied them to the Laguerre-
Gaussian modes as a specific case. Our treatment is based on general
expressions describing the radially-polarised electric and magnetic fields
which incorporate the longitudinal component and which were subject to
verification of Maxwell’s consistency conditions, namely that the electric
field in Cartesian coordinates follows from a derived expression of the
magnetic field using the Maxwell curl equation, also in Cartesian coordinates,
and the magnetic field follows from the electric field using the other Maxwell
curl equation. The final Cartesian expressions are then presented in
cylindrical polar coordinates. We have checked that the formalism presented is
verifiable to leading order in the paraxial approximation.

We have found that in general, the radially-polarised twisted light modes
exhibit only cycle-averaged transverse SAM components, which arise entirely
from the electric field part, are $\phi$-dependent, while the magnetic field
contribution to the SAM density is shown to be identically zero. The space
integral of the SAM density leading to the total SAM is therefore zero. Thus
we have confirmed that in general such radially-polarised modes have no SAM.
Next we evaluated the angular momentum, which, for paraxial modes, is always
the sum of the spin angular momentum and the orbital angular momentum and we
have found that the angular momentum density is proportional to the winding
number $\ell$, indicating that this contribution of the total angular momentum
density is purely orbital, but there are also transverse density components
which, like the SAM density case, depend on $\phi$ and so lead to zero on
spatial integration. The integrated total angular momentum is also
proportional to $\ell$ and so purely orbital in origin. This result is
consistent with the SAM result which was evaluated independently, that the
mode has no SAM.

The helicity density evaluations could also be carried out for a general
radially-polarised optical vortex mode characterised by an unspecified
${\tilde{\cal F}}$ with the density and the general results for the helicity
density and the total (integrated) helicity we arrived at are shown in Eq.(27)
for the helicity density and in Eq.(30) for the integrated helicity. We
proceeded to explore the helicity density variations for the spacial cases of
Laguerre-Gaussian modes and pointed out the special behaviour for $\ell=1$ in
that the helicity density does not vanish at $\ell=1$ at the core $\rho=0$,
while it does vanish for all $\ell\geq 2$. We explained this behaviour by
inspecting the $\rho$ variations of the two terms in the helicity density,
confirming that for $\ell=1$ the helicity density has overall dependence
$\rho^{0}$ and so the helicity density does not vanish at $\rho=0$ for
$\ell=1$.

Finally, we evaluated the total helicity of the radially-polarised Laguerre-
Gaussian optical vortex and found it equal to $(\ell/|\ell|){\cal Q}={\cal
N}{\cal Q}$, where ${\cal N}=\pm 1$ is interpreted as a Hopf index and $Q$ is
the action constant. The $\pm 1$ is reminiscent of $\sigma=\pm 1$ for circular
polarisation. This result holds for any radially-polarised LG mode however
large the magnitude of its winding number $\ell$ is, but the action constant
${\cal Q}$ and consequently the helicity are significant only for small beam
waist $w_{0}$ and diminish for all such LG modes of large $w_{0}$. Since this
type of helicity (and chirality) originate from the spatial structure of the
radially-polarised mode, they can be categorised as of the Kelvin type
Nechayev  _et al._ (2021).

###  V.1 Disclosures

The authors declare no conflicts of interest.

## Appendix A

We evaluate the integral in Eq.(37) which is

${\cal I}=\int_{0}^{\infty}x^{|\ell|-1}e^{-x}[L^{|\ell|}_{p}(x)]^{2}dx$ (39)

Consider the evaluation of following derivative of the function
$x^{|\ell|}e^{-x}[L^{|\ell|}_{p}(x)]^{2}$ with respect to $x$

$\frac{d}{dx}\\{x^{|\ell|}e^{-x}[L^{|\ell|}_{p}(x)]^{2}\\}=|\ell|x^{|\ell|-1}e^{-x}[L^{|\ell|}_{p}(x)]^{2}-x^{|\ell|}e^{-x}[L^{|\ell|}_{p}(x)]^{2}-2x^{|\ell|}e^{-x}L^{|\ell|}_{p}(x)L^{|\ell|+1}_{p-1}(x)$
(40)

So we may now write

$\displaystyle{\cal I}$ $\displaystyle=$
$\displaystyle\int_{0}^{\infty}x^{|\ell|-1}e^{-x}[L^{|\ell|}_{p}(x)]^{2}dx$
(41) $\displaystyle=$
$\displaystyle\frac{1}{|\ell|}\int_{0}^{\infty}\left(\frac{d}{dx}\\{x^{|\ell|}e^{-x}[L^{|\ell|}_{p}(x)]^{2}\\}+x^{|\ell|}e^{-x}[L^{|\ell|}_{p}(x)]^{2}+2x^{|\ell|}e^{-x}L^{|\ell|}_{p}(x)L^{|\ell|+1}_{p-1}(x)\right)dx$
$\displaystyle=$
$\displaystyle\frac{1}{|\ell|}\left([x^{|\ell|}e^{-x}(L^{|\ell|}_{p}(x))^{2}]_{0}^{\infty}+\int_{0}^{\infty}x^{|\ell|}e^{-x}[L^{|\ell|}_{p}(x)]^{2}dx+2\int_{0}^{\infty}x^{|\ell|}e^{-x}L^{|\ell|}_{p}(x)L^{|\ell+1|}_{p-1}(x)dx\right)$

where we have integrated by parts in the derivative term. The first term in
the last equality is zero on applying the integration limits, while the third
term is the standard orthogonality integral of the associated Laguerre
functions and so vanishes as well. We are thus left left with

${\cal
I}=\frac{1}{|\ell|}\int_{0}^{\infty}x^{|\ell|}e^{-x}[L^{|\ell|}_{p}(x)]^{2}dx=\frac{1}{|\ell|}\frac{(p+|\ell|)!}{p!}$
(42)

## References

  * Allen  _et al._ (1999) L. Allen, M. Padgett, and M. Babiker, Progress in optics 39, 291 (1999).
  * Allen  _et al._ (2003) L. Allen, S. M. Barnett, and M. J. Padgett,  _Optical Angular Momentum_ (IOP Publishing, 2003).
  * Zhan (2004) Q. Zhan, Opt. Express 12, 3377 (2004).
  * Franke-Arnold  _et al._ (2008) S. Franke-Arnold, L. Allen, and M. Padgett, Laser & Photonics Reviews 2 (2008).
  * Andrews and Babiker (2012) D. L. Andrews and M. Babiker, eds.,  _The Angular Momentum of Light_ (Cambridge University Press, Cambridge, 2012).
  * Torres and L. Torner (2011) J. P. Torres and E. L. Torner,  _Twisted photons_ (Wiley-VTH, 2011).
  * Babiker  _et al._ (2018) M. Babiker, D. L. Andrews, and V. E. Lembessis, Journal of Optics 21, 013001 (2018).
  * Holmes and Galvez (2019) B. M. Holmes and E. J. Galvez, Journal of Optics 21, 104001 (2019).
  * Mushiake  _et al._ (1972) Y. Mushiake, K. Matsumura, and N. Nakajima, Proceedings of the IEEE 60, 1107 (1972).
  * Kozawa and Sato (2005) Y. Kozawa and S. Sato, Opt. Lett. 30, 3063 (2005).
  * Zhan (2009) Q. Zhan, Adv. Opt. Photon. 1, 1 (2009).
  * Kozawa and Sato (2007) Y. Kozawa and S. Sato, J. Opt. Soc. Am. A 24, 1793 (2007).
  * Dorn  _et al._ (2003) R. Dorn, S. Quabis, and G. Leuchs, Physical Review Letters 91, 1 (2003).
  * Bashkansky  _et al._ (2010) M. Bashkansky, D. Park, and F. K. Fatemi, Opt. Express 18, 212 (2010).
  * Levy  _et al._ (ts) U. Levy, Y. Silberberg, and N. Davidson, Adv. Opt. Photon. 11, 828 (ts).
  * Nechayev  _et al._ (2021) S. Nechayev, J. S. Eismann, R. Alaee, E. Karimi, R. W. Boyd, and P. Banzer, Phys. Rev. A 103, L031501 (2021).
  * Koksal  _et al._ (2022a) K. Koksal, M. Babiker, V. E. Lembessis, and J. Yuan, J. Opt. Soc. Am. B 39, 459 (2022a).
  * Rañada (1989) A. F. Rañada, letters in mathematical physics 18, 97 (1989).
  * Rañada and Trueba (1995) A. F. Rañada and J. Trueba, Physics Letters A 202, 337 (1995).
  * Ranada and Trueba (1996) A. Ranada and J. Trueba, Nature 384, 124 (1996).
  * Tang and Cohen (2010) Y. Tang and A. E. Cohen, Phys. Rev. Lett. 104, 163901 (2010).
  * Cameron  _et al._ (2012) R. P. Cameron, S. M. Barnett, and A. M. Yao, New Journal of Physics 14, 053050 (2012).
  * Fernandez-Corbaton  _et al._ (2012) I. Fernandez-Corbaton, X. Zambrana-Puyalto, and G. Molina-Terriza, Physical Review A 86, 042103 (2012).
  * Afanasiev and Stepanovsky (1996) G. Afanasiev and Y. P. Stepanovsky, Il Nuovo Cimento A (1965-1970) 109, 271 (1996).
  * Bliokh  _et al._ (2017) K. Y. Bliokh, A. Y. Bekshaev, and F. Nori, Phys. Rev. Lett. 119, 073901 (2017).
  * Crimin  _et al._ (2019a) F. Crimin, N. Mackinnon, J. Götte, and S. Barnett, Applied Sciences 9, 828 (2019a).
  * Crimin  _et al._ (2019b) F. Crimin, N. Mackinnon, J. Götte, and S. Barnett, Journal of Optics 21, 094003 (2019b).
  * Woźniak  _et al._ (2019) P. Woźniak, I. De Leon, K. Höflich, G. Leuchs, and P. Banzer, Optica 6, 961 (2019).
  * Lembessis  _et al._ (2021) V. E. Lembessis, K. Koksal, J. Yuan, and M. Babiker, Phys. Rev. A 103, 013106 (2021).
  * Forbes and Jones (2021) K. A. Forbes and G. A. Jones, Journal of Optics 23, 115401 (2021).
  * Koksal  _et al._ (2022b) K. Koksal, M. Babiker, V. E. Lembessis, and J. Yuan, J. Opt. Soc. Am. B 39, 459 (2022b).


# Szilard engines and information-based work extraction for active systems

Paolo Malgaretti<EMAIL_ADDRESS>Helmholtz Institute Erlangen-
Nürnberg for Renewable Energy (IEK-11), Forschungszentrum Jülich, Cauerstr. 1,
91058 Erlangen, Germany  Holger Stark  Institut für Theoretische Physik,
Technische Universität Berlin, Hardenbergstr. 36, 10623 Berlin, Germany

###### Abstract

The out of equilibrium nature of active systems can be exploited for the
design of information-based engines. We design two types of an active Szilard
engine that use a Maxwell demon to extract work from an active bath composed
of non-interacting Active Brownian Particles (ABPs). The two engines exploit
either the quasi-static active pressure of ABPs or the long correlation time
of their velocities. For both engines the active bath allows to overcome the
Landauer principle and to extract larger work compared to conventional Szilard
engines operating in quasi thermal equilibrium. For both of our engines,we
identify the optimal regimes at which the work extracted and the efficiency
are maximized. Finally, we discuss them in the context of synthetic and
biological active systems.

Since Maxwell proposed his demon as an attempt to circumvent the second law of
thermodynamics tai (1867); Maxwell (1871), several attempts have been put
forward to rectify thermal fluctuations Feynman (1963); Reimann (2002). An
interesting twist to the idea of the demon has been provided by Szilard
Szilard (1929); Leff and Rex (2002). He proposed an engine that can extract
work from a single thermal reservoir using information obtained from the
system when measuring some observable Brillouin (1960); Parrondo  _et al._
(2015). Recently, several works have addressed the realization of Szilard
engines in experiments and theory spanning from the single electron Koski  _et
al._ (2014) and boson systems Bengtsson  _et al._ (2018); Aydiner (2021) up to
the macromolecular Ribezzi-Crivellari and Ritort (2019) and colloidal Paneru
and Pak (2020) scale.

Up to now the “fluctuations” in such Szilard engines have been regarded as
thermal. This imposes constraints on the dynamics of the bath, which, among
others, limit the extracted work by the Landauer principle. However, in many
situations the fluctuations are not of thermal origin. Rather, they are
induced by the motion of active agents, such as swimming bacteria, a school of
fish, or a swarm of drones Marchetti  _et al._ (2013); Ebbens and Howse
(2010); Bechinger  _et al._ (2016); Zöttl and Stark (2016); Doostmohammadi
_et al._ (2018); Gompper  _et al._ (2020). This implies that Landauer’s
principle does not apply anymore raising the question how well Szilard engines
perform in contact with such non-thermal baths. The striking feature of active
agents is their persistent motion, i.e., even within the overdamped regime
their (active) velocity has a finite correlation time Bechinger  _et al._
(2016); Zöttl and Stark (2016), in contrast to the delta-correlated velocities
of their passive counterparts. This is captured by the Active Brownian
Particle (ABP) model Romanczuk  _et al._ (2012) characterized by a velocity
correlation time $\tau_{M}$ and speed $v_{0}$ Zöttl and Stark (2016), which
gives a persistence length $\lambda=v_{0}\tau_{M}$. For $t\gg\tau_{M}$ the
motion of ABPs is diffusive with effective diffusion constant
$D_{eff}=D+v_{0}^{2}\tau_{M}/d$ in $d$ dimensions Zöttl and Stark (2016),
whereas for $t\lesssim\tau_{M}$ their motion is persistent and it differs from
their equilibrium counterparts. Despite its simplicity, such a simple model
has been exploited to explain many experimentally observed phenomena spanning
from wall accumulationRotschild (1963); Elgeti  _et al._ (2015) to Motility
Induced Phase Separation (MIPS)Cates and Tailleur (2013, 2015). Indeed, recent
works have addressed the role of a bath of active entities on the performance
of heat-like Zheng  _et al._ (2011); Pietzonka  _et al._ (2019); Holubec  _et
al._ (2020); Malgaretti  _et al._ (2021); Fodor and Cates (2021); Gronchi and
Puglisi (2021); Speck (2022); Ion and Urna (2022) and ratchet-like engines
Sokolov  _et al._ (2010); Di Leonardo  _et al._ (2010); Malgaretti  _et al._
(2013, 2014); Yang and Ripoll (2014); Michelin  _et al._ (2015); Malgaretti
and Stark (2017); Holubec  _et al._ (2017).

In this letter we show that the persistent motion of active agents, modeled as
ABPs, enables two generic realizations of an active Szilard engine. First, in
the traditional quasi-static regime, the performance of the Szilard engine in
contact with an active bath differs from its passive counterpart due to the
active pressure Malgaretti  _et al._ (2021). Second, we suggest a dynamic
Szilard engine that does not operate in the quasi-static regime but exploits
the finite correlation time of the active velocity. In order to keep the model
as general as possible, we rely on the minimal model of an ABP to capture the
effect of persistent motion. Indeed, such a simple (yet general) model allows
for closed expressions of the work extracted and its efficiency depending on
the key parameters controlling the activity of the system. In particular, our
active Szilard engines exhibit a remarkable improvement of the performance
compared to their passive counterpart operating in thermal equilibrium, in
line with recent experimental results Zanin  _et al._ (2021).

Figure 1: Cartoon of the active Szilard engine. In the quasi-static regime,
the demon measures in which half of the box the active particle is located and
then places a wall (in green) with an attached weight at $x=L/2$. After a full
quasi-static expansion performing work $W_{exp}$, the green wall hits the box
and is removed by the demon. This puts the system in its initial state. In the
dynamic regime, the demon measures the position of the particle with precision
$\delta$ and the sign of the velocity. It puts a wall (in green) with
precision $\delta$ ahead of the particle and let it push against the wall for
a time $\tau$. Then, the wall is removed and the system is back in its initial
state.

Model. In the following, we focus on the motion of an overdamped ABP along the
direction perpendicular to the wall and model its perpendicular velocity
component as a random variable with finite correlation time $\tau_{M}$ Zöttl
and Stark (2016),

$\displaystyle\left\langle v(t)\right\rangle=0\,,\quad\left\langle
v(t)v(t^{\prime})\right\rangle=v_{0}^{2}e^{-\frac{|t-t^{\prime}|}{\tau_{M}}}+2D\delta(t-t^{\prime})$
(1)

where $D$ is the (translational) thermal diffusion coefficient. Equations (1)
have been shown to properly reproduce the statistical properties of ABPs (as
modeled by Eqs. (S1) and (S2) in the Suppl. Mat.) as well as of Run-and-Tumble
bacteria Marini Bettolo Marconi and Maggi (2015). The case of passive
particles in thermal equilibrium is characterized by $v_{0}=0$. In the active
case and in the quasi-static regime ($t\gg\tau_{M}$) the model under study is
equivalent to that of a particle moving in one dimension and hopping between
two states, namely moving right or left Malgaretti  _et al._ (2021). In the
following, we analyze the performance of the two possible driving protocols of
active Szilard engines: the quasi-static regime and the dynamic regime.

Quasi-static active Szilard engine. First we focus on the case in which the
engine moves on time scales larger than $\tau_{M}$ and hence the time
evolution of the density of the ABPs can be regarded as quasi-static. The
Maxwell demon in this type of Szilard engine operates as follows. A single ABP
is confined in a box with lateral area $A$ (see Fig.5). The demon detects the
position of the particle and inserts a wall in the middle of the box
($x=L/2$). Due to the reduction of volume available to the particle, the
pressure in the occupied chamber increases and the volume of the chamber
expands till the movable wall hits the box (at $x=0,L$). Then, the wall is
removed and the system assumes its initial state again.

Figure 2: Work extracted during the expansion, $W_{exp}$, normalized by the
equilibrium value of a Szilard engine, as function of Pe (a) or $\Gamma$ (b),
with $L/R=100$.

In the following, we consider just a single active particle111The same results
holds also for the case of non-interacting active particles by simply
multiplying the pressure by the number $N$ of active particles. placed in the
left half of the box as in Fig. 5. The density is $\bar{\rho}=1/(Ax)$ and,
under a mean-field approximation, the pressure is given by Malgaretti  _et
al._ (2021):

$\displaystyle\beta\Pi(x)=\frac{R^{2}}{A}\frac{\kappa^{3}}{\operatorname{Pe}^{2}+2\Gamma\kappa
x}$ (2)

where

$\displaystyle\kappa=\frac{\sqrt{\operatorname{Pe}^{2}+2\Gamma}}{R},\,\,\Gamma=\frac{R^{2}}{\tau_{M}D},\,\,\operatorname{Pe}=\frac{v_{0}R}{D}\,.$
(3)

Here, $\kappa$ is the inverse of the length that characterizes the
accumulation of ABPs at the walls Malgaretti  _et al._ (2021) with $R$ the
radius of the particle, $\Gamma$ is the dimensionless hopping rate, and
$\operatorname{Pe}$ the Péclet number, respectively. Accordingly, the work
extracted during the expansion becomes

$\displaystyle\beta W_{exp}=A\int_{\frac{L}{2}}^{L}\beta\Pi(x)dx=\beta
W_{0}\ln\left[\frac{1+\chi L}{1+\chi L/2}\right]\,,$ (4)

where $xA$ is the momentary volume available to the particle and we have
identified the effective length $\chi^{-1}$ and strength $W_{0}$ of the work
cycle,

$\displaystyle\chi^{-1}=\frac{\operatorname{Pe}^{2}R}{2\Gamma\sqrt{\operatorname{Pe}^{2}+2\Gamma}},$
$\displaystyle\quad\beta W_{0}=1+\frac{\operatorname{Pe}^{2}}{2\Gamma}\,.$ (5)

Figure 2 shows the dependence of the work extracted from the system as
function of $\operatorname{Pe}$ and $\Gamma$. As expected, in the limit of
passive systems, i.e., for $\operatorname{Pe}\rightarrow 0$ or
$\Gamma\rightarrow\infty$, $\beta W_{exp}$ approaches the work of the Szilard
engine for a thermal bath $\beta W^{th}_{exp}=\ln 2$ Brillouin (1960);
Parrondo  _et al._ (2015). At variance, for active systems
($\operatorname{Pe}>0$ and $\Gamma\ll\infty$), the work extracted in the case
of an active bath exceeds that of the thermal bath. In particular, we can
identify an intermediate regime, for $\sqrt{2\Gamma}<\operatorname{Pe}<2\Gamma
L/R$ (see Ref.Malgaretti  _et al._ (2021)), in which $\beta
W_{exp}\propto\operatorname{Pe}^{2}$ whereas $\beta
W_{exp}\propto\operatorname{Pe}$ in the asymptotic regime,
$\operatorname{Pe}>2\Gamma L/R$. Note that $W_{exp}\gg W_{exp}^{th}$ is not in
contradiction with the Landauer principle Landauer (1961); Bérut  _et al._
(2012); Parrondo  _et al._ (2015) since the latter only applies to a thermal
bath. In the case of active baths the work is done at the expense of the
energy injected in the system to keep the bath in the active state. The
features we have described so far are generic for the active Szilard engine,
namely the considerable excess work during expansion as compared to a thermal
bath and its correct asymptotic behavior for $\operatorname{Pe}\rightarrow 0$
or $\Gamma\rightarrow\infty$. However, the total extracted work during one
cycle $W_{cyc}$ depends on the very realization of the active engine. In
Suppl. Mat. Sec. B, we present an alternative protocol. It includes a
compression step which results in a different dependence of $W_{cyc}$ on
$\operatorname{Pe}$ and $\Gamma$ including local maxima.

According to our model the energy input in the system over the cycling time
$\tau$ is due to the measurement $\mathcal{M}$ and due to the energy
$\mathcal{W}\tau$ consumed, on average, by the ABP during the cycle period
$\tau$ to keep itself in the active state. The input energy is partially
dissipated and partially converted into the work, $W_{exp}$ Kawai  _et al._
(2007). Therefore we have $\mathcal{W}\tau+\mathcal{M}=W_{exp}+W_{diss}$ with

$\displaystyle
W_{diss}=\mathcal{P}\tau+\frac{\gamma{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}_{w}}L^{2}}{\tau}+\mathcal{M}\,.$
(6)

Note, since the measurement does not contribute to the work, we add it to the
dissipated energy. Furthermore, $\mathcal{P}\tau$ is the energy dissipated, on
average, by the active particle during the cycle period $\tau$ and
$\gamma_{w}L^{2}/\tau$ is the (approximated) energy dissipated by the
piston222For heat engines that are in contact with thermal baths at different
temperatures the heat exchange with the thermal baths should be added to
$\langle W_{diss}\rangle$., where $\gamma_{w}$ is its effective friction
coefficient.

Next, we define the efficiency Schmiedl and Seifert (2007); Kjelstrup  _et
al._ (2010); Malgaretti  _et al._ (2021) as

$\eta=\frac{W_{exp}}{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}W_{exp}}+W_{diss}}=\left[{1+\dfrac{W_{diss}}{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}W_{exp}}}}\right]^{-1}\,.$
(7)

We remark that for the quasi-static processes assumed to derive Eq. (4),
$W_{exp}$ does not depend on $\tau$ and the dependence of $\eta$ on $\tau$ is
determined solely by the dissipated energy $W_{diss}$. Interestingly,
$W_{diss}$ diverges for both $\tau\rightarrow 0$ and $\tau\rightarrow\infty$
and it has a minimum at Malgaretti  _et al._ (2021)
$\tau_{opt}=\sqrt{\gamma_{w}L^{2}/\mathcal{P}}$ which, via Eq. (7), implies a
maximum for the efficiency. In order to be consistent with the quasi-static
assumption, $\tau_{opt}$ should exceed the typical relaxation time of the
density profile that, for large Pe, scales as $\tau_{rlx}\simeq L/v_{0}$.
Hence, $\tau_{rlx}<\tau_{opt}$ implies

$\displaystyle\gamma_{w}v_{0}^{2}>\mathcal{P}\,.$ (8)

For microswimmers it is well known that only a small part of the consumed mean
power $\mathcal{W}$ can be converted into work such as pushing a piston, while
the rest is needed for the swim mechanism. Estimating this work power as
$\gamma_{p}v_{0}^{2}$, where $\gamma_{p}$ is the friction coefficient of the
ABP, the efficiency  $\alpha$ of the ABP in a homogeneous and unbound fluid
becomes

$\displaystyle\alpha=\gamma_{p}v_{0}^{2}/\mathcal{W}\ll 1\,.$ (9)

Reported values for diffusiophoretic colloids are $\alpha\simeq
10^{-9}-10^{-5}$ Sabass and Seifert (2012); Shah  _et al._ (2020) whereas for
living organisms $\alpha\simeq 10^{-3}$ (Chlamydomonas Friedrich (2018)), and
$\alpha\simeq 10^{-2}$ (Paramecium Katsu-Kimura  _et al._ (2009) and
demembranated sperm flagella Chen  _et al._ (2015)). Thus, with
$\mathcal{W}\gtrsim\mathcal{P}$ and, using Eqs. (8)-(9), we obtain

$\displaystyle\gamma_{w}\alpha>\gamma_{p}\,.$ (10)

Since both, $\gamma_{p}$ and $\gamma_{w}$, scale with their linear size, the
constraint in Eq. (10) together with $\alpha\ll 1$ requires a significant
length scale separation between the ABP and the piston. This is also due to
the quasi-static process, we considered so far, where pressure instantaneously
relaxes to its stationary value. Therefore, now we formulate an engine where
this constraint is released.

Dynamic active Szilard engine The intrinsic out-of-equilibrium nature of
active systems provides additional means to control their dynamics. For
example, in the case of Active Brownian Particles (ABP) one can exploit the
finite correlation time of their velocity to design a specific “demon”. Within
such a scenario, the protocol of the demon works as follows. Every time $\tau$
the demon measures the position, with finite precision333The relation between
the precision of the position measurement, $\delta$, and the measurement cost,
$\mathcal{M}$, is not trivial due to the active nature of the ABP. We can
estimate it from the thermal case where $\mathcal{M}\simeq-
k_{B}T\ln(\delta/2L)$. The factor 2 comes from measuring the direction of
motionBrillouin (1960). $\delta$, and the direction of motion of the ABP and
it puts a wall ahead of the particle with precision $\delta$. The wall is
connected to a weight that opposes the motion of the particle with a force $F$
(see Fig.5). Thus, after colliding with the wall, the ABP pushes against it
and, hence, acts against both the force $F$ and the drag force of the wall,
$\gamma_{w}v_{w}$, where $\gamma_{w}$ is the friction coefficient of the wall
and $v_{w}$ its velocity. Accordingly, at contact with the wall the ABP moves
at a reduced velocity

$\displaystyle
v_{w}=\frac{\gamma_{p}v_{0}-F}{\gamma_{p}+\gamma_{w}}=v_{0}\frac{1-F/\gamma_{p}v_{0}}{1+\gamma_{w}/\gamma_{p}}\equiv
v_{0}\bar{v}_{w}$ (11)

with $\gamma_{p}$ the friction coefficient of the ABP. Now, at time $t=0$ the
velocity $v_{w}$ of the particle is measured while pushing against the wall.
Then, the work done against the conservative force $F$ in between two
measurements is

$W_{act}=F\int_{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\tau_{\delta}}}^{\tau}v(t)\big{|}_{v_{w}}dt\,,$
(12)

where $v(t)\big{|}_{v_{w}}$ means that at $t=0$ the velocity is $v_{w}$ and
$\tau_{\delta}\simeq\delta/v_{0}$ is the time the particle takes to reach the
wall placed with precision $\delta$. Accordingly, in the case of large
activity $v_{w}^{2}\gg D/\tau_{M}$ and relatively high measurement precision
$\tau_{\delta}\ll\min(\tau_{M},\tau)$, the average work is linear in
$\tau_{M}$ and becomes (see Suppl. Mat. Sec. B)

$\displaystyle\left\langle
W_{act}\right\rangle{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\simeq}Fv_{w}\tau_{M}\big{(}1-e^{-\bar{\tau}}\big{)}\,,$
(13)

where we introduced the reduced measurement time $\bar{\tau}=\tau/\tau_{M}$.
Introducing the rescaled work $\overline{W}_{act}=\langle
W_{act}\rangle/(\gamma_{p}v_{0}^{2}\tau_{M})$, we realize that it depends on
two dimensionless parameters, the reduced conservative force
$\bar{F}=F/\gamma_{p}v_{0}$ and the ratio of drag coefficients
$\bar{\gamma}=\gamma_{w}/\gamma_{p}$.

Figure 3:  Normalized work $\overline{W}_{act}=\langle
W_{act}\rangle/(\gamma_{p}v_{0}^{2}\tau_{M})$ extracted during a measurement
interval of rescaled duration $\bar{\tau}$ as function of $\bar{\tau}$ with
$\gamma_{w}/\gamma_{p}=1$ (a) or $\bar{F}$ with $\bar{\tau}=1$ (b).

Figure 3(a) plots the dependence of the rescaled $\left\langle
W_{act}\right\rangle$ on $\bar{\tau}$. As expected, in the limit of thermal
equilibrium, $\tau_{M}\rightarrow 0$, $\left\langle W_{act}\right\rangle$
vanishes, whereas for $\tau\gg\tau_{M}$, $\langle\mathcal{W}_{act}\rangle$
attains its maximum value $Fv_{w}\tau_{M}$, which in rescaled units becomes
$\bar{F}\bar{v}_{w}$. Due to Eq. (11) for $\bar{v}_{w}$, $\left\langle
W_{act}\right\rangle$ has a quadratic dependence on $\bar{F}$ [see Fig. 3(b)]
and it attains a maximum at $\bar{F}=1/2$, _i.e._ , when the force exerted by
the piston is half of the stall force of the ABP. The magnitude of
$\beta\left\langle W_{act}\right\rangle$ can be estimated for an ABP with
radius $R=\bar{R}\cdot 10^{-6}\mu\text{m}$, velocity $v_{0}=\bar{v}_{0}\cdot
10^{-6}\text{m/sec}$, and moving in water with $\tau_{M}=1/D_{r}=8\pi\eta
R^{3}/k_{B}T$ as $\beta\gamma_{p}v_{0}^{2}\tau_{M}\simeq
25\bar{R}^{4}\bar{v}^{2}_{0}$. Accordingly, for $\bar{R}=2$ and $\bar{v}=2$
the work extracted at $\bar{F}=1/2$ is $\beta\left\langle
W_{act}\right\rangle\simeq 240$ whereas for biological swimmers it amounts to
$\beta\left\langle W_{act}\right\rangle\simeq 10^{8}$ for Clamydomonas
($R=10\mu\text{m}$, $v=50\mu\text{m/sec}$, see Ref.Friedrich (2018)) and
$\beta\left\langle W_{act}\right\rangle\simeq 10^{12}$ for Paramecium
($R=25\mu\text{m}$, $v=10^{3}\mu\text{m/sec}$, see Ref.Katsu-Kimura  _et al._
(2009)).

Interestingly, in contrast to the quasi-static Szilard engine [see Eq. (4)],
in the current case the average work $\left\langle W_{act}\right\rangle$
retains an explicit dependence on $\tau$ that will be also visible in the
efficiency. In order to compute its mean value using the definition of Eq.
(7), we need to introduce the total dissipated energy per time step in full
analogy to Eq. (6),

$\displaystyle\left\langle
W_{diss}\right\rangle=\mathcal{P}\tau+\mathcal{M}+\left\langle
W^{pst}_{diss}\right\rangle\,,$ (14)

where

$\displaystyle\left\langle
W^{pst}_{diss}\right\rangle=\frac{\gamma_{w}v_{w}^{2}\tau_{M}}{2}\overline{W}^{pst}_{diss}\enspace\mathrm{with}\enspace\overline{W}^{pst}_{diss}=1-e^{-2\bar{\tau}}$
(15)

is the energy dissipated by the piston (see Suppl. Mat. Sec.C). In the limit
of small measurement times, $\bar{\tau}\ll 1$, it approaches $\left\langle
W^{pst}_{diss}\right\rangle=\gamma_{w}v_{w}^{2}\tau$, as expected, while for
long measurement times, $\bar{\tau}\gg 1$, it plateaus at $\left\langle
W^{pst}_{diss}\right\rangle=\gamma_{w}v_{w}^{2}\tau_{M}/2$.

Figure 4: Efficiency $\left\langle\tilde{\eta}\right\rangle$ as function of
$\bar{\tau}$, for $\overline{\mathcal{M}}=0.1,1,10$ lighter colors standing
for larger values of $\overline{\mathcal{M}}$. The black dot-dashed line is
$\propto 1/\bar{\tau}$, the black dashed line is $\propto\bar{\tau}$, and the
grey dashed line is the loci of the maxima for
$\overline{\mathcal{M}}\in[0.01,100]$.

Finally, the mean efficiency after some algebra (see Suppl. Mat. Sec.D),
becomes

$\displaystyle\left\langle\eta\right\rangle$
$\displaystyle\simeq\left[1+\dfrac{1}{\alpha}\dfrac{\bar{\tau}+\overline{\mathcal{M}}+\frac{1}{2}\alpha\bar{\gamma}\bar{v}_{w}^{2}\overline{W}^{pst}_{\\!diss}}{\overline{W}_{\\!act}}\right]^{-1}$
(16)

where we approximated $\gamma_{p}v_{0}^{2}/\mathcal{P}\simeq\alpha$ and we
have introduced the dimensionless energy cost of measurement
$\overline{\mathcal{M}}=\frac{\mathcal{M}}{\mathcal{P}\tau_{M}}$. Equation
(16) can be further simplified since $\alpha\ll 1$ and
$\bar{\gamma}\bar{v}_{w}^{2}\overline{W}^{pst}_{\\!diss}\lesssim\bar{\tau}$,
which follows from Eq. (15) and the definition of $\bar{\gamma}$ and
$\bar{v}_{wall}$ (for more details, see Suppl. Mat. Sec. E). Hence, we obtain

$\displaystyle\langle\eta\rangle\simeq\alpha\dfrac{\overline{W}_{act}}{\bar{\tau}+\overline{\mathcal{M}}}\,.$
(17)

Figure 4 shows that the efficiency of the dynamic Szilard engine,
$\langle\eta\rangle$ is bound by that of the ABP, $\alpha$. Moreover,
$\langle\eta\rangle$ has a non-monotonous dependence on $\bar{\tau}$
independently whether the measurement cost is much smaller then
($\overline{\mathcal{M}}\ll 1$), similar to ($\overline{\mathcal{M}}\simeq
1$), or larger then ($\overline{\mathcal{M}}\gg 1$) the energy dissipated by
the ABP to keep in the active state. This is at variance to the monotonic
increase of $\langle W_{act}\rangle$ in $\bar{\tau}$ [see Fig. 3(a)] and is
due to the linear increase of the power dissipated by the ABP with
$\bar{\tau}$. Interestingly, the value, $\bar{\tau}_{opt}$, at which the
efficiency is maximized grows with $\overline{\mathcal{M}}$ as
$\bar{\tau}_{opt}\simeq\ln\overline{\mathcal{M}}$ (see the grey-dashed line in
Fig. 4 and Fig. S3), whereas the maximum value of $\langle\eta\rangle$ can
reach up to $10\%$ of the efficiency of the ABP. Finally, $\langle\eta\rangle$
depends on $\bar{F}$ solely through $\langle W_{act}\rangle$ and hence it
retains the non-monotonous dependence on $\bar{F}$ [see Fig.3(b)].

In this article we have discussed the general features of Szilard engines
operating in contact with active baths. Importantly, the active bath differs
from a thermal bath by the persistent motion of its constituents. The
resulting enhanced active pressure and the exponential time correlations of
the velocity enable two designs of an active Szilard engine that exploit
different regimes.

In the quasi-static regime the slow expansion of the active Szilard engine on
times much larger than the typical relaxation time of the ABP’s density is
driven by the active pressure. Its magnitude depends on both the Péclet number
as well as the box size Malgaretti  _et al._ (2021). The  work extracted from
the active bath during expansion outnumbers the equilibrium counterpart. Thus,
the Landauer limit Landauer (1961); Bérut  _et al._ (2012); Parrondo  _et al._
(2015) can be overcome at the expense of the additional energy pumped into the
system by the active bath. Since the active bath undergoes a quasi-static work
cycle, $W_{exp}$ is independent of the cycling time $\tau$. The efficiency
depends on $\tau$ solely through the energy dissipated by the ABP and the
piston such that it becomes maximum at a finite value $\tau_{opt}$. For longer
cycling times, the dissipation is essentially due to the energy spent to keep
the system active and the efficiency scales $\propto 1/\tau$.

In the dynamic regime the active Szilard engine probes the velocity
correlation time, $\tau_{M}$, of the ABP. This implementation of the active
Szilard engine has no practical counterpart for thermal baths. In fact, for an
ABP with radius $R=1\mu\text{m}$ suspended in water, we have $\tau_{M}\simeq
D_{r}^{-1}=8\pi\eta R^{3}/k_{B}T\simeq 6\text{sec}$. In contrast, for a
passive colloid of the same size and also suspended in water the velocity
correlation time is $\tau_{th}\simeq m/\gamma_{p}\simeq 10^{-7}\text{sec}$,
where $m$ is the mass of the colloidal particle. Accordingly, while an ABP can
perform work on a relatively slow time scale of seconds, for a passive
particle the demon would need to monitor velocity correlations on a very fast
time scale of $10^{-7}\sec$. Interestingly, the work extracted,
$\langle{W}_{act}\rangle$, becomes maximum when $F$ is half of the stall force
$F_{opt}=\gamma_{p}v_{0}/2$. Since the extracted work scales as
$\beta\langle{W}_{act}\rangle\propto R^{4}v_{0}^{2}$, it varies in a large
range from $1$ for micron-sized colloids up to $10^{12}$ for a _Paramecium_.
The dynamic Szilard engine works with an efficiency $\langle\eta\rangle$
bounded from above by the efficiency $\alpha$ of an ABP. Interestingly,
$\langle\eta\rangle$ shows a non-monotonous dependence on both the
dimensionless measuring time $\bar{\tau}=\tau/\tau_{M}$ with a maximum at
$\bar{\tau}_{opt}\approx 1$, as well as on the conservative force with the
maximum attained at $F_{opt}=\gamma_{p}v_{0}/2$. For biological swimmers
$\alpha\in[10^{-4},10^{-2}]$, hence according to Fig. 4 the efficiency can be
as high as $\langle\eta\rangle\simeq 10^{-3}$, which is quite larger than
typical efficiencies of other micrometric engines with values smaller than
$\sim 10^{-8}$ (see the diverse systems discussed in Ref. Maggi  _et al._
(2015)).

Our results highlight how the persistent motion of active bath particles and
the resulting features of enhanced active pressure and long-time velocity or
orientational correlations, determine the dynamics of information-based
engines. These features can be exploited to design novel micro- and nano-
engines that outperform those relying on equilibrium baths. Using the method
of optical video microscopy and real-time image analysis (see, _e.g._ , Tobias
Bäuerle and Bechinger (2020)), we envisage the possibility to mimic the
Maxwell demon and to ultimately realize an active Szilard engine.

## Aknowledgments

P.M. acknoledges Wendong Wang, Gaurav Gardi, and Vimal Kishore for useful
discussions.

## References

  * tai (1867) _Life and Scientific Work of Peter Guthrie Tait_ (Cambridge University Press, 1867) p. 213–215.
  * Maxwell (1871) J. Maxwell,  _Theory of Heat_ (Longmans, Green and Co., London, 1871) Chap. 12.
  * Feynman (1963) R. Feynman,  _The Feynman Lectures on Physics, Vol. 1_ (1963).
  * Reimann (2002) P. Reimann, Phys. Rep. 361, 57 (2002).
  * Szilard (1929) L. Szilard, (1929).
  * Leff and Rex (2002) H. Leff and A. F. Rex,  _Maxwell’s Demon 2_ (Tyalor & Francis, Bristol, 2002).
  * Brillouin (1960) L. Brillouin,  _Science and Information Theory, 2nd Ed._ (Dover Publications, New York, 1960).
  * Parrondo  _et al._ (2015) J. M. R. Parrondo, J. M. Horowitz, and T. Sagawa, Nature Physics 11, 131 (2015).
  * Koski  _et al._ (2014) J. V. Koski, V. F. Maisi, J. P. Pekola, and D. V. Averin, Proceedings of the National Academy of Sciences 111, 13786 (2014).
  * Bengtsson  _et al._ (2018) J. Bengtsson, M. N. Tengstrand, A. Wacker, P. Samuelsson, M. Ueda, H. Linke, and S. M. Reimann, Phys. Rev. Lett. 120, 100601 (2018).
  * Aydiner (2021) E. Aydiner, Scientific Reports 11, 1576 (2021).
  * Ribezzi-Crivellari and Ritort (2019) M. Ribezzi-Crivellari and F. Ritort, Nature Physics 15, 660 (2019).
  * Paneru and Pak (2020) G. Paneru and H. K. Pak, Advances in Physics: X 5, 1823880 (2020).
  * Marchetti  _et al._ (2013) M. C. Marchetti, J. F. Joanny, S. Ramaswamy, T. B. Liverpool, J. Prost, M. Rao, and R. A. Simha, Rev. Mod. Phys. 85, 1143 (2013).
  * Ebbens and Howse (2010) S. J. Ebbens and J. R. Howse, Soft Matter 6, 726 (2010).
  * Bechinger  _et al._ (2016) C. Bechinger, R. Di Leonardo, H. Löwen, C. Reichhardt, G. Volpe, and G. Volpe, Rev. Mod. Phys. 88, 045006 (2016).
  * Zöttl and Stark (2016) A. Zöttl and H. Stark, Journal of Physics: Condensed Matter 28, 253001 (2016).
  * Doostmohammadi  _et al._ (2018) A. Doostmohammadi, J. Ignés-Mullol, J. M. Yeomans, and F. Sagués, Nature Communications 9, 3246 (2018).
  * Gompper  _et al._ (2020) G. Gompper, R. G. Winkler, T. Speck, A. Solon, C. Nardini, F. Peruani, H. Löwen, R. Golestanian, U. B. Kaupp, L. Alvarez, T. Kiørboe, E. Lauga, W. C. K. Poon, A. DeSimone, S. Muiños-Landin, A. Fischer, N. A. Söker, F. Cichos, R. Kapral, P. Gaspard, M. Ripoll, F. Sagues, A. Doostmohammadi, J. M. Yeomans, I. S. Aranson, C. Bechinger, H. Stark, C. K. Hemelrijk, F. J. Nedelec, T. Sarkar, T. Aryaksama, M. Lacroix, G. Duclos, V. Yashunsky, P. Silberzan, M. Arroyo, and S. Kale, Journal of Physics: Condensed Matter 32, 193001 (2020).
  * Romanczuk  _et al._ (2012) P. Romanczuk, M. Bär, W. Ebeling, B. Lindner, and L. Schimansky-Geier, Eur. Phys. J. Spec. Top. 202, 1 (2012).
  * Rotschild (1963) L. Rotschild, Nature , 198 (1963).
  * Elgeti  _et al._ (2015) J. Elgeti, R. G. Winkler, and G. Gompper, Rep. Prog. Phys. 78, 056601 (2015).
  * Cates and Tailleur (2013) M. E. Cates and J. Tailleur, EPL (Europhys. Lett.) 101, 20010 (2013).
  * Cates and Tailleur (2015) M. E. Cates and J. Tailleur, Annu. Rev. Condens. Matter Phys. 6, 219 (2015).
  * Zheng  _et al._ (2011) M. C. Zheng, F. M. Ellis, T. Kottos, R. Fleischmann, T. Geisel, and T. c. v. Prosen, Phys. Rev. E 84, 021119 (2011).
  * Pietzonka  _et al._ (2019) P. Pietzonka, É. Fodor, C. Lohrmann, M. E. Cates, and U. Seifert, Phys. Rev. X 9, 041032 (2019).
  * Holubec  _et al._ (2020) V. Holubec, S. Steffenoni, G. Falasco, and K. Kroy, Phys. Rev. Res. 2, 043262 (2020).
  * Malgaretti  _et al._ (2021) P. Malgaretti, P. Nowakowski, and H. Stark, EPL (Europhysics Letters) 134, 20002 (2021).
  * Fodor and Cates (2021) É. Fodor and M. E. Cates, EPL 134, 10003 (2021).
  * Gronchi and Puglisi (2021) G. Gronchi and A. Puglisi, Phys. Rev. E 103, 052134 (2021).
  * Speck (2022) T. Speck, Phys. Rev. E 105, L012601 (2022).
  * Ion and Urna (2022) S. Ion and B. Urna, Arxiv , 2201.00796 (2022).
  * Sokolov  _et al._ (2010) A. Sokolov, M. M. Apodaca, B. A. Grzybowski, and I. S. Aranson, Proc. Natl. Acad. Sci. 107, 969 (2010).
  * Di Leonardo  _et al._ (2010) R. Di Leonardo, L. Angelani, D. Dell’Arciprete, G. Ruocco, V. Iebba, S. Schippa, M. P. Conte, F. Mecarini, F. De Angelis, and E. Di Fabrizio, Proc. Natl. Acad. Sci. 107, 9541 (2010).
  * Malgaretti  _et al._ (2013) P. Malgaretti, I. Pagonabarraga, and J. M. Rubi, J. Chem. Phys. 138, 194906 (2013).
  * Malgaretti  _et al._ (2014) P. Malgaretti, I. Pagonabarraga, and J. Rubi, European Physics Journal Special Topics 223, 3295 (2014).
  * Yang and Ripoll (2014) M. Yang and M. Ripoll, Soft Matter 10, 1006 (2014).
  * Michelin  _et al._ (2015) S. Michelin, T. D. Montenegro-Johnson, G. De Canio, N. Lobato-Dauzier, and E. Lauga, Soft Matter 11, 5804 (2015).
  * Malgaretti and Stark (2017) P. Malgaretti and H. Stark, J. Chem. Phys. 146, 174901 (2017).
  * Holubec  _et al._ (2017) V. Holubec, A. Ryabov, M. H. Yaghoubi, M. Varga, A. Khodaee, M. E. Foulaadvand, and P. Chvosta, Entropy 19 (2017), 10.3390/e19040119.
  * Zanin  _et al._ (2021) G. L. Zanin, M. Antesberger, M. J. Jacquet, P. H. S. Ribeiro, L. A. Rozema, and P. Walther, (2021), arXiv:2107.09686 [quant-ph] .
  * Marini Bettolo Marconi and Maggi (2015) U. Marini Bettolo Marconi and C. Maggi, Soft Matter 11, 8768 (2015).
  * Landauer (1961) R. Landauer, IBM Journal of Research and Development 5, 183 (1961).
  * Bérut  _et al._ (2012) A. Bérut, A. Arakelyan, A. Petrosyan, S. Ciliberto, R. Dillenschneider, and E. Lutz, Nature 483, 187 (2012).
  * Kawai  _et al._ (2007) R. Kawai, J. M. R. Parrondo, and C. V. den Broeck, Phys. Rev. Lett. 98, 080602 (2007).
  * Schmiedl and Seifert (2007) T. Schmiedl and U. Seifert, EPL (Europhys. Lett.) 81, 20003 (2007).
  * Kjelstrup  _et al._ (2010) S. Kjelstrup, D. Bedeaux, E. Johannessen, and J. Gross,  _Non-Equilibrium Thermodynamics for Engineers_ (World Scientific Publishing Company, Singapore, 2010).
  * Sabass and Seifert (2012) B. Sabass and U. Seifert, The Journal of Chemical Physics 136, 064508 (2012).
  * Shah  _et al._ (2020) Z. H. Shah, S. Wang, L. Xian, X. Zhou, Y. Chen, G. Lin, and Y. Gao, Chem. Commun. 56, 15301 (2020).
  * Friedrich (2018) B. M. Friedrich, Phys. Rev. E 97, 042416 (2018).
  * Katsu-Kimura  _et al._ (2009) Y. Katsu-Kimura, F. Nakaya, S. A. Baba, and Y. Mogami, Journal of Experimental Biology 212, 1819 (2009).
  * Chen  _et al._ (2015) D. T. Chen, M. Heymann, S. Fraden, D. Nicastro, and Z. Dogic, Biophysical Journal 109, 2562 (2015).
  * Maggi  _et al._ (2015) C. Maggi, F. Saglimbeni, M. Dipalo, F. De Angelis, and R. Di Leonardo, Nature Communications 6, 7855 (2015).
  * Tobias Bäuerle and Bechinger (2020) R. C. L. Tobias Bäuerle and C. Bechinger, , 2547 (2020).

##  I Supplemental Material

###  I.1 Active Brownian Particle (ABP)

The motion of an active Brownian particle is governed by (see Ref. Romanczuk
_et al._ (2012); Zöttl and Stark (2016))

$\displaystyle\gamma_{t}\dot{\mathbf{x}}(t)=\,v_{0}\bm{\theta}(t)+\bm{\xi}(t)\,,\quad\dot{\bm{\theta}}=\bm{\eta}(t)\times\bm{\theta}$
(18)

where, $v_{0}$ is the magnitude of the active velocity, $\bm{\theta}$ its
direction and, introducing the translational and rotational diffusion
coefficients $D_{t,r}$, $\bm{\eta}$ and $\bm{\xi}$ are delta-correlated white
noise

$\displaystyle\langle\bm{\xi}(t)\rangle=$ $\displaystyle
0\,,\quad\langle\bm{\xi}(t)\bm{\xi}(t^{\prime})\rangle=2D_{t}\delta(t-t^{\prime})$
(19a) $\displaystyle\langle\bm{\eta}(t)\rangle=$ $\displaystyle
0\,,\quad\langle\bm{\eta}(t)\bm{\eta}(t^{\prime})\rangle=2D_{r}\delta(t-t^{\prime})$
(19b)

This model leads to a finite correlation time of the active component of
velocity (see Ref. Romanczuk  _et al._ (2012); Zöttl and Stark (2016))

$\displaystyle
v_{0}^{2}\langle\bm{\theta}(t)\cdot\bm{\theta}(t^{\prime})\rangle=v_{0}^{2}e^{-\frac{|t-t^{\prime}|}{\tau_{M}}}$
(20)

with $\tau_{M}=\frac{1}{2D_{r}}$.

###  I.2 Alternative realization of the quasi-static active Szilard engine

Figure 5: Cartoon of an alternative realization of the quasi-static active
Szilard engine. The daemon divides the cylinder by an additional wall (in
green) and at the same time fixes the piston (in green) of the empty part of
the cylinder. After the expansion of the right piston, the daemon removes the
wall, let the system compress to the original size, and then releases the left
piston.

The Maxwell daemon in this type of Szilard engine operates as follows. A
single particle confined between two pistons of area $A$ exerts a pressure
$\beta\Pi_{0}$ on the pistons such that the stationary distance between them
is $L$444For a single ideal gas particle $\beta\Pi_{0}=(AL)^{-1}$. (see
Fig.5). The daemon detects the position of the particle, inserts a wall in the
middle of the box ($x=0$), and at the same time fixes the position of the
piston of the empty chamber (depicted in green in Fig.5). Due to the reduction
of volume available to the particle, the pressure in the occupied chamber
increases and the volume of the chamber expands till $x=L$, i.e., when the
pressure matches the outer pressure $\Pi_{0}$. Once $x=L$ is reached, the
expansion stops since, by definition, $\Pi=\Pi_{0}$. Then, the wall is
released and the piston returns to $L/2$. The work done by this engine during
the expansion mathces $W_{exp}$ as reported in the main text. However, in the
current case a (negative) work is performed during compression. Accordingly,
using Eq. (2) we calculate the work during compression

$\displaystyle\beta$ $\displaystyle
W_{comp}=A\int_{\frac{3}{2}L}^{L}\\!\\!\\!\beta\Pi\left(x\right)dx=\beta
W_{0}\ln\left[\frac{1+\chi L}{1+3/2\chi L}\right]\,.$ (21)

Figure 6: Work extracted during a work cycle, $W_{cyc}$, normalized by the
value for thermal equilibrium, as function of Pe (upper panel) or $\Gamma$
(lower panel) with $L/R=100$.

We remark that during compression the internal wall has been removed and hence
the pressure of the ABP should be calculated taking into account that now the
volume available to the particle is $A(x+L/2)$. Accordingly, the work over a
cycle, $\beta W_{cyc}=\beta(W_{exp}+W_{comp})$, reads:

$\displaystyle\beta W_{cyc}=-\beta W_{0}\ln\left[1-\frac{1}{4}\left(\frac{\chi
L}{1+\chi L}\right)^{2}\right]\,.$ (22)

As shown in Fig. 6, the work extracted over a cycle attains a plateau at
$W_{cyc}\simeq\ln(4/3)$ for small values of the Péclet number,
$\operatorname{Pe}<\sqrt{2\Gamma}$, and at
$W_{cyc}\simeq\frac{\Gamma}{2}\frac{L^{2}}{R^{2}}$ for large values,
$\operatorname{Pe}>2\Gamma L/R$. For intermediate values of
$\operatorname{Pe}$, $W_{cyc}$ scales as $\operatorname{Pe}^{2}$.
Interestingly, the extracted work has a non monotonous dependence on $\Gamma$
and $W_{cyc}$ is maximized for $\Gamma_{max}\sim\operatorname{Pe}$ (see Fig.S1
in Suppl. Mat.). Finally, for very large values of $\Gamma$, the active bath
reduces to a thermal one, and the extracted work reaches a plateau $\beta
W_{cyc}\simeq\ln(4/3)$ i.e., the work extracted from a thermal bath.
Finally, we report the values of $\Gamma_{max}$ obtained by the numerical
solution of $\partial_{\Gamma}\beta W_{cyc}=0$ where $\beta W_{cyc}$ is
defined in Eq. (8) of the main text.

Figure 7: $\Gamma_{max}$, as function of $\operatorname{Pe}$.

It is interesting to note that within this realization of the Szilard engine,
the dependence of the extracted work on $\Gamma$ in non-monotonous, in
contrast to the monotonic decay shown by the Szilard engine reported in the
main text.

###  I.3 Derivation of Eq. (17) of the main text

In the following we report the derivation of Eq. (17). When the ABP is pushing
against the wall its nominal velocity is reduced to

$\displaystyle
v_{w}=\frac{\gamma_{p}v_{0}-F}{\gamma_{p}+\gamma_{w}}=v_{0}\frac{1-F/\gamma_{p}v_{0}}{1+\gamma_{w}/\gamma_{p}}\equiv
v_{0}\bar{v}_{w}$ (23)

where $\gamma_{p}$ is the friction coefficient of the ABP. Accordingly the
work reads

$\displaystyle\left\langle W_{act}\right\rangle$
$\displaystyle=F\int_{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\tau_{\delta}}}^{\tau}\left\langle
v(t)\big{|}_{v(0)}\right\rangle dt$ (24)
$\displaystyle=\frac{F}{v(0)}\int_{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\tau_{\delta}}}^{\tau}\\!\\!\left\langle
v(t)v(0)\right\rangle dt$ (25)
$\displaystyle=\frac{F}{v_{w}}\int_{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\tau_{\delta}}}^{\tau}\\!\\!\left\langle
v(t)v(0)\right\rangle dt$ $\displaystyle\text{Using the definition of }v(0)$
(26)
$\displaystyle=F\left(v_{w}+\frac{D}{v_{w}\tau_{M}}\right)\int_{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\tau_{\delta}}}^{\tau}e^{-\frac{t}{\tau_{M}}}dt$
Using Eq. (1) (27)
$\displaystyle=F\left(v_{w}+\frac{D}{v_{w}\tau_{M}}\right)\tau_{M}\left({\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}e^{-\frac{\tau_{\delta}}{\tau_{M}}}}-e^{-\frac{\tau}{\tau_{M}}}\right)$
integrating (28)
$\displaystyle=Fv_{w}\tau_{M}\left(1+\frac{D}{v_{w}^{2}\tau_{M}}\right)\left({\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}e^{-\frac{\tau_{\delta}}{\tau_{M}}}}-e^{-\frac{\tau}{\tau_{M}}}\right)$
$\displaystyle\text{grouping }v_{w}$ (29)

Finally, for large activity $\frac{D}{v_{w}^{2}\tau_{M}}\ll 1$ and high
measurement precision $\tau_{\delta}\ll\min(\tau_{M},\tau)$ the last
expression reduces to

$\displaystyle\left\langle
W_{act}\right\rangle\\!=Fv_{w}\tau_{M}\left(1-e^{-\bar{\tau}}\right)$ (30)

where we have introduced the time ratio $\bar{\tau}=\tau/\tau_{M}$. For later
use we introduce the the rescaled work

$\displaystyle\overline{W}_{\\!act}=\frac{\left\langle
W_{act}\right\rangle}{\gamma_{p}v_{0}^{2}\tau_{M}}=\frac{Fv_{w}}{\gamma_{p}v_{0}^{2}}\left(1-e^{-\bar{\tau}}\right)=\bar{F}\bar{v}_{w}\left(1-e^{-\bar{\tau}}\right)$
(31)

where we introduced the dimensionless external force,
$\bar{F}=F/v_{0}\gamma_{p}$ and we used Eq. (23). Then, we have

$\displaystyle\left\langle
W_{act}\right\rangle\\!=\gamma_{p}v_{0}^{2}\tau_{M}\overline{W}_{\\!act}\,.$
(32)

###  I.4 Derivation of Eq. (20) of the main text

Within a mean-field approximation, the energy dissipated by the piston during
the time $\tau$ can be expressed as:

$\displaystyle\langle W^{pst}_{diss}\rangle$
$\displaystyle=\gamma_{w}\int_{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\tau_{\delta}}}^{\tau}\\!\\!\left\langle
v(t)\big{|}_{v(0)}\right\rangle^{2}dt$ (33)
$\displaystyle=\frac{\gamma_{w}}{v^{2}(0)}\int_{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\tau_{\delta}}}^{\tau}\\!\\!\left\langle
v(t)v(0)\right\rangle^{2}dt$ (34)

In the large activity regime we can approximate $v(0)=v_{w}^{2}+D/\tau_{M}\sim
v_{w}^{2}$. Accordingly the last expression reduces to

$\displaystyle\langle W^{pst}_{diss}\rangle$
$\displaystyle=\gamma_{w}v_{w}^{2}\int_{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\tau_{\delta}}}^{\tau}e^{-2\frac{t}{\tau_{M}}}dt$
(36)
$\displaystyle=\frac{\tau_{M}}{2}\gamma_{w}v_{w}^{2}\left({\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}e^{-2\frac{\tau_{\delta}}{\tau_{M}}}}-e^{-2\frac{\tau}{\tau_{M}}}\right)$
(37)

As for the work, in the case of high measurement precision
$\tau_{\delta}\ll\min(\tau_{M},\tau)$ the previous expression reduces to

$\displaystyle\langle
W^{pst}_{diss}\rangle=\frac{\tau_{M}}{2}\gamma_{w}v_{w}^{2}\left(1-e^{-2\frac{\tau}{\tau_{M}}}\right)$
(38)

and the overall dissipated energy hence reads:

$\displaystyle\langle W_{diss}\rangle$
$\displaystyle=\mathcal{P}\tau+\mathcal{M}+\frac{\tau_{M}}{2}\gamma_{w}v_{w}^{2}\left(1-e^{-2\frac{\tau}{\tau_{M}}}\right)$
$\displaystyle\simeq\mathcal{P}\tau_{M}\left[\bar{\tau}+\overline{\mathcal{M}}+\frac{1}{2}\alpha\bar{\gamma}\bar{v}_{w}^{2}\overline{W}^{pst}_{\\!diss}\right]$
(39)

where we used

$\displaystyle\bar{\gamma}=$ $\displaystyle\frac{\gamma_{w}}{\gamma_{p}}$ (40)
$\displaystyle\overline{\mathcal{M}}=$
$\displaystyle\frac{\mathcal{M}}{\mathcal{P}\tau_{M}}$ (41)
$\displaystyle\frac{\gamma_{w}v_{0}^{2}}{\mathcal{P}}=\frac{\gamma_{p}v_{0}^{2}}{\mathcal{P}}\bar{\gamma}\simeq$
$\displaystyle\alpha\bar{\gamma}$ (42)
$\displaystyle\overline{W}^{pst}_{\\!diss}=$ $\displaystyle
1-e^{-2\bar{\tau}}$ (43)

Accordingly, the efficiency, Eq. (20), reads

$\displaystyle\langle\eta\rangle$ $\displaystyle=\dfrac{1}{1+\dfrac{\langle
W_{diss}\rangle}{\langle W_{act}\rangle}}$
$\displaystyle\simeq\dfrac{1}{1+\dfrac{1}{\alpha}\dfrac{\bar{\tau}+\overline{\mathcal{M}}+\frac{1}{2}\alpha\bar{\gamma}\bar{v}_{w}^{2}\overline{W}^{pst}_{\\!diss}}{\overline{W}_{\\!act}}}$
(44)

Finally that is Eq. (20) of the main text.

###  I.5 Magnitude of
$\frac{\gamma_{w}}{\gamma_{p}}\bar{v}_{w}^{2}\overline{W}^{pst}_{\\!diss}$

We comment on the magnitude of
$\frac{\gamma_{w}}{\gamma_{p}}\bar{v}_{w}^{2}\overline{W}^{pst}_{\\!diss}$.
For clarity sake, we express it as

$\displaystyle\frac{\gamma_{w}}{\gamma_{p}}\bar{v}_{w}^{2}\overline{W}^{pst}_{\\!diss}=\frac{\gamma_{w}}{\gamma_{p}}\frac{1-F/\gamma_{p}v_{0}}{1+\gamma_{w}/\gamma_{p}}\left(1-e^{-2\bar{\tau}}\right)$
(45)

Clearly the last expression is maximized at $\bar{\tau}\rightarrow\infty$ and
$F=0$. Accordingly, in the left panel of Fig.8 we report the values of
$\frac{\gamma_{w}}{\gamma_{p}}\bar{v}_{w}^{2}\overline{W}^{pst}_{\\!diss}$ for
$F=0$ and $\bar{\tau}=10$ while in the right panel the values of
$\frac{\gamma_{w}}{\gamma_{p}}\bar{v}_{w}^{2}\overline{W}^{pst}_{\\!diss}$
normalized by $\bar{\tau}$ as function of $\bar{\tau}$ for
$\gamma_{w}/\gamma_{p}=1$ and $\bar{F}=0$ i.e. at the maximum value of the
left panel.

Figure 8: Left: magnitude of
$\frac{\gamma_{w}}{\gamma_{p}}\bar{v}_{w}^{2}\overline{W}^{pst}_{\\!diss}$ as
function of $\frac{\gamma_{w}}{\gamma_{p}}$ with $F=0$ and
$\bar{\tau}=0.1,0.5,1,10$ where lighter colors stand for larger values of
$\bar{\tau}$. Right: magnitude of
$\frac{\gamma_{w}}{\gamma_{p}}\bar{v}_{w}^{2}\overline{W}^{pst}_{\\!diss}\big{/}\bar{\tau}$
as function of $\bar{\tau}$ for $\gamma_{w}/\gamma_{p}=1$ and $\bar{F}=0$ i.e.
at the maximum value of the left panel.

###  I.6 Dependence of $\bar{\tau}_{opt}$ on $\overline{\mathcal{M}}$

We report the values of $\tau_{opt}$ obtained by the numerical solution of
$\partial_{\bar{\tau}}\langle\tilde{\eta}\rangle=0$ where $\tilde{\eta}$ is
defined in Eq. (21) of the main text.

Figure 9: Dependence of $\bar{\tau}_{opt}$ on $\overline{\mathcal{M}}$ for any
values of $\bar{F}$ and $\bar{\gamma}$. The dot-dashed lines stands for
$\ln\left(\overline{\mathcal{M}}\right)$.