# Dynamics of Uniaxial-to-Biaxial Nematics Switching in Suspensions of Hard
Cuboids

Effran Mirzad Rafael  Department of Chemical Engineering and Analytical
Science, The University of Manchester, Manchester, M13 9PL, United Kingdom
Luca Tonti  Department of Chemical Engineering and Analytical Science, The
University of Manchester, Manchester, M13 9PL, United Kingdom  Daniel Corbett
Department of Chemical Engineering and Analytical Science, The University of
Manchester, Manchester, M13 9PL, United Kingdom  Alejandro Cuetos  Department
of Physical, Chemical and Natural Systems, Pablo de Olavide University, 41013
Sevilla, Spain  Alessandro Patti<EMAIL_ADDRESS>Department
of Chemical Engineering and Analytical Science, The University of Manchester,
Manchester, M13 9PL, United Kingdom

###### Abstract

Field-induced reorientation of colloidal particles is especially relevant to
manipulate the optical properties of a nanomaterial for target applications.
We have recently shown that surprisingly feeble external stimuli are able to
transform uniaxial nematic liquid crystals (LCs) of cuboidal particles into
biaxial nematic LCs. In the light of these results, here we apply an external
field that forces the reorientation of colloidal cuboids in nematic LCs and
sparks a uniaxial-to-biaxial texture switching. By Dynamic Monte Carlo
simulation, we investigate the unsteady-state reorientation dynamics at the
particle scale when the field is applied (uniaxial-to-biaxial switching) and
then removed (biaxial-to-uniaxial switching). We detect a strong correlation
between the response time, being the time taken for the system to reorient,
and particle anisotropy, which spans from rod-like to plate-like geometries.
Interestingly, self-dual shaped cuboids, theoretically considered as the most
suitable to promote phase biaxiality for being exactly in between prolate and
oblate particles, exhibit surprisingly slow response times, especially if
compared to prolate cuboids.

##  I Introduction

Colloids are biphasic systems comprising particles homogeneously dispersed in
a medium. In colloidal suspensions, the dispersed phase consists of solid
particles, while the continuous phase is a liquid. In particular, the
dispersed particles should have at least in one direction a dimension roughly
between 1 nm and 1 $\mu$m, so that gravitational and thermal forces compensate
each other Mormann  _et al._ (2017). This balance allows the dispersed
particles to remain suspended and to diffuse randomly via Brownian motion,
named after the Scottish botanist Robert Brown who, in 1827, described the
persistent and casual jumpy moves of organelles suspended in water Brown
(2015). If the suspended particles are anisotropic, under certain conditions,
they can self-organise into liquid-crystalline phases. Liquid crystals (LCs)
are mesophases that flow like liquids but, exhibit a significant degree of
internal ordering like crystals. A common LC morphology is the uniaxial
nematic ($\rm N_{U}$) phase where particles have one axis pointing
collectively in the same direction, but their centres of mass are randomly
distributed. This merely orientational ordering allows nematic LCs to exhibit
optical birefringence while maintaining mechanical fluidity. Currently, LCs
are deployed in a multitude of optical technology, including (but not limited
to) commercial displays Chen  _et al._ (2018), displays for virtual augmented
reality Huang  _et al._ (2018) and even smart windows with memory displays
Dabrowski  _et al._ (2018). The roll-out of these high-tech products are
coupled alongside further advances to understand and control the morphology of
LCs Meyer  _et al._ (2020) as well as techniques to optimise how they are
manufactured Cousins  _et al._ (2019), highlighting the relevance of LCs in
today’s research landscape. Very recently, there has been reignited interest
in the biaxial nematic phase ($\rm N_{B}$) and its potential to be
incorporated into display technology. In contrast to the $\rm N_{U}$ phase,
the $\rm N_{B}$ phase possesses two optical axes due to the alignment of the
three directors, making it very appealing for the design of nanomaterials with
novel optical properties Luckhurst and Sluckin (2015). Equally important, the
$\rm N_{B}$ phase is also foreseen to realise faster switching through its
minor axis switching mode, an aspect that could potentially improve refresh
rates in displays Lee  _et al._ (2007); Berardi and Zannoni (2008); Ricci  _et
al._ (2015). Despite these promising features, the existence of stable
molecular $\rm N_{B}$ phases is still an ongoing debate within the LC
community. While a biaxial geometry is indeed necessary to observe $\rm N_{B}$
phases, it has been shown that this requirement might not be sufficient as the
$\rm N_{B}$ phase tends to be metastable with respect to other phases,
including $\rm N_{U}$ and smectic (Sm) LCs Taylor and Herzfeld (1991). This is
for instance the case of colloidal cuboids, which can only form $\rm N_{B}$
phases at sufficiently large size dispersity Belli  _et al._ (2011); Rafael
_et al._ (2020), extreme anisotropy Dussi  _et al._ (2018), in the presence of
depletants Belli  _et al._ (2012) or upon application of an external field
Leferink op Reinink  _et al._ (2014); Cuetos  _et al._ (2019). Unless at least
one of these conditions are met, systems of monodisperse or bidisperse cuboids
cannot form $\rm N_{B}$ phases Cuetos  _et al._ (2017); Patti and Cuetos
(2018). For the monodisperse case, simulation results showed excellent
qualitative agreement with sedimentation experiments of highly uniform
colloidal cuboids synthesised by Yang and co-workers, which also preclude the
existence of stable $\rm N_{B}$ phases at equilibrium Yang  _et al._ (2018).
Research efforts have also been made to ascertain, by theory, simulation and
experiments, the existence of $\rm N_{B}$ phases in systems of other biaxial
particles Chiappini  _et al._ (2019); Tasios and Dijkstra (2017); van den Pol
_et al._ (2009); Peroukidis  _et al._ (2013); Peroukidis and Vanakaras (2013);
Peroukidis (2014); Skutnik  _et al._ (2020); Van den Pol  _et al._ (2010a, b);
Peroukidis  _et al._ (2020) or mixtures of uniaxial particles Cuetos  _et al._
(2008).

Despite such a widespread interest in mapping the phase behaviour of colloidal
suspensions of biaxial particles, including cuboids, the study of their
dynamics is still at an embryonic stage, especially for the difficulty of
finding suitable interaction potentials that could describe exotic geometries
and still be framed within a simulation technique. Exotic shapes are commonly
described by hard-core potentials, but these cannot be directly employed in
Brownian dynamics (BD) or Molecular Dynamics (MD) simulations. Nevertheless,
it is only by assessing the dynamics that one will be able to draw relevant
conclusions on the potential use of nematic or other LC phases in specific
applications. With this in mind, over the last years, our group has developed
a stochastic method that can qualitatively and quantitatively mimic the
Brownian motion of colloids as obtained by BD simulations Patti and Cuetos
(2012); Cuetos and Patti (2015); Corbett  _et al._ (2018); Daza  _et al._
(2020). This method, referred to as Dynamic Monte Carlo (DMC), has become an
established simulation technique not only for the study of the dynamics of
biaxial particles, such as cuboids and curved rods Chiappini  _et al._
(2020a); Cuetos and Patti (2020), but also for the dynamics of uniaxial
particles, like rods, for which soft potentials are indeed available Lebovka
_et al._ (2019); Chiappini  _et al._ (2020b).

With regards to the equilibrium dynamics of cuboids, we found that the system
long-time relaxation dramatically depends on particle anisotropy, being slower
at the self-dual shape, the geometry that would preferentially stabilise
biaxial nematics Cuetos and Patti (2020). By definition, the self-dual shape
is an intermediate geometry between prolate and oblate, where length ($L$),
width ($W$) and thickness ($T$) are such that $W=\sqrt{LT}$. Our simulations
also confirmed the occurrence of a Fickian and Gaussian dynamics at both short
and long times, thus providing an alternative picture to the claimed
universality of Fickian yet non-Gaussian dynamics in soft-matter systems
Cuetos  _et al._ (2018). For its potential impact in nanotechnology, equally
intriguing is the out-of-equilibrium dynamics of cuboids, especially because
it can spark phase switching and new material properties. In general, the
reorientation dynamics of biaxial particles induced by an external stimulus
has received very limited attention. At the molecular scale, Lee and co-
workers studied the reorientation dynamics of $\rm N_{B}$ phases of bent-core
mesogens and measured primary and secondary axis switching, finding the latter
either 3 or 100 times faster than the former depending on the mesogen Lee  _et
al._ (2007). Although this study was met with some scepticism Stannarius
(2008), Zannoni and co-workers later on performed MD simulations on $\rm
N_{B}$ phases of biaxial Gay-Berne ellipsoids and confirmed that the rotation
of minor axes is indeed faster, although only up to one order of magnitude,
than the rotation of the main axis, both in the bulk Berardi and Zannoni
(2008) and under confinement Ricci  _et al._ (2015). Following our recent
findings on the field-induced stability of the $\rm N_{B}$ phase Cuetos  _et
al._ (2019), here we explore the field-induced dynamics of switching from
uniaxial to biaxial nematics of colloidal cuboids, with special interest in
the particle reorientation dynamics and associated response time. More
specifically, we are interested to study the kinetics of reorientation of LCs
transitioning between two different nematic textures, namely an $\rm N_{U}$
$\rightarrow$ $\rm N_{B}$ transition under an external field, and an $\rm
N_{B}$ $\rightarrow$ $\rm N_{U}$ relaxation when the field is switched off,
and estimate the associated response times. To gain an insight into the impact
of particle anisotropy on the dynamics of phase switching, we consider
monodisperse systems of prolate, oblate and self-dual shaped cuboids. Their
ability of reorienting under the effect of an external field is assessed by
employing a DMC algorithm specifically designed to track the dynamics of out-
of-equilibrium colloidal systems Corbett  _et al._ (2018).

This paper is organised as follows. We first introduce the methodology to
simulate our systems and characterise the dynamics. We then discuss the
results of our simulations by analysing the effect of particle anisotropy on
the out-of-equilibrium dynamics in $\rm N_{U}$ $\rightarrow$ $\rm N_{B}$ and
$\rm N_{B}$ $\rightarrow$ $\rm N_{U}$ switching before finally drawing our
conclusions.

##  II Model and Simulation Methodology

We modelled monodisperse colloidal cuboids as hard board-like particles (HBPs)
constrained in a cubic box with periodic boundaries. The behaviour of hard-
core systems is basically determined by the packing fraction, which is given
by:

$\eta\equiv\frac{Nv_{o}}{V}$ (1)

where $N$ is the number of particles, $v_{o}$ the volume of an individual HBP
and $V$ the volume of the simulation box. The particle thickness, $T$, is set
as the system unit length. Consequently, particle length and width are given
in units of $T$, and read $L^{*}\equiv L/T$ and $W^{*}\equiv W/T$,
respectively. In particular, $L^{*}=12$ for all systems studied, while $W^{*}$
assumed values between 1 (rod-like HBPs) and 12 (plate-like HBPs) and included
$W^{*}=\sqrt{L^{*}}$ at the self-dual shape. Similar to our previous work
Cuetos  _et al._ (2019), we apply an external field that promotes alignment of
the particle intermediate axis, defined by:

$U_{\rm
ext}=-\frac{\varepsilon_{f}}{2}\Big{[}3\cdot(\rm\hat{\textbf{x}}_{i}\cdot\rm\hat{\textbf{e}})^{2}-1\Big{]}$
(2)

where $\varepsilon_{f}$ is the field strength, $\rm\hat{\textbf{x}}$ is the
unit vector associated with the width of particle $i$, while
$\rm\hat{\textbf{e}}$ is the field direction. We have set the reduced field
strength, $\varepsilon^{*}_{f}\equiv\varepsilon_{f}\beta=3$, with $\beta$ the
inverse temperature, which provides a measure of the relative strength of the
field applied with respect to thermal energy. The unit vectors
$\rm\hat{\textbf{y}}$ and $\rm\hat{\textbf{z}}$ are associated with thickness
and length, respectively. Similar to our previous work, we remind that this
external field model does not intend to mimic a real electric or magnetic
field, but rather, we focus on its effect to reorient particles. The
orientation of the particle unit vectors before and after application of the
field is schematically displayed in Fig. 1.

Figure 1: Reorientation of a prolate (a) and oblate (b) HBP due to external
field $\rm\hat{\textbf{e}}$ coupled to the particle intermediate axis
$\rm\hat{\textbf{x}}$. The $\rm N_{U}$ phases of self-dual shaped particles in
this work always have a director along $\hat{\textbf{n}}$, and therefore
rotate around $\hat{\textbf{z}}$ (similar to rod-like particles) when an
external field is applied.

The focus of this work is on the reorientation dynamics of HBPs in an external
field. We first performed standard Monte Carlo (MC) simulations in the
canonical ensemble in a cubic box containing $N=2000$ HBPs to equilibrate $\rm
N_{U}^{+}$ and $\rm N_{U}^{-}$ phases at $\eta=0.34$, where $\rm N_{U}^{+}$
and $\rm N_{U}^{-}$ refer, respectively to prolate and oblate nematic LCs. At
this packing fraction, the $\rm N_{U}$ phases (either prolate or oblate) are
stable across all anisotropies Cuetos  _et al._ (2017). Each MC cycle consists
of $N$ attempts to displace and/or rotate HBPs, which are accepted if no
overlaps are detected. To determine the occurrence of overlaps between pairs
of HBPs, we implemented the separating axes theorem by Gottschalk et al.
Gottschalk  _et al._ (1996), adapted by John and Escobedo to study tetragonal
parallelepipeds John and Escobedo (2005); John  _et al._ (2008). To quantify
the system long-range orientational order, we calculated the nematic order
parameter and director associated to each particle axis. To this end, we
performed the diagonalisation of a second-rank symmetric tensor of the form:

$\textbf{Q}^{\lambda\lambda}=\frac{1}{2N}\Biggl{\langle}\sum_{i=1}^{N}\big{(}3\hat{\lambda}_{i}\cdot\hat{\lambda}_{i}-\textbf{I}\big{)}\Biggr{\rangle}$
(3)

where $\hat{\lambda}_{i}$ = $\hat{\textbf{x}}$, $\hat{\textbf{y}}$,
$\hat{\textbf{z}}$ are unit vectors aligned with W, T and L respectively,
while I is the identity tensor. The diagonalisation of
$\textbf{Q}^{\lambda\lambda}$ results in three eigenvalues ($S_{2,W}$,
$S_{2,T}$, $S_{2,L}$) and their corresponding eigenvectors
($\rm\hat{\textbf{m}}$, $\rm\hat{\textbf{l}}$, $\rm\hat{\textbf{n}}$). The
nematic director of an $\rm N_{U}$ phase is the eigenvector with the largest
eigenvalue. For instance, if the largest positive eigenvalue is $S_{2,L}$,
then the nematic director is $\hat{\textbf{n}}$, indicating high degree of
orientational order along the $\hat{\textbf{z}}$ axis of the particles. The
biaxial order parameters can also be evaluated using the same symmetric
tensor. For example, the biaxial order parameter that quantifies the
fluctuations of particles’ axes $\hat{\textbf{x}}$ and $\hat{\textbf{y}}$
respectively along the directors $\hat{\textbf{m}}$ and $\rm\hat{\textbf{l}}$
reads Cuetos  _et al._ (2017)

$B_{2,L}=\frac{1}{3}\;(\hat{\bf m}\cdot\textbf{Q}^{xx}\cdot\hat{\bf
m}+\hat{\bf l}\cdot\textbf{Q}^{yy}\cdot\hat{\bf l}-\hat{\bf
m}\cdot\textbf{Q}^{yy}\cdot\hat{\bf m}-\hat{\bf
l}\cdot\textbf{Q}^{xx}\cdot\hat{\bf l})$ (4)

The values of $B_{2,W}$ and $B_{2,T}$ can be calculated using similar
expressions. Following the definition introduced in our former work, a phase
is considered to be biaxial if $B_{2}\geq 0.35$, although weak biaxial phases
can already be observed for $B_{2}\geq 0.20$ Cuetos  _et al._ (2019); Rafael
_et al._ (2020). We monitor the evolution of uniaxial and biaxial order
parameters until they have plateaued and fluctuate in a bounded range. The
equilibrated configurations are then used for external field application in
DMC simulations.

To study the dynamics, we performed DMC simulations in the canonical ensemble.
Because our goal is producing realistic time trajectories, unphysical moves,
such as cluster moves, swaps, jumps and changes in box dimension (which would
result in centres of mass rescaling) are not implemented. The position of the
particle centre of mass is updated by decoupling the displacement
$\delta\textbf{r}_{i}$ into three contributions, with
$\delta\textbf{r}_{i}=X_{W}\hat{\textbf{x}}+X_{T}\hat{\textbf{y}}+X_{L}\hat{\textbf{z}}$.
Rotational moves are performed by three consecutive reorientations around
$\hat{z}$, $\hat{x}$ and $\hat{y}$, with maximum rotations of $Y_{L}$, $Y_{W}$
and $Y_{T}$, respectively. The extent of particle displacement and rotation
are picked from uniform distributions that depend on the particle
translational, $D^{tra}_{\alpha,i}$, and rotational, $D^{rot}_{\alpha,i}$,
diffusion coefficients at infinite dilution, with $\alpha={L,W,T}$. Maximum
displacements and rotations are given by:

$\absolutevalue{X_{\alpha}}\leq\sqrt{2D^{tra}_{\alpha,i}\delta t_{MC}}$ (5)

$\absolutevalue{Y_{\alpha}}\leq\sqrt{2D^{rot}_{\alpha,i}\delta t_{MC}}$ (6)

where $\delta t_{MC}$ is the DMC timescale for one cycle, and is set to
$\delta t_{MC}=10^{-2}\tau$ for all simulations, with $\tau$ the time unit.
The coefficients $D^{tra}_{\alpha,i}$ and $D^{rot}_{\alpha,i}$ have been
estimated by using the open-source software HYDRO++ Carrasco and de la Torre
(1999); García de la Torre  _et al._ (2007). The interested reader is referred
to our previous work Cuetos and Patti (2020) for the specific values of these
translational and rotational diffusivities in units of $D_{o}\equiv
T^{2}\tau^{-1}$ and $D_{r}\equiv$ rad${}^{2}\tau^{-1}$ respectively. For a
monodisperse out-of-equilibrium system, the Brownian dynamics timescale,
$\delta t_{BD}$, can be obtained by rescaling the MC time scale as follows

$\delta t_{BD}=\frac{\mathcal{A}_{c}}{3}\delta t_{MC}$ (7)

where $\mathcal{A}_{c}$ is the time-dependent acceptance rate calculated at
the $c$th MC cycle over the transitory unsteady state Corbett  _et al._
(2018). We determine $\mathcal{A}_{c}$ by performing an MC cycle at a fixed
$\delta t_{MC}$ and integrating the above equation numerically:

$t_{\text{BD}}(\mathcal{C}_{\text{MC}})=\delta
t_{\text{MC}}\sum_{c=0}^{\mathcal{C}_{\text{MC}}}\frac{\mathcal{A}_{c}}{3}$
(8)

where $t_{\rm BD}(\mathcal{C}_{MC})$ is the Brownian time after
$\rm\mathcal{C}_{MC}$ MC cycles. It should be noted that the DMC method does
not consider solvent mediated hydrodynamic interactions, which are expected to
become especially relevant at strong external fields or large packing
fractions.

Figure 2: (a) Uniaxial and biaxial order parameters of a system with HBPs of
$W^{*}=1.0$ undergoing equilibration with an external field with
$\varepsilon^{*}_{f}=3$. The field is switched on at $t/\tau\approx 280$ and
switched off at $t/\tau\approx 550$. The shaded orange area corresponds to
$t_{\rm ON}$, while the green shaded area corresponds to $t_{\rm OFF}$. (b)
Response times ($t_{\rm ON}$ and $t_{\rm OFF}$) as a function of $W^{*}$ with
$\varepsilon^{*}_{f}=3$. The dashed vertical line in (b) represents the self-
dual shape that separates the prolate and oblate geometries.

To characterise the dynamics, we estimated (i) the response times, (ii) the
mean square angular displacement (MSAD), and (iii) the angular self-part of
the van-Hove function (s-VHF). We refer to the field-on ($t_{\rm ON}$) and
field-off ($t_{\rm OFF}$) response times as the time taken for the biaxial
order parameter to reach, respectively, 95$\%$ (field-on) and 105$\%$ (field-
off) of its equilibrium value. In particular, when a field is applied to an
$\rm N_{U}^{+}$ phase, $t_{\rm ON}$ is the time taken for $B_{2,L}$ to reach
95$\%$ of its equilibrium value in the field-on steady state. A schematic
illustration of how we performed this evaluation is reported in Fig. 2(a).
Both sets of response times have been calculated from an average over 50
independent trajectories per system. Approximately 2$\%$ of these trajectories
have given response times that were very different from those generally
observed. Since these anomalies tend to distort averages and give misleadingly
large error bars, we have considered them as outliers and excluded them from
the average. To this end, we employed the Modified Z-score method, a multiple
outlier rejection technique to identify statistical anomalies Iglewicz and
Hoaglin (1993). In particular, the Modified Z-score, $M_{j}$, is given by the
expression:

$M_{j}=\frac{0.6745\times\big{(}t_{j}-\bar{t}\big{)}}{MAD}$ (9)

where $t_{j}$ is the response time of trajectory $j$, $\bar{t}$ is the median
response time of the 50 trajectories and $MAD$ stands for median absolute
deviation. An observable is considered an outlier only if $M_{j}>3.5$ Iglewicz
and Hoaglin (1993).

The MSAD provides the ensemble average of the particle angular displacements
over time. To compute the MSAD, we employ a definition of an unbounded MSAD
akin to the translational mean square displacement. To this end, we first
introduce the definition of a rotational displacement vector which takes the
form Mazza  _et al._ (2006, 2007):

$\overrightarrow{\varphi}(t)=\int_{0}^{t}\delta\overrightarrow{\varphi}(t^{\prime})dt^{\prime}$
(10)

where $\delta\overrightarrow{\varphi}(t^{\prime})$ is a vector with direction
${\hat{\lambda}_{i}}(t^{\prime})\times\hat{\lambda}_{i}(t^{\prime}+dt^{\prime})$
and magnitude
$|\delta\overrightarrow{\varphi}(t^{\prime})|=\cos^{-1}[{\hat{\lambda}_{i}}(t^{\prime})\cdot\hat{\lambda}_{i}(t^{\prime}+dt^{\prime})]$.
From this, we can define the MSAD, which is mathematically expressed as:

$\langle\varphi^{2}(t)\rangle=\frac{1}{N}\Bigg{\langle}\sum\limits_{i=1}^{N}|{\overrightarrow{\varphi}}_{i}(t)-{\overrightarrow{\varphi}}_{i}(0)|^{2}\Bigg{\rangle}$
(11)

where $\overrightarrow{\varphi}_{i}$ is the rotational displacement vector of
particle $i$ defined in Eq. 10. Angular brackets denote average over different
trajectories. Finally, the so-defined rotational displacements are employed to
compute the angular s-VHF Mazza  _et al._ (2006, 2007):

$G(\varphi,t)=\frac{1}{N}\Biggl{\langle}\sum\limits_{i=1}^{N}\delta(\varphi-|\overrightarrow{\varphi}_{i}(t+t_{0})-\overrightarrow{\varphi}_{i}(t_{0})|)\Biggr{\rangle}$
(12)

where the symbol $\delta$ is the Dirac delta function. Basically,
$G(\varphi,t)$ provides the probability distribution of angular displacements
of particles within a time $t+t_{0}$ given their position at time $t_{0}$.

##  III Results

Figure 3: Schematic illustration of a field-induced $\rm N_{U}$ $\rightarrow$
$\rm N_{B}$ and a free $\rm N_{B}$ $\rightarrow$ $\rm N_{U}$ switching.
Different colours represent different orientations. Readers interested on the
detailed structure of the $\rm N_{U}$ and $\rm N_{B}$ phases are referred to
our past works for snapshots with higher resolution Cuetos  _et al._ (2019);
Rafael  _et al._ (2020).

Upon application of a sufficiently strong external field, an $\rm N_{U}$ phase
can be transformed into an $\rm N_{B}$ phase Cuetos  _et al._ (2019). This
transformation is not permanent, and, when the field is removed, uniaxiality
is restored. The time taken by the particles to reorient along the field
director measures the system’s ability of switching to a more ordered
configuration. Vice-versa, when the field is removed, the particles are left
free to rotate and the system recovers its original uniaxial state. A
schematic illustration of both transitory states is given in Fig. 3. We have
measured the response time associated to both $\rm N_{U}$ $\rightarrow$ $\rm
N_{B}$ and $\rm N_{B}$ $\rightarrow$ $\rm N_{U}$ transitions upon application
of the field $U_{\rm ext}$ with $\varepsilon^{*}_{f}=3$. The effect of
changing field intensity between $\varepsilon^{*}_{f}=1.5$ and 3 on the $\rm
N_{U}$ $\rightarrow$ $\rm N_{B}$ response time has also been assessed and is
available, for the interested reader, in Appendix A. The resulting response
times, $t_{\rm ON}$ and $t_{\rm OFF}$, are reported in Fig 2(b). Since the
formation of a field-induced $\rm N_{B}$ phase is dependent on the alignment
of the particle intermediate axis $\hat{x}$ with the external field, the
discussion that follows is relative to this axis, unless otherwise stated. To
start with, we notice that $t_{\rm ON}<t_{\rm OFF}$ across the complete set of
anisotropies (see Fig. 2(b)). In other words, at a given particle width, the
$\rm N_{U}$ $\rightarrow$ $\rm N_{B}$ switching is faster than the $\rm N_{B}$
$\rightarrow$ $\rm N_{U}$ switching. To understand the origin of this
behaviour, we calculated the MSADs of our systems and compared the field-on
and field-off profiles for each anisotropy.

The MSAD of systems with $W^{*}=2.5$ and 6 are shown, respectively, in the top
and bottom frames of Fig. 4. At very short times, the field-on and field-off
MSADs are very similar to each other, with the former becoming larger
immediately after and up to relatively long time scales. Over this period of
time, field-induced rotation is faster than free rotation. However, on time
scales comparable to $t_{\rm ON}$, a crossover between the two MSADs is
observed. On these time scales and beyond, free rotation grows significantly
much faster with time than field-induced rotation. We therefore conclude that
the presence of the external field accelerates the system orientational
dynamics by forcing the reorientation of the particle $\hat{x}$ axis along the
field director. As more and more HBPs are oriented, the field-on MSAD grows
less and less with time and would eventually saturate to a plateau if the
field strength was sufficiently large to offset and overcome thermal forces.
By contrast, the field-off MSAD practically shows the same behaviour with time
over the full time scale, as expected in free rotational diffusion. As for the
effect of anisotropy on the response time, we first discuss the case of the
field-induced uniaxial-to-biaxial transitory state.

Figure 4: MSAD in field-on and field-off scenarios of a system of HBPs with
reduced width (a) $W^{*}=2.5$ and (b) $W^{*}=6$. The field-on simulations
apply an external field of strength $\varepsilon^{*}_{f}=3$. The dashed
vertical lines indicate $t_{\rm ON}$ of each systems ($t_{\rm OFF}$ is out of
scale and not shown). The insets in (a) and (b) show the MSAD at shorter
timescales. Figure 5: MSAD of the $\hat{x}$, $\hat{y}$ and $\hat{z}$ axes of
HBPs for systems undergoing $\rm N_{U}\rightarrow N_{B}$ transition with
$\varepsilon^{*}_{f}=3$ for (a) prolate; (b) self-dual and (c) oblate HBPs,
and $\rm N_{B}\rightarrow N_{U}$ transition when the field is switched off for
a (d) prolate; (e) self-dual and (f) oblate HBPs. The dashed vertical lines in
each figure represent the $t_{\rm ON}$ for (a)-(c) and $t_{\rm OFF}$ for
(d)-(f).

In Fig. 2(b), $t_{\rm ON}$ increases with $W^{*}$, implying that the
reorientation is slower for oblate than for prolate particles. More
specifically, for rod-like HBPs ($W^{*}=1$), we observe a rapid switching with
$t_{\rm ON}/\tau\approx 41$, whereas for plate-like HBPs ($W^{*}=12$), it is
significantly slower, with $t_{\rm ON}/\tau\approx 4200$. Consequently, making
HBPs more oblate leads to a slower field-induced $\rm N_{U}\rightarrow N_{B}$
transition. To confirm these preliminary tendencies, we compare the MSADs of
the field-on regimes of each anisotropy along the three axes. The top frames
of Fig. 5 display the field-on MSADs of systems containing HBPs with
$W^{*}=1$, 3.46 and 12. We notice that the MSAD of the particle axis oriented
as the nematic director of the original $\rm N_{U}$ phase is the smallest
across all the geometries. More specifically, the MSAD of rod-like particles
in Fig. 5(a) exhibits a strong rotational coupling between $\hat{x}$ and
$\hat{y}$ particle axes, while $\hat{z}$ is practically unaffected by the
application of the field. Such a strong angular correlation between $\hat{x}$
and $\hat{y}$, with
$\langle\varphi^{2}_{W}\rangle=\langle\varphi^{2}_{T}\rangle$ over time, is
due to the square cross-sectional area of this specific set of HBPs, where
$W=T$. For similar reasons, plate-like HBPs with $W=L$ exhibit strong
rotational correlations between their axes $\hat{x}$ and $\hat{z}$, with
$\langle\varphi^{2}_{W}\rangle=\langle\varphi^{2}_{L}\rangle$ (see Fig. 5(c)),
while $\langle\varphi^{2}_{T}\rangle$, slightly increasing over time for mere
thermal fluctuations, remains very small, practically insensible to the
external field. In systems of self-dual shaped HBPs ($W^{*}=3.46$), we observe
that the MSADs of $W$ and $T$ are initially coupled, but then diverge over
time. This behaviour is observed for all anisotropies that are not perfectly
rod-like or plate-like and agrees very well with the tendencies reported in
our recent work on the equilibrium dynamics of HBPs Cuetos and Patti (2020).

When analysing the field-on MSADs of the particle axes perpendicular to the
original nematic director, we also notice an initially linear, rather steep
dependence on time, followed by an intermediate non-linear behaviour and
subsequently by a second linear regime at times comparable to $t_{\rm ON}$.
Such a long-time linear regime suggests that HBPs’ angular displacements are
gradually reducing, due to the system approaching a new equilibrium state.
Under these conditions, further rotations of the particle intermediate axis
$\hat{x}$, which is already aligned with the field, are suppressed, and only
small angular fluctuations are detected. At much larger field strengths, with
thermal fluctuations completely inhibited, we expect this second linear regime
to plateau at long times. In agreement with Fig. 2(b), we also observe that
systems with rod-like HBPs ($W^{*}=1$) only take $t/\tau\approx 17$ to reach
$<\varphi^{2}_{W}>=0.6$ rad2, whereas systems with self-dual shaped
($W^{*}=3.46$) or plate-like ($W^{*}=12$) HBPs take, respectively,
$t/\tau\approx 250$ and $t/\tau\approx 1600$ to achieve the same MSAD value.
This suggests that prolate HBPs tend to reorient significantly faster when an
external field is applied, leading to a relatively rapid equilibration.

Figure 6: (a) Angular s-VHFs of (a) field-on reorientation at time $t/\tau=35$
and (b) field-off reorientation at time $t/\tau=100$.

To gain a better insight into the dynamics of reorientation during this first
transitory unsteady state, we calculated the s-VHFs of all anisotropies at
$t/\tau=35$, corresponding to a time when all field-on cases are still
undergoing equilibration. As the MSAD for field-off scenarios are linear, we
arbitrarily picked $t/\tau=100$ to show the field-off s-VHFs. The s-VHFs shown
in Fig. 6 refer to the intermediate axis and have been normalised such that
$\int_{0}^{\infty}4\pi\varphi^{2}G(\varphi,t)d\varphi=1$. The first evident
conclusion, confirming the results discussed so far, is that prolate HBPs
rotate faster than oblate HBPs. This can be appreciated in Fig. 6(a) by
pinpointing the location of the peak of $G(\varphi_{W},t)$, which indicates
the most probable rotation achieved by particles of a given geometry at
$t/\tau=35$. In particular, the peak of $G(\varphi_{W=T},t)$ and
$G(\varphi_{W=L},t)$ suggests that rod-like and plate-like particles have
rotated, respectively, by $\varphi_{W}\approx 1.4$ rad and $\varphi_{W}\approx
0.2$ rad. By increasing particle width from $W^{*}=1$ to 12, the peak of the
angular s-VHFs gradually displaces towards lower rotations. Not only does the
particle anisotropy determine the location of the peak of these distributions
at a given time, but also their broadness. In other words, the angular s-VHFs
provide relevant information on the most probable rotation performed by HBPs
and on the existence of HBPs that rotate faster or slower than the average. In
particular, the presence of fast- and slow-responsive HBPs is evinced by the
tails of the distributions in Fig. 6(a), especially broad at $W^{*}=1$ and
then narrower and narrower up to $W^{*}=12$. Therefore, rod-like HBPs rotate
relatively fast, but heterogeneously (broad $G(\varphi_{W},t)$ peaked at large
distances), whereas plate-like HBPs are significantly slower, but rotate much
more homogeneously (narrow $G(\varphi_{W},t)$ peaked at short distances).

Figure 7: (a) Angular s-VHFs at various times, expressed as percentage of
$t_{\rm OFF}$, in systems of HBPs with (a) $W^{*}=3.46$ and (b) $W^{*}=4$.

When the field is switched off, the system recovers its original uniaxial
symmetry with the particles free to reorient under the mere effect of thermal
fluctuations. The time $t_{\rm OFF}$ taken by this field-off reorientation to
re-establish the $\rm N_{U}$ phase is again much shorter in nematics of
prolate HBPs (see Fig. 2(b)). In particular, $t_{\rm OFF}$ shows a tendency to
increase with particle width up to $W^{*}=4$, where $t_{\rm OFF}/\tau\approx
15\cdot 10^{3}$. When HBPs acquire a modest oblate geometry, such as from
$W^{*}=6$ to $W^{*}=8$, $t_{\rm OFF}$ drastically decreases to around $t_{\rm
OFF}/\tau\approx 6\cdot 10^{3}$, before increasing again at $W^{*}=12$. We
also notice that the free reorientation at $W^{*}=3.46$ and 4 is particularly
slower than that of other anisotropies and deserves an explanation. The MSADs
in the field-off scenarios, shown in Fig. 5(d)-(f), exhibit a linear profile
throughout the simulation due to the absence of an external field and decrease
upon increasing $W^{*}$. This tendency is also detected in the field-off
s-VHFs of Fig. 6(b), where, similarly to the field-on case, the angular
displacement decreases at increasing particle width. These elements would
suggest a scenario where $t_{\rm OFF}$ increases and free rotation becomes
slower upon increasing $W^{*}$. Because $t_{\rm OFF}$ is peaked for
approximately self-dual shaped particles and then decreases (see Fig. 2(b)),
there must be an additional element contributing to the field-free
reorientation from the $\rm N_{B}$ to the $\rm N_{U}$ phase. We believe that
this element is related to the ability of self-dual shaped HBPs of retaining
phase biaxiality when the field is switched off. In our recent work on the
field-induced phase behaviour of HBPs, we found that the self-dual shape
requires a surprisingly weak external field, compared to prolate and oblate
geometries, to spark an $\rm N_{U}\rightarrow N_{B}$ transition Cuetos  _et
al._ (2019). In particular, the minimum field strength to stabilise $\rm
N_{B}$ LCs was found to be $\varepsilon^{*}_{f}=0.1$ and 0.25 at $W^{*}=3.46$
and 4, respectively, and then increasing to $\varepsilon^{*}_{f}=0.5$ at
$W^{*}=3$ and to $\varepsilon^{*}_{f}=1$ at $W^{*}=6$. Therefore, we believe
that the field free $\rm N_{U}\rightarrow N_{B}$ transition at or very close
to the self-dual shape is affected by a metastability of the $\rm N_{B}$ phase
in off-field case. The underlying metastability allows the system to retain
the induced biaxiality over a longer time when compared with nematics of
oblate of prolate HBPs. In Fig. 7, we show the evolution of the s-VHFs of
$W^{*}=3.46$ and $W^{*}=4$ at different times up to $t_{\rm OFF}$. At short to
intermediate time scales, after having switched the field off, these s-VHFs
exhibit a double peak that suggests the presence of two populations of HBPs.
Because these two populations rotate at sufficiently different rates, we can
label them as slow and fast. The first peak survives over a relatively long
period of time, between 0.2$t_{\rm OFF}$ and 0.7$t_{\rm OFF}$, turning
gradually into a shoulder that disappears at longer times. These peaks and
subsequent shoulders are especially pronounced in the case of $W^{*}=4$ (Fig.
7(b)), explaining why $t_{\rm OFF}$ at $W^{*}=4$ is significantly slower than
$t_{\rm OFF}$ at $W^{*}=3.46$. Double peaks and shoulders are not observed at
other anisotropies or in field-on transitions (not shown here), indicating
that these tendencies are especially relevant only in the field-off relaxation
of self-dual shaped particles. In addition to their propensity towards
biaxiality retention, HBPs with $W^{*}=4$ rotate more slowly than perfectly
self-dual shaped particles, as shown in Fig. 6(b) and in agreement with the
tendencies observed in $\rm N_{U}$ phase in the absence of external fields
Cuetos and Patti (2020). The resulting large value of $t_{\rm OFF}$ is
therefore determined by the interplay between the particle’s ability to rotate
and the system’s tendency of retaining phase biaxiality. This interplay
explains the non-monotonic trend of $t_{\rm OFF}$ with particle shape in Fig.
2(b)) and provides, along with $t_{\rm ON}$, a useful guideline to select the
most suitable particle anisotropy for the design of field-responsive
nanomaterials.

##  IV Conclusions

In summary, by Dynamic Monte Carlo simulation, we studied the field-induced
dynamics in uniaxial nematic LCs of colloidal HBPs. By forcing the particles
to reorient around the nematic director, the external field induces an $\rm
N_{U}\rightarrow N_{B}$ phase transition that takes the system to a new steady
state. When the field is switched off, the biaxiality is gradually lost and
the $\rm N_{U}$ phase is restored. The time taken for the system to reorient,
also referred to as response time, strongly depends on the particle
anisotropy. The response times in $\rm N_{U}\rightarrow N_{B}$ and $\rm
N_{B}\rightarrow N_{U}$ switching were calculated and compared across all
anisotropies studied. Despite being the optimal shape to promote phase
biaxiality, the switching dynamics of self-dual shape HBPs is less
satisfactory compared to prolate HBPs. In particular, rod-like HBPs with
$W^{*}=1$ exhibit the fastest reorientation times in both the field-on and
field-off cases. The analysis of MSADs and s-VHFs show that the response time
is a result of a trade-off between particle rotational diffusion and phase
biaxiality retention, being both determined by shape anisotropy. Prolate HBPs
were found to rotate faster than self-dual shaped or oblate HBPs, allowing
rapid phase switching between the two nematic phases. Systems of HBPs with
geometry equal or very close to the self-dual shape exhibit a particularly
slow field-free reorientation, most likely due to relatively low field
strength required to transform $\rm N_{U}$ into $\rm N_{B}$ phases and
favouring the former in absence of an external field Cuetos  _et al._ (2019).
The ability of retaining biaxiality over a longer period of time is
corroborated by the existence of a double peak in the angular s-VHFs of $\rm
N_{B}\rightarrow N_{U}$ transition at short-to-intermediate time scales. This
double peak suggests the existence of two populations of (quasi) self-dual
shaped HBPs whose reorientation is not uniform and delays the system
relaxation. While prolate HBPs are especially field-responsive and exhibit a
rapid field-free reorientation, when one analyses the distribution of their
angular displacements over time, this appears to be very broad, with particles
exhibiting a very heterogeneous ability of rotating. By contrast, oblate HBPs,
while significantly less responsive, are characterised by a very narrow
distribution of angular displacements. All these elements offer a fundamental
understanding of the impact of shape anisotropy on the dynamics of uniaxial-
to-biaxial switching and a guidance to formulate nanomaterials with specific
switching dynamics for target applications.

##  V Acknowledgements

EMR would like to thank the Malaysian Government Agency Majlis Amanah Rakyat
for funding his PhD at the University of Manchester. AP and LT acknowledges
the financial support from the Leverhulme Trust Research Project Grant
RPG-2018-415. AC acknowledge the Spanish Ministerio de Ciencia, Innovación y
Universidades and FEDER for funding (project PGC2018-097151-B-I00). EMR, LT,
DC and AP acknowledge the assistance given by IT Services and the use of
Computational Shared Facility at the University of Manchester. Finally, we
thank Gerardo Campos-Villalobos (Utrecht University) for his assistance in
generating the snapshots in Fig. 3.

##  VI Data Availability

The data that support the findings of this study are available from the
corresponding author upon reasonable request.

##  Appendix A Effect of Field Strength

In this appendix, we briefly discuss the effect of altering field intensity on
the field-on response times, $t_{\rm ON}$ of HBPs. Here, we report the
response times for the field-on case for all anisotropies studied at field
strengths from $\varepsilon^{*}_{f}=1.5$ to $\varepsilon^{*}_{f}=3$. These
field intensities result in the formation of strong $\rm N_{B}$ phases with
$B_{2}\geq 0.35$ Cuetos  _et al._ (2019). The results are shown in Fig. 8.

Figure 8: Changes in field-on response time, $t_{\rm ON}$ as a function of
$\varepsilon^{*}_{f}$ across different anisotropies.

At constant $\varepsilon^{*}_{f}$, we observe that $t_{\rm ON}$ generally
increases with $W^{*}$ and this increment is significant. For instance, at
$\varepsilon^{*}_{f}=2$, the response time increases by two orders of
magnitude from $t_{\rm ON}\approx 57$ at $W^{*}=1$ to $t_{\rm ON}\approx 5700$
at $W^{*}=12$. We conclude that prolate particles tend to rotate faster than
oblate particles, regardless the field strength. Upon increasing field
strength, the reorientation becomes faster as suggested by the gradual
decrease of $t_{\rm ON}$ with $\varepsilon^{*}_{f}$. In addition, we note that
the statistical errors in $t_{\rm ON}$ decrease with $\varepsilon^{*}_{f}$,
most likely due to a stronger suppression of rotational fluctuations. This
tendency is consistent with the works by Zannoni and co-workers Berardi and
Zannoni (2008). For rod-like HBPs ($W^{*}=1$), increasing $\varepsilon_{f}$
does not significantly affect $t_{\rm ON}$, probably because the reorientation
capability of these HBPs is very close to its saturation value.

##  Appendix B Response Times

$W^{*}$ |  $t_{\rm ON}$ | $t_{\rm OFF}$
---|---|---
| $\varepsilon^{*}_{f}=1.5$ | $\varepsilon^{*}_{f}=2$ | $\varepsilon^{*}_{f}=2.5$ | $\varepsilon^{*}_{f}=3$ | $\varepsilon^{*}_{f}=3$
1 | $60.7\pm 10.0$ | $56.8\pm 2.9$ | $52.1\pm 5.2$ | 40.7 $\pm$ 3.6 | 75.0 $\pm$ 12.4
2.5 | $545.2\pm 100.9$ | $422.0\pm 64.4$ | $332.9\pm 40.7$ | 258.3 $\pm$ 29.4 | 878.6 $\pm$ 160.1
3 | $1057.5\pm 243.0$ | $746.4\pm 118.2$ | $548.4\pm 106.8$ | 420.9 $\pm$ 59.2 | 2181.0 $\pm$ 434.8
3.46 | $1764.9\pm 365.6$ | $1416.6\pm 247.5$ | $1130.0\pm 205.6$ | 837.6 $\pm$ 160.0 | 7552.2 $\pm$ 2758.4
4 | $4247.8\pm 683.3$ | $3012.6\pm 403.6$ | $2461.8\pm 273.4$ | 1948.8 $\pm$ 200.5 | 15004.6 $\pm$ 2011.0
6 | $3844.7\pm 641.8$ | $3380.6\pm 546.0$ | $2681.7\pm 345.7$ | 2272.2 $\pm$ 310.4 | 5975.5 $\pm$ 1200.2
8 | $4151.4\pm 846.8$ | $3728.3\pm 636.6$ | $3055.4\pm 272.5$ | 2700.0 $\pm$ 340.2 | 5602.0 $\pm$ 877.4
12 | $6084.4\pm 1376.9$ | $5700.9\pm 777.0$ | $4849.0\pm 598.6$ | 4175.5 $\pm$ 487.5 | 6906.0 $\pm$ 733.8
Table 1: Table of response times for $t_{\rm ON}$ ($\varepsilon^{*}_{f}=1.5$
to $3$) and $t_{\rm OFF}$ ($\varepsilon^{*}_{f}=3)$ with associated
statistical errors.

In Table I, we report $t_{\rm ON}$ for $\varepsilon^{*}_{f}=1.5$, 2, 2.5 and
3, and $t_{\rm OFF}$ for $\varepsilon^{*}_{f}=3$.

## References

  * Mormann  _et al._ (2017) W. Mormann, K.-H. Hellwich, J. Chen, and E. S. Wilks, Pure Appl. Chem. 89, 1695 (2017).
  * Brown (2015) R. Brown,  _The Miscellaneous Botanical Works of Robert Brown_, edited by J. J. Bennett, Cambridge Library Collection - Botany and Horticulture, Vol. 1 (Cambridge University Press, 2015).
  * Chen  _et al._ (2018) H.-W. Chen, J.-H. Lee, B.-Y. Lin, S. Chen, and S.-T. Wu, Light Sci. Appl. 7, 17168 (2018).
  * Huang  _et al._ (2018) Y. Huang, E. Liao, R. Chen, and S.-T. Wu, Appl. Sci. 8, 2366 (2018).
  * Dabrowski  _et al._ (2018) R. Dabrowski, J. Dziaduszek, J. Bozetka, W. Piecek, R. Mazur, M. Chrunik, P. Perkowski, M. Mrukiewicz, M. Żurowska, and D. Weglowska, J. Mol. Liq. 267, 415 (2018).
  * Meyer  _et al._ (2020) C. Meyer, C. Blanc, G. R. Luckhurst, P. Davidson, and I. Dozov, Sci. Adv. 6, eabb8212 (2020).
  * Cousins  _et al._ (2019) J. Cousins, S. Wilson, N. Mottram, D. Wilkes, and L. Weegels, Phys. Fluids 31, 083107 (2019).
  * Luckhurst and Sluckin (2015) G. R. Luckhurst and T. J. Sluckin,  _Biaxial nematic liquid crystals: theory, simulation and experiment_ (Wiley, 2015).
  * Lee  _et al._ (2007) J.-H. Lee, T.-K. Lim, W.-T. Kim, and J.-I. Jin, J. Appl. Phys. 101, 034105 (2007).
  * Berardi and Zannoni (2008) R. Berardi and C. Zannoni, J. Chem. Phys. 113, 5971 (2008).
  * Ricci  _et al._ (2015) M. Ricci, R. Berardi, and C. Zannoni, J. Chem. Phys. 143, 084705 (2015).
  * Taylor and Herzfeld (1991) M. P. Taylor and J. Herzfeld, Phys. Rev. A. 44, 3742 (1991).
  * Belli  _et al._ (2011) S. Belli, A. Patti, M. Dijkstra, and R. van Roij, Phys. Rev. Lett. 107, 148303 (2011).
  * Rafael  _et al._ (2020) E. M. Rafael, A. Cuetos, D. Corbett, and A. Patti, Soft Matter 16, 5565 (2020).
  * Dussi  _et al._ (2018) S. Dussi, N. Tasios, T. Drwenski, R. van Roij, and M. Dijkstra, Phys. Rev. Lett. 120, 177801 (2018).
  * Belli  _et al._ (2012) S. Belli, M. Dijkstra, and R. van Roij, J. Phys.: Condens. Matter 24, 284128 (2012).
  * Leferink op Reinink  _et al._ (2014) A. B. G. M. Leferink op Reinink, S. Belli, R. van Roij, M. Dijkstra, A. V. Petukhov, and G. J. Vroege, Soft Matter 10, 446 (2014).
  * Cuetos  _et al._ (2019) A. Cuetos, E. M. Rafael, D. Corbett, and A. Patti, Soft Matter 15, 1922 (2019).
  * Cuetos  _et al._ (2017) A. Cuetos, M. Dennison, A. Masters, and A. Patti, Soft Matter 13, 4720 (2017).
  * Patti and Cuetos (2018) A. Patti and A. Cuetos, Mol. Simul. 44, 516 (2018).
  * Yang  _et al._ (2018) Y. Yang, G. Chen, S. Thanneeru, J. He, K. Liu, and Z. Nie, Nat. Commun. 9, 4513 (2018).
  * Chiappini  _et al._ (2019) M. Chiappini, T. Drwenski, R. Van Roij, and M. Dijkstra, Phys. Rev. Lett. 123, 068001 (2019).
  * Tasios and Dijkstra (2017) N. Tasios and M. Dijkstra, J. Chem. Phys. 146, 144901 (2017).
  * van den Pol  _et al._ (2009) E. van den Pol, A. Petukhov, D. Thies-Weesie, D. Byelov, and G. Vroege, Phys. Rev. Lett. 103, 258301 (2009).
  * Peroukidis  _et al._ (2013) S. D. Peroukidis, A. G. Vanakaras, and D. J. Photinos, Phys. Rev. E 88, 062508 (2013).
  * Peroukidis and Vanakaras (2013) S. D. Peroukidis and A. G. Vanakaras, Soft Matter 9, 7419 (2013).
  * Peroukidis (2014) S. D. Peroukidis, Soft Matter 10, 4199 (2014).
  * Skutnik  _et al._ (2020) R. A. Skutnik, I. S. Geier, and M. Schoen, Mol. Phys. , 1 (2020).
  * Van den Pol  _et al._ (2010a) E. Van den Pol, A. Lupascu, P. Davidson, and G. Vroege, J. Chem. Phys. 133, 164504 (2010a).
  * Van den Pol  _et al._ (2010b) E. Van den Pol, A. Lupascu, M. Diaconeasa, A. Petukhov, D. Byelov, and G. Vroege, Chem. Phys. Lett. 1, 2174 (2010b).
  * Peroukidis  _et al._ (2020) S. D. Peroukidis, S. H. Klapp, and A. G. Vanakaras, Soft Matter 16, 10667 (2020).
  * Cuetos  _et al._ (2008) A. Cuetos, A. Galindo, and G. Jackson, Phys. Rev. Lett. 101, 237802 (2008).
  * Patti and Cuetos (2012) A. Patti and A. Cuetos, Phys. Rev. E 86, 011403 (2012).
  * Cuetos and Patti (2015) A. Cuetos and A. Patti, Phys. Rev. E 92, 022302 (2015).
  * Corbett  _et al._ (2018) D. Corbett, A. Cuetos, M. Dennison, and A. Patti, Phys. Chem. Chem. Phys. 20, 15118 (2018).
  * Daza  _et al._ (2020) F. A. G. Daza, A. Cuetos, and A. Patti, Phys. Rev. E 102, 013302 (2020).
  * Chiappini  _et al._ (2020a) M. Chiappini, A. Patti, and M. Dijkstra, Phys. Rev. E 102, 040601 (2020a).
  * Cuetos and Patti (2020) A. Cuetos and A. Patti, Phys. Rev. E 101, 052702 (2020).
  * Lebovka  _et al._ (2019) N. I. Lebovka, N. V. Vygornitskii, and Y. Y. Tarasevich, Phys. Rev. E 100, 042139 (2019).
  * Chiappini  _et al._ (2020b) M. Chiappini, E. Grelet, and M. Dijkstra, Phys. Rev. Lett. 124, 087801 (2020b).
  * Cuetos  _et al._ (2018) A. Cuetos, N. Morillo, and A. Patti, Phys. Rev. E 98, 042129 (2018).
  * Stannarius (2008) R. Stannarius, J. Appl. Phys. 104, 034105 (2008).
  * Gottschalk  _et al._ (1996) S. Gottschalk, M. C. Lin, and D. Manocha, Comp. Graph. 30, 171 (1996).
  * John and Escobedo (2005) B. S. John and F. A. Escobedo, J. Phys. Chem. B 109, 23008 (2005).
  * John  _et al._ (2008) B. S. John, C. Juhlin, and F. A. Escobedo, J. Chem. Phys. 128, 044909 (2008).
  * Carrasco and de la Torre (1999) B. Carrasco and J. G. de la Torre, Biophys. J. 76, 3044 (1999).
  * García de la Torre  _et al._ (2007) J. García de la Torre, G. del Rio Echenique, and A. Ortega, J. Phys. Chem. B 111, 955 (2007).
  * Iglewicz and Hoaglin (1993) B. Iglewicz and D. C. Hoaglin,  _How to detect and handle outliers_ , Vol. 16 (Asq Press, 1993).
  * Mazza  _et al._ (2006) M. G. Mazza, N. Giovambattista, F. W. Starr, and H. E. Stanley, Phys. Rev. Lett. 96, 057803 (2006).
  * Mazza  _et al._ (2007) M. G. Mazza, N. Giovambattista, H. E. Stanley, and F. W. Starr, Phys. Rev. E 76, 031203 (2007).