# In search of a precursor for crystal nucleation of hard and charged colloids

Marjolein de Jager<EMAIL_ADDRESS>Soft Condensed Matter, Debye Institute
of Nanomaterials Science, Utrecht University, Utrecht, Netherlands  Frank
Smallenburg  Université Paris-Saclay, CNRS, Laboratoire de Physique des
Solides, 91405 Orsay, France  Laura Filion  Soft Condensed Matter, Debye
Institute of Nanomaterials Science, Utrecht University, Utrecht, Netherlands

###### Abstract

The interplay between crystal nucleation and the structure of the metastable
fluid has been a topic of significant debate over recent years. In particular,
it has been suggested that even in simple model systems such as hard or
charged colloids, crystal nucleation might be foreshadowed by significant
fluctuations in local structure around the location where the first nucleus
arises. We investigate this using computer simulations of spontaneous
nucleation events in both hard and charged colloidal particles. To detect
local structural variations, we use both standard and unsupervised machine
learning methods capable of finding hidden structures in the metastable fluid
phase. We track numerous nucleation events for the face-centered cubic and
body-centered cubic crystal on a local level, and demonstrate that all signs
of crystallinity emerge simultaneously from the very start of the nucleation
process. We thus conclude that there is no precursor for the nucleation of
charged colloids.

††preprint: APS/123-QED

##  I Introduction

Crystal nucleation plays an important role in fields ranging from colloidal
self-assembly, to protein crystallization, and even polymorph selection in
pharmaceuticals Karthika, Radhakrishnan, and Kalaichelvi (2016); Sosso  _et
al._ (2016). However, despite its importance in a number of essential fields,
the detailed mechanism of forming a crystal nucleus still remains a topic of
continuous debate.

The simplest theory which addresses crystal nucleation is classical nucleation
theory (CNT). In CNT, the metastable fluid is continuously undergoing thermal
fluctuations, where small, solid clusters form and dissolve until one appears
which is large enough (critically large) to grow out into a macroscopic
crystal. The size of such a critical cluster is given simply by balancing the
bulk free-energy gain associated with transitioning into the more stable solid
phase, with the surface free-energy cost of having a finite crystal cluster
immersed in the fluid. This picture, however, becomes significantly more
complicated when one considers the possibility of multiple competing crystal
structures, typically referred to as polymorphs. In systems with crystal
polymorphs, the crystalline phase that first nucleates in the metastable fluid
is not necessarily the stable phase. Theories to address such situations, such
as the Ostwald step rule Ostwald (1897), and the Alexander-McTague theory
Alexander and McTague (1978) have proven unreliable in explaining polymorph
selection (see e.g. Refs. Russo and Tanaka, 2012a; Taffs and Patrick Royall,
2016; Ouyang  _et al._ , 2016).

One complication when studying such questions is the close interplay between
local structural motifs that occur naturally in the fluid, and the ones that
might emerge when the crystal forms. It has been suggested that motifs hiding
in the fluid are predictive of, or even responsible for, the location or
polymorph of the nucleus that forms Russo and Tanaka (2012b); Gispen  _et al._
(2023). To investigate this possibility, one avenue forward could be to
explore just how much information the metastable liquid is hiding regarding
the nucleation process. Over the last two decades a plethora of studies have
appeared presenting contradictory observations ten Wolde and Frenkel (1999);
Schilling  _et al._ (2010); Russo and Tanaka (2016); Tan, Xu, and Xu (2014);
Li  _et al._ (2020); Russo and Tanaka (2012b); Hu and Tanaka (2022); Lu  _et
al._ (2015); Berryman  _et al._ (2016); Tanaka  _et al._ (2019); Lechner,
Dellago, and Bolhuis (2011). In particular, in simple systems such as charged
colloids which nucleate into either the face-centered cubic (FCC) or body-
centered cubic (BCC) crystal, some studies have argued that local structural
order develops before the local density increases Russo and Tanaka (2012b); Hu
and Tanaka (2022); Lu  _et al._ (2015), while other authors found evidence
that the two processes happen simultaneously Berryman  _et al._ (2016).

To address this issue, some recent, intriguing studies have explored how
modifying (via biasing) the structure of the fluid – either enhancing or
suppressing specific local motifs – affects the nucleation process Hu and
Tanaka (2022); Taffs and Patrick Royall (2016). In principle, such studies
might be able to give one direct evidence that a specific local structure
either enhances or suppresses the nucleation process. Unfortunately, however,
biasing the structure of the fluid modifies not only its local structure but
also its thermodynamics, meaning that comparisons with the unbiased case are
inconclusive.

The more direct route to trying to explore how various kinds of local ordering
interplay in crystal nucleation is simply to simulate the nucleation event,
and follow the various structural and density features as nucleation happens.
At first glance this would appear to be a straightforward approach. However,
the challenge in this case lies in the difficulty in creating local order
parameters that are unbiased. For example, order parameters that are tuned to
recognize the crystalline regions from fluid might struggle at the boundary
between the fluid and crystal – a highly important aspect at the beginning of
nucleation. Similar issues exist for other order parameters, making it very
difficult to pinpoint the start of the nucleation process and hence to
determine whether structural order emerges before, during, or after
densification. Hence, in some cases instead of capturing accurately whether
local structure exists in the highly fluctuating metastable fluid, one ends up
examining the properties of the order parameter instead of the properties of
the fluid.

While this problem is never fully avoidable, one option to try and avoid
accidental biases is to exploit multiple different measures for local order –
for instance measures associated with symmetries like bond-order parameters
and order associated with the topological connections between neighboring
particles – such as topological cluster classification (TCC). Interestingly,
new unsupervised machine learning (UML) algorithms also give new avenues to
probe structure (see e.g. Refs. Becker  _et al._ , 2022; Reinhart  _et al._ ,
2017; Boattini, Dijkstra, and Filion, 2019; Boattini  _et al._ , 2020; Coli
and Dijkstra, 2021; van Damme  _et al._ , 2020; Gardin  _et al._ , 2021;
Coslovich, Jack, and Paret, 2022; Paret, Jack, and Coslovich, 2020; Adorf  _et
al._ , 2019) . Recent studies have even demonstrated that simple, UML-based
approaches are able to extract variations in disorder in the structure of
supercooled fluids from e.g. a simple vector of bond order parameters Boattini
_et al._ (2020); Paret, Jack, and Coslovich (2020); Coslovich, Jack, and Paret
(2022). Intriguingly, this includes identifying variations in local structure
that are not easily extracted by looking at each element of the vector
individually.

In this paper, we attempt to take the utmost care in identifying local
signatures of the fluid and and revisit the question: are there hidden local
structures present in metastable fluid that foreshadow the location of the
imminent formation of a crystal nucleus? Specifically, we apply both classical
and UML-based methods to the nucleation of hard and charged colloids in both
the regime of strong screening and weak screening, for which respectively the
FCC and BCC crystals nucleate. To this end, we simulate numerous spontaneous
nucleation events, and closely follow all nucleation events as a function of
time. In particular, similar to Ref. Berryman  _et al._ , 2016 we zoom in on
the regions where the nuclei are born and analyze the local fluctuations in
density and structure of the metastable fluid. By doing this we can locally
track whether there is a delay between the increase in local structural
ordering and local density prior to the start of nucleation. Such a delay
would indicate the presence of a precursor. However, within the limits of this
study, we find no evidence of such a precursor in the systems we studied.

##  II Model

We consider a system of $N$ like-charged hard spheres of diameter $\sigma$
suspended in a solvent containing salt. The effective interaction potential
between these colloids is given by the repulsive hard-core Yukawa potential

$\beta\phi(r)=\begin{dcases}\beta\epsilon\;\frac{e^{-\kappa\sigma(r/\sigma-1)}}{r/\sigma}&\quad\text{for
}r\geq\sigma,\\\ \infty&\quad\text{for }r<\sigma,\end{dcases}$ (1)

with contact value
$\beta\epsilon=Z^{2}\lambda_{B}/\sigma(1+\kappa\sigma/2)^{2}$, where $Z$ is
the charge of the colloids in electron charge, $\lambda_{B}$ is the Bjerrum
length, $\kappa$ is the inverse Debye screening length, and $\beta=1/k_{B}T$,
with $k_{B}$ the Boltzmann constant and $T$ the temperature. Note that in the
limit of zero charge ($Z\to 0$) or infinite screening
($\kappa\sigma\to\infty$), this potential reduces to the hard-sphere
potential. The interaction potential was truncated and shifted such that the
shift was never more than $10^{-5}k_{B}T$.

Nucleation of both the BCC and FCC phases in this system has been studied in
the past (see e.g. Refs. Auer and Frenkel, 2005; Desgranges and Delhommelle,
2007; Browning, Doherty, and Fredrickson, 2008; Gispen and Dijkstra, 2022; de
Jager and Filion, 2022). In a previous study de Jager and Filion (2022), we
used umbrella sampling to calculate the nucleation barriers and rates of
highly screened charged particles. In this paper, we will study the nucleation
of some of these (nearly-)hard systems, as well as the nucleation of weakly
screened charged particles. To be able to compare the nucleation processes of
different systems, we select state points with approximately equal barrier
heights. In particular, we will simulate brute-force nucleation events of
systems with barrier heights around 15-18$k_{B}T$. Information on the
nucleation barriers of the systems studied is given in Tab. 1. Note that
systems with a Debye screening length of $1/\kappa\sigma=0.01$ were found to
behave essentially as “hard” spheres when mapped with an effective hard-sphere
diameter de Jager and Filion (2022). A brief explanation of the methods used
for computing the nucleation barriers as well as some additional information
on these systems can be found in the Supplemental Materials (SM).

Table 1:  For each system studied, the packing fraction of the supersaturated fluid $\eta^{*}$ at which the brute force nucleation is performed together with the corresponding supersaturation $\beta|\Delta\mu|$. The last columns give the critical nucleus size $n^{*}$ and barrier height $\beta\Delta G^{*}$ obtained using umbrella sampling. The error in $\beta\Delta G^{*}$ is no more than 1.  | $\beta\epsilon$ | $1/\kappa\sigma$ | $\eta^{*}$ | $\beta|\Delta\mu|$ | $n^{*}$ | $\beta\Delta G^{*}$
---|---|---|---|---|---|---
FCC | hard spheres | 0.5385 | 0.585 | 75 | 16.5
81 | 0.01 | 0.4681 | 0.584 | 84 | 16.3
8 | 0.04 | 0.4400 | 0.541 | 69 | 14.8
BCC | 81 | 0.40 | 0.1305 | 0.321 | 122 | 18.0

##  III Methods

To explain the methods we use for studying the nucleation events, we need to
discuss two things: i) how we identify local structure, and ii) how we track
nucleation events locally.

###  III.1 Identifying local structure

We use three different methods to classify the local structure. The first
method considers just the averaged bond-orientational order parameters (BOPs)
of Lechner and Dellago Lechner and Dellago (2008). For this, we first
calculate for each particle $i$ the complex quantities

$q_{lm}(i)=\frac{1}{N_{b}(i)}\sum_{j\in\mathcal{N}_{b}(i)}Y_{l}^{m}(\theta_{ij},\phi_{ij}),$
(2)

where $\mathcal{N}_{b}(i)$ is the set of the $N_{b}(i)$ nearest neighbors of
particle $i$, $\,Y_{lm}\left(\theta,\phi\right)$ are the spherical harmonics
with $m\in[-l,l]$, and $\theta_{ij}$ and $\phi_{ij}$ are the polar and
azimuthal angles of the vector $\mathbf{r}_{ij}=\mathbf{r}(j)-\mathbf{r}(i)$
connecting particles $i$ and $j$. We use the SANN algorithm van Meel  _et al._
(2012) to determine the nearest neighbors. Next, we average these complex
quantities over the set of nearest neighbors as well as the particle itself

$\bar{q}_{lm}(i)=\frac{1}{N_{b}(i)+1}\sum_{j\in\\{i,\mathcal{N}_{b}(i)\\}}q_{lm}(j).$
(3)

Finally, we compute the rotationally invariant averaged BOPs

$\bar{q}_{l}(i)=\sqrt{\frac{4\pi}{2l+1}\sum_{m=-l}^{l}|\bar{q}_{lm}(i)|^{2}},$
(4)

and

$\bar{w}_{l}(i)=\frac{w_{l}(i)}{\left(\sum_{m=-l}^{l}|\bar{q}_{lm}(i)|^{2}\right)^{3/2}},$
(5)

with

$w_{l}(i)=\\!\\!\\!\\!\\!\\!\sum_{\begin{matrix}\scriptstyle
m_{1},m_{2},m_{3}\\\ \scriptstyle
m_{1}+m_{2}+m_{3}=0\end{matrix}}\\!\\!\\!\\!\\!\\!\begin{pmatrix}l&l&l\\\
m_{1}&m_{2}&m_{3}\end{pmatrix}q_{lm_{1}}(i)q_{lm_{2}}(i)q_{lm_{3}}(i),$ (6)

where the term in brackets is the Wigner $3j$ symbol, which is only non-zero
when $m_{1}+m_{2}+m_{3}=0$. Note that $w_{l}=0$ when $l$ is odd. Depending on
the choice of $l$, these BOPs are sensitive to different (crystal) symmetries.
For example, $\bar{q}_{6}$ is very helpful in distinguishing more fluid-like
environments from more solid-like environments such as FCC and BCCLechner and
Dellago (2008).

While the BOPs are extremely useful for detecting specific symmetries in the
local structure of the fluid, they are not necessarily optimal for detecting
the most important structural variations in disordered systems such as the
metastable fluid. Recent work has shown that BOPs in combination with
unsupervised machine learning algorithms is highly effective at autonomously
detecting variations that might be difficult to see by studying the individual
BOPs Boattini, Dijkstra, and Filion (2019); Boattini  _et al._ (2020). Hence,
for the second method, we use an unsupervised algorithm to autonomously detect
local structural fluctuations in the metastable fluid. In Ref. Boattini,
Dijkstra, and Filion, 2019, Boattini et al. showed that a neural-network-based
autoencoder, which is given $\bar{q}_{l}$ with $l\in[1,8]$ as input, does a
good job in distinguishing a whole range of different local structures.
Similarly, in Ref. van Damme  _et al._ , 2020 van Damme et al. showed that
using principal component analysis (PCA) as dimensionality reduction method
also does a good job in distinguishing the sizable assortment of crystal
structures formed by rounded tetrahedra. For our system we found that PCA and
an autoencoder preformed equally well in distinguishing order, and hence chose
to use the simpler PCA algorithm in our analysis. To specifically focus on
finding local fluctuations or signatures in the metastable fluid, we train the
PCA model on $\bar{q}_{l}$ with $l\in[1,8]$ of configurations containing only
fluid particles and no (significant) solid nuclei. We then use this trained
PCA model to analyze the entire nucleation trajectory. Note that in this way,
the first principal component corresponds to the largest BOPs-related
structural variation in the metastable fluid.

For the third and last method, we use an altogether different approach for
classifying local structure. In particular we use the topological cluster
classification (TCC) algorithm developed by Malins et al. Malins  _et al._
(2013) to detect any local motifs that the BOPs might have overlooked. More
specifically, we use TCC to calculate the population of certain types of
clusters, as well as the number of clusters of a certain type a particle is
involved in.

Note that for computing the nucleation barriers, we additionally need a binary
classification method which labels a particle as either fluid or solid. For
this we use the 6-fold Ten Wolde bonds ten Wolde, Ruiz-Montero, and Frenkel
(1996)

$d_{6}(i,j)=\frac{\sum_{m=-6}^{m=6}q_{6m}(i)q_{6m}^{*}(j)}{\sqrt{\left(\sum_{m=-6}^{m=6}|q_{6m}(i)|^{2}\right)\left(\sum_{m=-6}^{m=6}|q_{6m}(j)|^{2}\right)}},$
(7)

where ∗ indicates the complex conjugate and $q_{6m}$ are the bond-
orientational order parameters given by Eq. (2). Particle $i$ is classified as
solid if it has 6 or more neighboring particles $j$ with which it has a solid-
like bond, i.e. $d_{6}(i,j)>0.7$. We also use this fluid-solid classification
to initially locate and follow the nucleation event. However, we want to point
out that, although this provides a general overview of the nucleation event,
it is not an ideal order parameter to study the onset of nucleation. Its
binary nature with the thresholds for $d_{6}(i,j)$ and the number of solid-
like bonds causes it to overlook subtle increases in the local structural
ordering of the fluid and hence reacts more slowly to the nucleation than
other continuous order parameters.

###  III.2 Simulating and tracking nucleation events

To obtain the spontaneous nucleation events we use brute force MC, KMC, and MD
simulations in the $NVT$-ensemble, where we simulate $N=10976$ and $N=11664$
particles for the systems forming the FCC and BCC phase, respectively. The
difference between the MC and KMC simulations is the acceptance ratio of the
trial particle moves. For the MC simulations this acceptance ratio is around
30%, whilst for the KMC simulations it is around 85% resulting in dynamics
that mimic Brownian motion Sanz and Marenduzzo (2010). The MD simulations are
performed using LAMMPS with a Nose-Hoover thermostat Plimpton (1995). We only
perform MD simulations for the soft hard-core Yukawa system, i.e. the one with
$1/\kappa\sigma=0.40$. As this system does not feel its hard core, it can
simply be ignored. Once we have obtained the numerous nucleation events, we
want to track the local density and structure to determine if there is a
difference between the increase in local structural ordering and local density
at the start of nucleation. To this end, for each nucleation event we find the
position $\mathbf{r}_{0}$ that best captures the center of the nucleus at the
start of nucleation. For this we use the average center-of-mass of the
precritical nucleus as a starting point and, if needed, by eye adjust it to
best capture the birthplace of the crystal nucleus. Next, for each snapshot of
the nucleation trajectory, starting well before the start of nucleation, we
determine all particles inside a sphere of radius $R$ around $\mathbf{r}_{0}$,
and take the average of the local properties of these particles. This is
similar to what Berryman et al. did in Ref. Berryman  _et al._ , 2016. The
local structural properties that we consider are explained in the previous
subsection. Additionally, we define for each particle a local packing fraction
measured via the volume of its Voronoi cell. The volumes of the Voronoi cells
were obtained using voro++ Rycroft (2009). As we are searching for local
precursors, we choose $R$ such that the selected region contains around 30-40
particles. This size provides a good balance between being large enough to
obtain relatively stable averages of the local properties, and being small
enough to ensure that the averaged properties still represent the local
situation.

##  IV Results

###  IV.1 Structure of the metastable fluid

Before we look into the actual crystal nucleation, we first characterize the
structural properties of the metastable fluid.

First, we examine the globally averaged values of the local BOPs, and plot the
results in Fig. 1 for the metastable fluids of essentially hard spheres and of
soft spheres as a function of the supersaturation. We see that $\bar{q}_{6}$,
$\bar{q}_{8}$, and $\bar{w}_{8}$ are most prominent in both metastable fluids,
and that all BOPs are only marginally affected by the increase in
supersaturation. Furthermore, notice that the values in both systems are
surprisingly similar, even though the metastable fluid of essentially hard
spheres later forms an FCC crystal, whereas the fluid of the soft spheres will
form a BCC crystal. The most prevalent difference between the two systems can
be found in $\bar{w}_{6}$, which is smaller for the nearly-hard spheres than
for soft spheres, and for high supersaturation even becomes on average
negative for the nearly-hard spheres whereas it stays positive for the soft
spheres. See the SM for more analysis on the $\bar{w}_{l}$’s. We, thus,
conclude that the fluid’s “knowledge” about which crystal phase it should
nucleate into is difficult to distinguish from the global values of the BOPs.

In addition to the BOPs, we take a look at the presence of the different TCC
clusters in the metastable fluids. Figure 2 shows the population of various
TCC clusters in the metastable fluid of hard spheres, and soft spheres. Even
though we see some small deviations in the populations of the different
metastable fluids – e.g. clusters 6A, 8A, 8K, 9K, and BCC_9 have a slightly
higher population in the fluid of soft spheres and 9B, 10B, 11C, 11E, and 12D
have a slightly higher population in the fluid of nearly-hard spheres – the
values are again surprisingly similar. This indicates once more that it is
difficult to determine which crystal phase will nucleate from the metastable
fluid for the systems studied here.

a) |  | b)
---|---|---
|  |
c) |  | d)
|  |
Figure 1:  The mean of a,b) the first eight $\bar{q}_{l}$’s and c,d) the first four even $\bar{w}_{l}$’s as a function of the supersaturation for the metastable fluids of hard-core Yukawa with a,c) $\beta\epsilon=81$ and $1/\kappa\sigma=0.01$, and b,d) $\beta\epsilon=81$ and $1/\kappa\sigma=0.40$.  Figure 2:  The population of various TCC clusters in the metastable fluids of hard spheres ($\eta=0.5385$), nearly-hard hard-core Yukawa particles ($\beta\epsilon=8$, $1/\kappa\sigma=0.04$, $\eta=0.4400$), and soft hard-core Yukawa particles ($\beta\epsilon=81$, $1/\kappa\sigma=0.40$, $\eta=0.1305$).  Figure 3:  PCA on the metastable fluids of a,b) hard spheres ($\eta=0.5385$) and c,d) soft hard-core Yukawa particles ($\beta\epsilon=81$, $1/\kappa\sigma=0.40$, $\eta=0.1305$). a,c) give the weight of each BOP in the first and second principal component. b,d) show the distribution of the fluid particles in the PC1-PC2 plane, with the addition of the corresponding bulk FCC phase of hard spheres ($\eta=0.5981$) and bulk BCC phase of soft hard-core Yukawa ($\beta\epsilon=81$, $1/\kappa\sigma=0.40$, $\eta=0.1311$). e) Pearson correlation between PC1 and the number of TCC clusters a particle is involved in for the hard spheres. The color of the bars indicate clusters that consist of one or more tetrahedral subclusters (blue), one or more square pyramidal subclusters (yellow), or both/neither (gray).  a) $\eta$ | b) PC1 | c) 9B | d) 11F
---|---|---|---
|  |  |
|  |  |
Figure 4:  Snapshot of the metastable fluid of hard spheres ($\eta=0.5385$)
colored by a) the local packing fraction, b) the first principal component, c)
the number of 9B clusters per particle, and d) the number of 11F clusters per
particle.

Next, we characterize the local ordering of the metastable fluid on the
single-particle level. As explained in the Methods, we train a PCA model using
only configurations of the metastable fluid. In all cases the first principal
component (PC1) explains around 70% of the total variance of the input. To
illustrate what kind of fluctuations PCA picks up in the fluid, Figs. 3a,c)
show for the metastable fluids of hard spheres and of soft spheres the weight
of each BOP in the first and second principal component. We see that PC1 is
mostly made up of $\bar{q}_{6}$ and $\bar{q}_{8}$. Furthermore, Figs. 3b,d)
show the distribution of these metastable fluid particles in the PC1-PC2
plane, as well as the distributions of the corresponding FCC and BCC phase.
(Recall that the data for the crystalline particles was not used in training
the PCA models.) In this scatter plot we can neatly see that the crystal
phases lie in the region of large PC1. To get a better understanding of the
real-space distribution of these particles with above or below average PC1, we
take a look at a single snapshot of the metastable fluid of hard spheres and
color the particles according to their local packing fraction, PC1, and the
number of 9B and 11F clusters a particles is involved in, see Fig. 4. Even
though the spatial correlations in the local packing fraction are not clearly
visible, we can clearly distinguish by eye large spatial regions of above or
below average PC1. The autocorrelation functions of these spatial correlations
can be found in the SM. Notice that the regions with above average PC1
correspond to an absence of 9B clusters and a high presence of 11F clusters,
while regions with below average PC1 correspond to a high presence of 9B
clusters and an absence of 11F clusters. Thus, there is a negative correlation
between PC1 and 9B clusters and positive correlation between PC1 and 11F. The
precise correlations of these two and other TCC clusters with PC1 are shown in
Fig. 3e). Analogous to what was found in Ref. Boattini  _et al._ , 2020, the
TCC clusters can be roughly divided into two groups: those with a negative
correlation, which essentially are all clusters consisting of one or more
tetrahedral subclusters, and those with a positive correlation, which contain
the clusters consisting of one or more square pyramidal subclusters.

a) |  | b)
---|---|---
|  |
Figure 5:  Autocorrelation function of the first and second principal components and the local packing fraction in the metastable fluids of a) hard spheres ($\eta=0.5385$) and b) soft hard-core Yukawa particles ($\beta\epsilon=81$, $1/\kappa\sigma=0.40$, $\eta=0.1305$). The different dashing and darkness of the color indicate the simulation method, and the time is in terms of the long-time diffusion time $\tau_{d}$.  a) $\;\quad(t-t_{0})/\tau_{d}=-0.18$ | b) $\;\quad(t-t_{0})/\tau_{d}=0.40$ | c) $\;\quad(t-t_{0})/\tau_{d}=0.96$ | d) $\;\quad(t-t_{0})/\tau_{d}=1.43$
---|---|---|---
|  |  |
Figure 6:  Four snapshots of a typical nucleation event of hard spheres
($\eta=0.5385$). Here $t_{0}$ indicates the start of nucleation. Fluid
particles are displayed at a quarter of their actual size to make the nucleus
visible, and red indicates the particles inside the studied region, i.e. those
inside the sphere of radius $R$ around the center of nucleation
$\mathbf{r}_{0}$.

Combining the observation of large spatial regions of above average PC1 with
the proximity of these particles to the crystal phases in the PC1-PC2 scatter
plot (Fig. 3), we can conclude that regions of above average PC1 form a good
candidate for harboring a precursor for crystal nucleation. In the next
section, we will investigate these regions whilst tracking the nucleation
events. However, before we turn our attention to that, we need to determine
the temporal correlations of the local structure such that we know the time
window before the start of nucleation during which we can search for a
precursor. Figure 5 shows, for multiple simulation methods, the
autocorrelation functions (ACFs) of the first two principal components and the
local packing fraction in the metastable fluids of hard spheres and of soft
spheres. Here, we give the time in terms of the long-time diffusion time
$\tau_{d}=\sigma^{2}/6D_{l}$, where $D_{l}$ is the long-time diffusion
coefficient obtained from the mean-squared displacement. Notice that the ACFs
are essentially independent of the simulation method, which confirms that the
dynamics are also independent of the choice of simulation method. We see in
both systems that PC2 and the local packing fraction decay extremely fast in
time. Although PC1 decays more slowly, i.e. within half a diffusion time for
the hard spheres and one diffusion time for the soft spheres, this is still
relatively fast. The decay time of PC1 provides a good estimate for the time
window before the start of nucleation in which we can search for a precursor.

a) |  | b)
---|---|---
|  |
c) |  | d)
|  |
e) |  | f)
|  |
g) |  | h)
|  |
Figure 7:  Left: typical nucleation event of hard spheres ($\eta=0.5385$), same as in Fig. 6. Right: typical nucleation event of soft hard-core Yukawa particles ($\beta\epsilon=81$, $1/\kappa\sigma=0.40$, $\eta=0.1305$). Both events were obtained using MC simulations. The vertical dashed line in each figure indicates the start of nucleation $t_{0}$, and the shaded area indicates the time window before $t_{0}$ for which the ACF of PC1 $>0.05$ (see Fig. 5). In a-f) the black line (right axis) gives the size of the biggest nucleus present in the studied region. The other lines give the average value of a-b) $\bar{q}_{l}$ for $l\in[3,8]$, c-d) $\bar{w}_{l}$ for $l\in[4,6,8]$, e-f) PC1. In g-h) the blue line (right axis) gives PC1, while the yellow line (left axis) gives the local packing fraction $\eta_{\text{local}}$. Note that in e-g) the horizontal dashed lines give the reference value of PC1 in the fluid and solid phase. In g) the right axis is scaled in such a way that the reference values in the fluid and solid of PC1 and $\eta_{\text{local}}$ lie on top of each other. In h) the horizontal dashed lines give the reference value of $\eta_{\text{local}}$ in the fluid and solid phase.  a) |  | b)
---|---|---
|  |
c) |  | d)
|  |
e) |  | f)
|  |
g) |  | h)
|  |
Figure 8:  For the same events as in Fig. 7, i.e. hard spheres (left) and soft
hard-core Yukawa particles (right), the average number of clusters per
particle (left axis, yellow) for a couple of TCC clusters together with the
average value of PC1 (right axis, blue). The horizontal dashed lines indicate
the reference values of the number of clusters per particle and PC1 in the
fluid. In a-f) the left axis is scaled in such a way that these lines lie on
top of each other. In g-h) this was not possible without inverting one of the
axes.

###  IV.2 Nucleation study

We now turn our attention to crystal nucleation. As explained in the methods,
we simulate numerous spontaneous nucleation events using (K)MC and MD
simulations, and track for all these events the local properties of the region
where nucleation starts. Here, we discuss our observations using two typical
nucleation events: one of the hard-spheres system and one of soft spheres
($\beta\epsilon=81$, $1/\kappa\sigma=0.40$). Both these nucleation events were
obtained using MC simulations. More nucleation events, where we either used
other simulation methods or studied the other systems mentioned in Tab. 1, can
be found in the SM. To better illustrate which region we study while tracking
a nucleation event, Fig. 6 shows a couple of snapshots of the nucleation event
of hard spheres where the particles inside the studied region are colored red.
Figure 7 shows for this event and the nucleation event of soft spheres the
average properties of the particles in this studied region. Before we discuss
what we see, let us again point out that the nucleus size (black line) is not
an ideal order parameter for tracking nucleation since its binary nature
causes it to overlook subtle increases in the local structural ordering at the
onset of nucleation. It does, however, provide a general overview of the
nucleation event, such as when the nucleus reaches its critical size (see Tab.
1). That being said, let us first discuss the BOPs of the studied region. We
observe no notable change in the behavior of the BOPs before the start of
nucleation, but, as soon as nucleation starts, we see a sharp increase in the
values of $\bar{q}_{6}$ and $\bar{q}_{8}$ for both systems. Furthermore, for
the hard spheres, we see that, once nucleation starts, $\bar{q}_{4}$
increases, $\bar{w}_{6}$ stays negative, and $\bar{w}_{4}$ decreases. This all
indicates that indeed the FCC phase nucleates. On the other hand, for the soft
spheres, we see that as nucleation starts $\bar{q}_{4}$ barely increases,
$\bar{w}_{6}$ stays positive, and $\bar{w}_{4}$ keeps fluctuating around zero,
which all indicates that indeed the BCC phase nucleates. Similar to the
behavior of the BOPs, we observe no notable change in the behavior of PC1
before the start of nucleation, but see a sharp increase in its value once
nucleation starts. Note that this increase in PC1 is visible before the number
of particles classified as crystalline starts to rise (black line),
demonstrating that PC1 is a better order parameter for tracking the start of
nucleation than the nucleus size according to our definition. Lastly, for the
hard spheres, we see that the local packing fraction increases simultaneously
with PC1 as soon as nucleation starts, and that no notable behavior can be
observed before the start of nucleation. This strongly indicates that increase
in structural ordering and local density go together and, thus, that there is
no apparent precursor. Unfortunately, as the difference between the packing
fraction of the fluid and solid phases is extremely small for the soft
spheres, i.e. less than 0.001, it is not possible to observe any increase in
the local packing fraction on top of the normal fluctuations. Hence, we cannot
draw any conclusions on the local packing fraction of the soft spheres.

Next, to show that we have missed no subtle changes in the local structure and
thus confirm that there is no precursor for nucleation, we further examine the
local structure using TCC. In Fig. 8 we show for four of the most relevant TCC
clusters the average number of clusters a particle is involved in and compare
it with PC1. These four clusters are: i) 6A, which has the strongest positive
correlation with PC1 and is present in both bulk FCC and bulk BCC, ii) 8A,
which has the second strongest positive correlation with PC1 and is present in
bulk BCC but not bulk FCC, iii) FCC, which is present in bulk FCC but not bulk
BCC, and iv) 7A, which has the strongest negative correlation with PC1 and can
neither be found in bulk FCC nor bulk BCC. Similar figures of other TCC
clusters can be found in the SM. For all clusters we see that there is no
significant change prior to the start of nucleation. Furthermore, we see that
the trends of the 6A cluster coincide almost perfectly with those of PC1.
Similarly, we see that the trends of the 8A cluster closely follow the trends
of PC1. However, for hard spheres the initial increase in 8A clusters is
followed by a decrease. Since 8A is a cluster that is usually found in bulk
BCC and not in bulk FCC, this initial increase might be surprising. This can
be explained via the observation that 8A clusters are found in high
concentrations near the surface of growing nucleiGispen  _et al._ (2023). As a
result, the number of these clusters decreases once the nucleus grows beyond
our averaging radius. For the FCC cluster, we observe a sharp increase during
the nucleation of hard spheres. Notice, however, that this increase starts
slightly later than the increase in PC1. This is not surprising as the FCC
cluster is a relatively large cluster, i.e. it contains 13 particles, and
consequently is not present in the first stages of nucleation. For the soft
spheres there is no significant increase in FCC clusters, as expected. Lastly,
we take a look at the 7A cluster. In contrast to the other three clusters,
this five-fold symmetric cluster has a strong negative correlation with PC1.
Moreover, it is strongly present in the metastable fluid phases, whereas its
presence in the FCC and BCC phases is negligible. It is, therefore, not
surprising that we observe an immediate and sharp decrease in 7A clusters as
soon as nucleation starts.

##  V Conclusions

To conclude, we have characterized the local structure of various metastable
fluids of charged colloids using multiple methods: bond-orientational order
parameters (BOPs), principal component analysis (PCA) on the BOPs, and
topological cluster classification (TCC). In doing this have attempted to
avoid artefacts due to biases in our chosen order parameters. For all systems
we have found that any local structural ordering has a relatively short
lifetime, resulting in a short time window prior to the start of nucleation in
which a precursor could exist. By tracking the local structure of the spatial
region coinciding with the birthplace of the crystal nucleus, we show that
inside this time window no atypical behavior in the local structural order is
observed using any of our structural order parameters. Furthermore, we
demonstrate that all structural characteristics that differ significantly
between the fluid and crystal phase start changing simultaneously as soon as
nucleation starts. Specifically in the case of FCC, this includes the local
density, which starts growing immediately as soon as structural order emerges.
We, thus, conclude that we find no evidence for a precursor for the crystal
nucleation of hard and charged colloids.

##  VI Acknowledgements

L.F. and M.d.J. acknowledge funding from the Vidi research program with
project number VI.VIDI.192.102 which is financed by the Dutch Research Council
(NWO).

## References

  * Karthika, Radhakrishnan, and Kalaichelvi (2016) S. Karthika, T. Radhakrishnan, and P. Kalaichelvi, Crystal Growth & Design 16, 6663 (2016).
  * Sosso  _et al._ (2016) G. C. Sosso, J. Chen, S. J. Cox, M. Fitzner, P. Pedevilla, A. Zen, and A. Michaelides, Chemical reviews 116, 7078 (2016).
  * Ostwald (1897) W. Ostwald, Zeitschrift für physikalische Chemie 22, 289 (1897).
  * Alexander and McTague (1978) S. Alexander and J. McTague, Phys. Rev. Lett. 41, 702 (1978).
  * Russo and Tanaka (2012a) J. Russo and H. Tanaka, Soft Matter 8, 4206 (2012a).
  * Taffs and Patrick Royall (2016) J. Taffs and C. Patrick Royall, Nat. Commun. 7, 13225 (2016).
  * Ouyang  _et al._ (2016) W. Ouyang, C. Fu, Z. Sun, and S. Xu, Phys. Rev. E 94, 042805 (2016).
  * Russo and Tanaka (2012b) J. Russo and H. Tanaka, Sci. Rep. 2, 1 (2012b).
  * Gispen  _et al._ (2023) W. Gispen, G. M. Coli, R. van Damme, C. P. Royall, and M. Dijkstra, ACS Nano 17, 8807–8814 (2023).
  * ten Wolde and Frenkel (1999) P. R. ten Wolde and D. Frenkel, Physical Chemistry Chemical Physics 1, 2191 (1999).
  * Schilling  _et al._ (2010) T. Schilling, H. J. Schöpe, M. Oettel, G. Opletal, and I. Snook, Phys. Rev. Lett. 105, 025701 (2010).
  * Russo and Tanaka (2016) J. Russo and H. Tanaka, J. Chem. Phys. 145, 211801 (2016).
  * Tan, Xu, and Xu (2014) P. Tan, N. Xu, and L. Xu, Nat. Phys. 10, 73 (2014).
  * Li  _et al._ (2020) M. Li, Y. Chen, H. Tanaka, and P. Tan, Science advances 6, eaaw8938 (2020).
  * Hu and Tanaka (2022) Y.-C. Hu and H. Tanaka, Nat. Commun. 13, 4519 (2022).
  * Lu  _et al._ (2015) Y. Lu, X. Lu, Z. Qin, and J. Shen, Solid State Communications 217, 13 (2015).
  * Berryman  _et al._ (2016) J. T. Berryman, M. Anwar, S. Dorosz, and T. Schilling, J. Chem. Phys. 145, 211901 (2016).
  * Tanaka  _et al._ (2019) H. Tanaka, H. Tong, R. Shi, and J. Russo, Nature Reviews Physics 1, 333 (2019).
  * Lechner, Dellago, and Bolhuis (2011) W. Lechner, C. Dellago, and P. G. Bolhuis, Phys. Rev. Lett. 106, 085701 (2011).
  * Becker  _et al._ (2022) S. Becker, E. Devijver, R. Molinier, and N. Jakse, Phys. Rev. E 105, 045304 (2022).
  * Reinhart  _et al._ (2017) W. F. Reinhart, A. W. Long, M. P. Howard, A. L. Ferguson, and A. Z. Panagiotopoulos, Soft Matter 13, 4733 (2017).
  * Boattini, Dijkstra, and Filion (2019) E. Boattini, M. Dijkstra, and L. Filion, J. Chem. Phys. 151, 154901 (2019).
  * Boattini  _et al._ (2020) E. Boattini, S. Marín-Aguilar, S. Mitra, G. Foffi, F. Smallenburg, and L. Filion, Nat. Commun. 11, 1 (2020).
  * Coli and Dijkstra (2021) G. M. Coli and M. Dijkstra, ACS nano 15, 4335 (2021).
  * van Damme  _et al._ (2020) R. van Damme, G. M. Coli, R. van Roij, and M. Dijkstra, ACS nano 14, 15144 (2020).
  * Gardin  _et al._ (2021) A. Gardin, C. Perego, G. Doni, and G. M. Pavan, arXiv (2021).
  * Coslovich, Jack, and Paret (2022) D. Coslovich, R. L. Jack, and J. Paret, J. Chem. Phys. 157, 204503 (2022).
  * Paret, Jack, and Coslovich (2020) J. Paret, R. L. Jack, and D. Coslovich, J. Chem. Phys. 152, 144502 (2020).
  * Adorf  _et al._ (2019) C. S. Adorf, T. C. Moore, Y. J. Melle, and S. C. Glotzer, J. Phys. Chem. B 124, 69 (2019).
  * Auer and Frenkel (2005) S. Auer and D. Frenkel, Adv. Comput. Simul. 173, 149 (2005).
  * Desgranges and Delhommelle (2007) C. Desgranges and J. Delhommelle, The Journal of chemical physics 126, 054501 (2007).
  * Browning, Doherty, and Fredrickson (2008) A. R. Browning, M. F. Doherty, and G. H. Fredrickson, Physical Review E 77, 041604 (2008).
  * Gispen and Dijkstra (2022) W. Gispen and M. Dijkstra, Physical Review Letters 129, 098002 (2022).
  * de Jager and Filion (2022) M. de Jager and L. Filion, J. Chem. Phys. 157, 154905 (2022).
  * Lechner and Dellago (2008) W. Lechner and C. Dellago, J. Chem. Phys. 129, 114707 (2008).
  * van Meel  _et al._ (2012) J. A. van Meel, L. Filion, C. Valeriani, and D. Frenkel, J. Chem. Phys. 136, 234107 (2012).
  * Malins  _et al._ (2013) A. Malins, S. R. Williams, J. Eggers, and C. P. Royall, J. Chem. Phys. 139, 234506 (2013).
  * ten Wolde, Ruiz-Montero, and Frenkel (1996) P. R. ten Wolde, M. J. Ruiz-Montero, and D. Frenkel, Faraday Discuss. 104, 93 (1996).
  * Sanz and Marenduzzo (2010) E. Sanz and D. Marenduzzo, J. Chem. Phys. 132, 194102 (2010).
  * Plimpton (1995) S. Plimpton, J. Comp. Phys. 117, 1 (1995).
  * Rycroft (2009) C. Rycroft, Chaos 19, 041111 (2009).


# Order-disorder charge density wave instability in the kagome metal
(Cs,Rb)V3Sb5

D. Subires  Donostia International Physics Center (DIPC), San Sebastián, Spain
A. Korshunov  European Synchrotron Radiation Facility (ESRF), BP 220, F-38043
Grenoble Cedex, France  A. H. Said  Advanced Photon Source, Argonne National
Laboratory, Lemont, IL 60439  L. Sánchez  Donostia International Physics
Center (DIPC), San Sebastián, Spain  Brenden R. Ortiz  Materials Department
and California Nanosystems Institute, university of California Santa Barbara,
Santa Barbara, California, 93106, USA  Stephen D. Wilson  Materials Department
and California Nanosystems Institute, university of California Santa Barbara,
Santa Barbara, California, 93106, USA  A. Bosak  European Synchrotron
Radiation Facility (ESRF), BP 220, F-38043 Grenoble Cedex, France  S. Blanco-
Canosa<EMAIL_ADDRESS>Donostia International Physics Center (DIPC), San
Sebastián, Spain  IKERBASQUE, Basque Foundation for Science, 48013 Bilbao,
Spain

(November 2022)

###### Abstract

##  I Abstract

The origin of the charge density wave phases in the kagome metal compound
AV3Sb5 is still under great scrutiny. Here, we combine diffuse and inelastic
x-ray scattering to identify a 3-dimensional precursor of the charge order at
the L point that condenses into a CDW through a first order phase transition.
The quasi-elastic critical scattering indicates that the dominant contribution
to the diffuse precursor is the elastic central peak without phonon softening.
However, the inelastic spectra show a small broadening of the Einstein-type
phonon mode on approaching TCDW. Our results point to the situation where the
Fermi surface instability at the L point is of order-disorder type with
critical growth of quasi-static domains. The experimental data indicate that
the CDW consists on an alternating Star of David and trihexagonal distortions
and its dynamics goes beyond the classical weak-coupling scenario and is
discussed within strong-electron phonon coupling and non-adiabatic models.

##  II Introduction

The intricate relationship between lattice geometry and topological electronic
behavior determines the ground state properties of materials. The non-trivial
band topology of the kagome lattice is being extensively explored as
candidates to engineer Dirac fermions Ye18 , topological flat bands Kan20 ,
magnetic Weyl semimetals Liu18  or quantum spin liquid behavior Yan11 ; Han12
. Besides, the combination between the non-trivial topology and strong
electronic correlations, originated from the interplay of orbital, charge and
spin degrees of freedom, introduces topological correlated physics Mazin2014
with more exotic phenomena, like topological Majorana modes.

Recently, a new family of non magnetic kagome metals AV3Sb5 (A=K, Rb, Cs) has
emerged as a fertile playground to investigate correlation-driven topological
phases Ort19 ; Ort20 ; Ort21 . Derived from the kagome structure, this family
of materials features a quasi 2D electronic band structure with van Hove
singularities (vHs) and Dirac crossings close to the Fermi level Kan22 ; Lou22
. Superconductivity sets in at low temperature and is intertwined with an
exotic charge density wave phase (CDW) Zhao2021 . Although the electronic
structure obtained from density functional theory (DFT) calculations is well
established, controversy surrounds the origin and stabilization of the CDW,
and the mechanism of superconductivity and its gap symmetry Zhao21 ; Tan21 ;
Wu21 .

The rich phase diagram of the kagome lattice was investigated with respect to
Fermi surface instabilities within the Hubbard model, revealing spin and CDW
orders and unconventional superconductivity emerging upon tuning the band
filling to the vHs Kie13 ; Wang13 . Diffraction experiments and scanning
tunneling microscopy (STM) unraveled multiple-q CDW modulations with an
inverse Star of David structure (SoD) and a trihexagonal (TrH) distortions
Miao22 ; Jian21 , following the phonon instabilities at the M and L points
predicted by DFT Tan21 ; Subedi22 . The q vector at M connects neighbouring
vHs in the vicinity of the Fermi surface of AV3Sb5, suggesting that the CDW is
triggered by the Peierls mechanism Kan22 ; Lou22 . On the other hand, the CDW
of CsV3Sb5 has been shown to break the time reversal symmetry Mie22 ,
unveiling an entanglement of superconductivity and charge order chirality
Jian21 . Nevertheless, the chiral CDW has been questioned by recent STM
experiments Li22 , pointing to a sample/surface dependent fermiology or the
presence of dynamical disorder. This is particularly appealing, since the
sample quality and A-site occupancy has a strong impact on the Tc Ort21  and
might also influence the possible origin and stabilization of the charge
modulations. Further, recent experiments reported a new unidirectional
modulation charge stripe order with 4a0 wavelength that breaks the rotational
symmetry Li22 , nematicity Nie22  and pair density wave order Chen21 .

Central to the interplay between CDW and $\mathbb{Z}_{2}$ band topology is the
nature of the CDW phase. From the thermodynamic point of view, specific heat
measurements Ort20 ; Ort21  show a small release of entropy through the CDW
transition, suggesting a first order phase transition, and tight binding
models describe the CDW as electronically mediated Den21 ; Fen21 . In
contrast, DFT predicts soft modes at M and L points of the Brillouin zone
Tan21 , for which, however, inelastic x-ray scattering (IXS) does not observe
any anomaly of the low energy acoustic branches Li22Miao . However, ARPES data
have reported deviations from the linear electronic dispersion Luo22 ,
providing that the electron-phonon coupling gives sufficient strength for the
CDW formation. Furthermore, optical spectroscopy ascribed the opening of the
CDW gap to the suppression of density of states at the M point Zhou21  with
strong phonon renormalizations Uyk22  without electronic anomalies at L.
Finally, in addition to the enigmatic driving force of the CDW, the multiple
qCDW and the different energy scales Nak21  hinder a detailed study of the
ground state of AV3Sb5, and the origin and stabilization of the CDW phase
transition remains unclear.

Here, we use a combination of diffuse scattering (DS) and inelastic x-ray
scattering (IXS) to identify frozen phonon fluctuations of the
2$\times$2$\times$2 (L point) 3D-CDW in AV3Sb5 at T$>$TCDW that condense into
a triple-QL CDW through a first order phase transition, preserving the six-
fold symmetry of the kagome lattice. The precursor of the charge order is
absent for the 2$\times$2 (M point) and 2$\times$2$\times$4 superstructures,
indicating that the 2$\times$2$\times$2 charge modulation is the leading Fermi
surface instability. Furthermore, the intensity of the quasi-elastic central
peak (CP) of IXS follows the diffuse scattering (DS) intensity, showing that
the CDW is of an order-disorder transformation type. This is corroborated by
the absence of a clear softening, but small broadening of the optical branch
at L. Our results solve the dispute about the primary order parameter and
highlight the critical role of the order-disorder phase transition to
understand the controversial experimental reports.

Figure 1: Charge density wave in CsV3Sb5. (a) Chemical structure of AV3Sb5.
Orange, red and cyan balls stand for A (Cs, Rb), V and Sb atoms. (b) Top view
of the kagome structure, showing the V-kagome net interlaced with the Sb
atoms. (c) Calculated phonon dispersions of fully-relaxed CsV3Sb5. Soft modes
are observed at M and L. The red ellipse marks the region where phonons are
measured in figure 4(d). (d) Sketch of the Brillouin zone along the main
symmetry directions. (e) Temperature dependence of the (2.5 0.5 0.5) CDW
reflection in CsV3Sb5, showing a sharp drop of intensity at the TCDW=95 K.
(Inset, rocking scans). Figure 2: Diffuse scattering maps of CsV3Sb5. (a) (h k
0.5) plane at 90 K (top), showing the CDW reflections (weak spots highlighted
within the yellow circle in the zoomed in area) and the precursor of the 3D
CDW at 100 K (bottom), marked with arrows. (b) and (c) (h k 0) and (h 0 l)
planes, respectively (no precursor). The DS maps of RbV3Sb5 are shown in the
Supplementary Figure 1.

##  III Results

###  III.1 x-ray diffraction

Figure 1(a) and (b) sketches the structure of CsV3Sb5 (P6/mmm space group).
The unit cell is described as individual sublattices of a kagome network of V
coordinated by Sb1 atoms and a hexagonal net of Sb2 above and below each
kagome layer. The lattice parameters obtained from the DS at room temperature
(a=b=5.49 Å, c=9.34 Å) are consistent with the reports in the literature Ort19
; Ort20 . The temperature dependence of the CDW instability at the L point is
displayed in the figure 1(e), showing a sharp onset at $\sim$95 K, indicative
of a first order-like phase transition. Bragg reflections with the same
temperature dependence (not shown) were also measured with propagation vectors
(0.5 0 0) and (0.5 0 0.25) r.l.u. Li22Miao .

###  III.2 Diffuse Scattering, DS

Besides the sharp and intense CDW diffraction spots in x-ray diffraction,
diffuse intensity is present at T$>$TCDW due to lattice vibrations and
correlated disorder. Figure 2 summarizes the processed DS maps for the CsV3Sb5
sample at 90 and 100 K for the (h k 0.5), (h k 0) and (h 0 l) planes. The
image reconstruction of the (h k 0.5) plane at 100 K (figure 2(a) bottom)
reveals the presence of DS corresponding to freezing of a transverse component
of phonons (streaks highlighted in the yellow square), with most pronounced
intensity around the (3 0 0.5) reciprocal lattice vector. Other satellites are
equivalent according to the symmetry operations of the P6/mmm space group. The
DS is visible below 150 K and condense into a 3QL order at TCDW (figure 2(a)
top). The fluctuating precursor observed in figure 2(a) along the $<$1 0 1$>$
direction is in agreement with the phonon calculations that predict a lattice
instability at the L point, figure 1(c) Tan21 ; Subedi22 . The correlation
length of the charge precursor ($\xi$) perpendicular to the streak of diffuse
intensity remains nearly constant as we cool down to TCDW and extends to 4-5
unit cells, while along the streak reduces to less than one unit cell, 3(b).
Remarkably, DS is absent for T$>$TCDW in the (h k 0) plane (figure 2(b)
bottom) and (h 0 l) plane (figure 2 (c) bottom) that maps the longitudinal
component of the (0.5 0 0.25) and (0.5 0 0.5) CDWs Li22Miao ; Li22 ; Jian21 .
However, CDW reflections are again present at T$<$TCDW, matching the hard
x-ray diffraction data of Fig. 1 (e). The presence of diffuse intensity as
satellites at only the (0.5 0 0.5) reciprocal lattice vector indicates that
the Fermi surface instability at the L position is the leading order parameter
and drives the CDW formation in the kagome metal AV3Sb5 (see Supplementary
Figure 1 for the DS maps of RbV3Sb5).

Figure 3: Temperature dependence of the diffuse scattering and IXS central
peak (a) 90 K, (b) 100 K and (c) 150 K around the (3 0 0.5) cut. The sketch in
(c) shows the reciprocal lattice vectors and the IXS scan direction through
the CDW, as indicated in the red line. (d) Profile of the diffuse scattering
at 120 K and its fitting to a sum of Bragg and Lorentzian profiles. (e)
Temperature dependence of the IXS spectra at the (2.5 0.5 0.5) reciprocal
lattice vector, highlighting the enhancement of the elastic line and drop of
intensity of the low-energy phonon at 96 K. (f) Temperature dependence of the
intensity of the DS signal, elastic central peak and optical phonon. Inset,
temperature dependence of the generalized susceptibility. T∗ denotes the
change of slope of the electronic susceptibility, $\sim 105$ K, presumably
related to percolation of CDW domains. The error bars are the fit to
uncertainty.

The reconstruction of the (h k 0.5) plane for CsV3Sb5 at different
temperatures around the (3 0 0.5) cut is shown in figure 3(a-c). Figure 3(d)
shows the profile of the diffuse precursor along the diffuse streak of the
figure 3(b). The critical scattering was modelled assuming a sum of Bragg and
a Lorentzian profile that accounts for the peak sharpening as T approaches
TCDW from the high temperature. The temperature dependence of the DS intensity
is extracted by integrating the diffuse precursor and plotted in figure 3(f).
It starts at temperatures below 150 K and increases down to 100-105 K, below
which it sharply diverges at the critical temperature, following the
temperature dependence of the energy integrated experiment of figure 1(e).

###  III.3 Inelastic x-ray scattering, IXS

To disentangle static and dynamic contributions and to get further insights
about the origin of the observed diffuse features, we have carried out IXS
scanning the (2.5 0.5 0.5) satellite. As shown in figure 3(e), the IXS spectra
consists on a (quasi)-elastic central peak (zero energy loss) evolving in
intensity below 120 K and a low energy non-dispersive phonon at 6 meV. The
intensity of the CP sharply drops upon warming from 95 K and follows the
intensity of the pretranslational fluctuations detected by DS above 105 K. The
CP of CsV3Sb5 shows parallels to the critical scattering observed in SrTiO3
Ris71 ; Sha74  and transition metal dichalcogenides Mon77  exhibiting a
critical divergence on approaching TCDW. On the other hand, directly
visualized in the raw data of figure 3(e), the intensity of the low-energy
mode increases smoothly upon cooling and drops below 100 K due to the
modification of the dynamic structure factor on approaching the CDW phase,
making the quasi-elastic contribution more important than the phonon intensity
to the DS. Following the applicability of the Ornstein-Zernike correlation
function Rou86 , the generalized susceptibility, $\chi$(q), proportional to
the q Fourier component of the displacement-displacement correlation function
$\langle\lvert$uq$\lvert^{2}\rangle$=kBT$\chi(q)$, displays a linear Curie-
Weiss behavior for the quasi-elastic component of the CP (inset of figure
3(f)) and vanishes at the critical temperature, indicating a temperature
dependence of the density of states. Nevertheless, $\chi(q)$ of the DS changes
slope at T*$\sim$ 105 K, coinciding with the drop of intensity of the low-
energy phonon and the deviation of the DS and CP intensities Note1.

Focusing on the inelastic part of the IXS spectra, figure 4(a) shows the
momentum dependence of the phonon dispersion along the (h h 0.5) direction
($\Delta$E=3 meV). The spectra consist on 3 weakly dispersive branches with
2.5 ($\gamma$), 6 ($\omega$) and 11 meV ($\delta$) at h=0. The lowest energy
mode merges with the 6 meV branch at 0.1 r.l.u. and, at the L point, only two
optical-like modes, $\omega$ and $\delta$, are observed. These two modes are
better visualized in figure 4(b) and its energy dispersion, obtained from the
fitting to damped harmonic oscillators (DHO) convoluted with the experimental
resolution, is displayed in the figure 4(c). They exhibit a weak dispersion in
reciprocal space and can be described as a system of independent quantum
oscillators or localized atomic vibrations (Einstein-type phonon modes). To
isolate the contribution of each mode, we increased the energy resolution down
to 1.5 meV and work at the (3.5 2.5 0.5) reciprocal lattice vector. Figure
4(d) shows the temperature dependence of the $\omega$ branch. Indeed, the
$\omega$ phonon in figure 4(b) consists on 2 modes, which we label here as
$\omega_{1}$ and $\omega_{2}$, in good agreement the DFPT calculations in
figure 1(c) (Supplementary Figure 2). From these data, the energy of
$\omega_{1}$ and $\omega_{2}$ remains constant as the temperature approaches
the CDW transition. This clearly demonstrates the absence of a Kohn anomaly
and discards a soft phonon nature of the CDW phase transition; i.e, it is not
phonon driven but of order-disorder type. This does not preclude the presence
of a soft mode below TCDW, which is, in fact, a fingerprint of both order-
disorder and displacive phase transitions. We need to stress here that it is
not the crystal quality, but the inherent dynamical disorder that stabilizes
the CDW in CsV3Sb5. On the other hand, a detailed analysis of the phonon
linewidth reveals an anomalous broadening at TCDW. At 150 K, $\omega_{1}$ is
not resolution limited and develops a small intrinsic broadening of $\sim$0.5
meV (a DHO with 0 meV intrinsic width corresponds to the instrumental
resolution of 1.5 meV). Upon cooling, the width of the phonon decreases and
becomes resolution limited at 105 K. However, its width increases again (not
resolution limited) below 100 K, coinciding with T*, a signature of a
proximity to a Fermi surface instability.

Figure 4: Einstein modes and phonon broadening (a) Color map of the low-energy
phonons of CsV3Sb5 along the (3-h h 0.5) direction at 120 K. The three modes
identified at 2.5 ($\gamma$, asterisks), 6 ($\omega$, triangles) and 11 meV
($\delta$, squares) are superimposed. (b) IXS spectra and (c) dispersion along
the (2+h 1-h 0.5) direction ($\Delta$E= 3 meV) of the $\delta$ and $\omega$
branches. (d) Temperature dependent high resolution IXS spectra ($\Delta$E=
1.5 meV) at the 3.5 2.5 0.5 reciprocal lattice vector. Inset, fitting details
of the scan at 150 K. Note that, $\omega$ in (b) is now resolved as 2 modes,
$\omega_{1}$ and $\omega_{2}$. (e) Temperature dependence of the linewidth of
$\omega_{1}$. The error bars represent the fit uncertainty. See Supplementary
for the fitting details and Supplementary Figure 3.

##  IV Discussion

The results presented here indicate that weak coupling limit theories,
$\omega_{2k_{F}}\tau_{eh}<$1 ($\omega_{2k_{F}}$ and $\tau_{eh}$ are the bare
phonon frequency and electron-hole pair lifetime, respectively), are not
adequate to interpret the absence of the Kohn anomaly and the major
contribution of the CP to diffuse scattering data. This locates the CDW
instability of AV3Sb5 in the strong coupling non-adiabatic or in the strong
electron-phonon interaction (EPI) regimes. Indeed, recent reflectivity
measurements have reported a much larger CDW gap than the weak-coupling BCS
value Zhou21 ; Uyk21 . The effect of the dynamical disorder was theoretically
modelled by Yu and Anderson Yu83  for the strong interaction of electrons and
single dispersionless (Einstein-type) optical modes and further extended by
Gor’kov Gor12  for any arbitrary ionic displacements. Within this scenario,
the CDW in AV3Sb5 realizes in two stages. First, an ionic displacement in the
kagome plane introduces a strong EPI that traps the electronic cloud near an
ion, changing the local electronic environment and breaking the adiabatic
approximation. This gives rise to a quasi-elastic response above the CDW
transition, an associated order-disorder dynamics and the absence of a Kohn
anomaly. In an order-disorder transformation, the harmonic potential well
deforms and the phase transition is achieved by the critical growth and
percolation of quasi-static domains on approaching the critical temperature
Cow80 . In each domain, the potential within the atom moves is quasi-harmonic
and phonon anomalies are expected. This naturally explains the change of slope
in the diffuse susceptibility of figure 3(f) and the small broadening of the
low-energy mode in figure 4(d) and (e). Strong EPI mechanisms have been
considered to explain the CDW formation in the one dimensional Peierls
instability of BaVS3 Ila21 ; Gir19 , transition metal dichalcogenides Bar75 ,
H-bond ferroelectrics Bli06  and the controversial charge order in high-Tc
cuprates LeTacon14 . On the other hand, in the non-adiabatic dynamics,
$\omega_{2k_{F}}$ fluctuates so rapidly that cannot couple with the
electron–hole condensate during the electron–hole pair lifetime. The rapidly
fluctuating phonon reduces the screening between atoms and prevents the
development of the Kohn anomaly at TCDW. Following this argument, the
suppression of the Peierls transition of BaVS3 under pressure was interpreted
as a reduction of the $\tau_{eh}$ Ber12 . This is a particularly appealing
scenario since the TCDW of CsV3Sb5 is completely suppressed at low pressures
and evolves towards a double superconducting dome Yu21 ; Zhu22 ; Yu22 ; Wang21
.

We would like to finish our discussion by placing our results within the
physics of the kagome metals AV3Sb5. The identification of the charge order
fluctuations at T$>$TCDW indicates that the condensation of the CDW order
parameter at the L point in (Cs,Rb)V3Sb5 is the leading instability, favoring
the alternating Star of David and trihexagonal (LLL) pattern Fer21 . The
absence of DS at the M point suggests that the strength of the coupling
parameter of the L and M modes in the Landau free energy expansion,
$\gamma_{ML}$, is rather small or the phase transition is more strongly first-
order Note2. Our DS results also indicate that the 6-fold symmetry is
preserved, discarding single QM/double QL orders (staggered SoD and staggered
TrH). Rather, our data seem to agree with the experiments that support a
combination of triple QM/triple QL order (LLL+MMM pattern) ComNatMat ; Luo22a
at low temperature.

On the other hand, the importance of the EPI has been highlighted recently by
neutron scattering Xie22  and, perhaps, acts concomitantly with nesting to
enhance the charge correlations. Coherent phonon spectroscopy reported a
condensation of M and L-point phonons, suggesting an order-disorder crossover
at lower temperature Rat21  to a different ordering pattern as observed here.
An order-disorder transition is expected to be of first order type, thus the
small release of entropy at TCDW Ort20 ; Ort21  fits within this scenario.
Moreover, the data presented here does not conflict with the possible chiral
charge density wave, but put constrains on the nematic order Nie22 ; Li22 .
Broken rotational symmetries can be affected by local dynamic disorder,
fluctuating domains and defects that introduce pinning potentials. These
pinning centers and local distortions introduce variations of the amplitudes
of the local atomic displacements, which can be detected by pair distribution
function (PDF) analysis Upr22 .

In conclusion, we have presented diffuse scattering and inelastic x-ray
scattering data to demonstrate that the CDW in the kagome metal AV3Sb5 is of
order-disorder type without phonon softening. A precursor of the CDW is
observed in diffuse scattering with propagation vector (0.5 0 0.5), indicating
that the CDW instability at the L point is the primary order parameter,
discarding a broken six-fold rotational symmetry at high temperature and
pointing to an LLL pattern. Furthermore, we demonstrate that the diffuse
signal contains intensity of the quasi-elastic central peak, likely due to the
growth of fluctuating domains, thus favouring strong coupling theories to
describe the CDW ground state of AV3Sb5. Our results present a step forward
and provide crucial information to understand the enigmatic charge order of
(Cs,Rb)V3Sb5.

##  V Methods

Single crystals of AV3Sb5 (A=Cs, Rb) with TCDW=95 and 104 K, were synthesized
by the self-flux growth method Ort19 . Diffuse scattering measurements were
performed with a Dectris PILATUS3 1M X area detector at the ID28 beamline at
the European Synchrotron Research Facility (ESRF) with an incoming energy
Ei=17.8 keV. Orientation matrix refinement and reciprocal space
reconstructions were done using the CrysAlis software package. Final DS images
were processed with an in-house code. Low energy phonons were measured by
inelastic x-ray scattering (IXS) at the ID28 IXS station at ESRF (Ei=17.8 keV,
$\Delta$E=3 meV) and at the HERIX beamline at the Argonne Photon Source (APS)
(Ei=23.72 keV, $\Delta$E=1.5 meV). The components ($h$ $k$ $l$) of the
scattering vector are expressed in reciprocal lattice units (r.l.u.), ($h$ $k$
$l$)= $h\mathbf{a}^{*}+k\mathbf{b}^{*}+l\mathbf{c}^{*}$, where
$\mathbf{a}^{*}$, $\mathbf{b}^{*}$, and $\mathbf{c}^{*}$ are the reciprocal
lattice vectors.

Lattice dynamic calculations were performed by Density Functional Perturbation
Theory (DFPT) as implemented in QUANTUM ESPRESSO package. The forces were
calculated in 8$\times$8$\times$4 supercells using the ultrasoft
parametrization of the PBE exchange-correlation functional. We set the energy
cutoffs of 60 and 600 Ry for the basis-set and treated the Van der Waals
corrections within the Grimme’s semiempirical approach. An 8$\times$8$\times$4
k-point grid and a Marzari-Vanderbilt smearing of 0.01 Ry was used for the
Brillouin zone integration. The convergence threshold was set to 10-19.

##  VI Data availability

Source data are provided with this paper. The DS and IXS data generated in
this study have been deposited in the Figshare database
https://doi.org/10.6084/m9.figshare.21990461.

## References

  * (1)

##  VII References

* (2) Ye, L., Kang, M., Liu, J., Cube, F., Wicker, C., Suzuki, T., Jozwiak, C., Bostwick, A., Rotenberg, E., Bell, D., Fu, L., Comin, R. & Checkelsky, J. Massive Dirac fermions in a ferromagnetic kagome metal. Nature. 555, 638-642 (2018).
* (3) Kang, M., Ye, L., Fang, S., You, J., Levitan, A., Han, M., Facio, J., Jozwiak, C., Bostwick, A., Rotenberg, E., Chan, M., McDonald, R., Graf, D., Kaznatcheev, K., Vescovo, E., Bell, D., Kaxiras, E., Brink, J., Richter, M., Prasad Ghimire, M., Checkelsky, J. & Comin, R. Dirac fermions and flat bands in the ideal kagome metal FeSn. Nature Materials. 19, 163-169 (2020).
* (4) Liu, E., Sun, Y., Kumar, N., Muechler, L., Sun, A., Jiao, L., Yang, S., Liu, D., Liang, A., Xu, Q., Kroder, J., Süß, V., Borrmann, H., Shekhar, C., Wang, Z., Xi, C., Wang, W., Schnelle, W., Wirth, S., Chen, Y., Goennenwein, S. & Felser, C. Giant anomalous Hall effect in a ferromagnetic kagome-lattice semimetal. Nature Physics. 14, 1125-1131 (2018).
* (5) Yan, S., Huse, D. & Steven R. White Spin-Liquid Ground State of the S = 1/2 Kagome Heisenberg Antiferromagnet. Science. 332, 1173-1176 (2011).
* (6) Han, T., Helton, J., Chu, S., Nocera, D., Rodriguez-Rivera, J., Broholm, C. & Lee, Y. Fractionalized excitations in the spin-liquid state of a kagome-lattice antiferromagnet. Nature. 492, 406-410 (2012).
* (7) Mazin, I., Jeschke, H., Lechermann, F., Lee, H., Fink, M., Thomale, R. & Valentí, R. Theoretical prediction of a strongly correlated Dirac metal. Nature Communications. 5, 4261 (2014).
* (8) Ortiz, B., Gomes, L., Morey, J., Winiarski, M., Bordelon, M., Mangum, J., Oswald, I., Rodriguez-Rivera, J., Neilson, J., Wilson, S., Ertekin, E., McQueen, T. & Toberer, E. New kagome prototype materials: discovery of KV3Sb5, RbV3Sb5, and CsV3Sb5. Phys. Rev. Materials. 3, 094407 (2019).
* (9) Ortiz, B., Teicher, S., Hu, Y., Zuo, J., Sarte, P., Schueller, E., Abeykoon, A., Krogstad, M., Rosenkranz, S., Osborn, R., Seshadri, R., Balents, L., He, J. & Wilson, S. CsV3Sb5: A $\mathbb{Z}_{2}$ Topological Kagome Metal with a Superconducting Ground State. Phys. Rev. Lett.. 125, 247002 (2020).
* (10) Ortiz, B., Sarte, P., Kenney, E., Graf, M., Teicher, S., Seshadri, R. & Wilson, S. Superconductivity in the $\mathbb{Z}_{2}$ kagome metal KV3Sb5. Phys. Rev. Materials. 5, 034801 (2021).
* (11) Kang, M., Fang, S., Kim, J., Ortiz, B., Ryu, S., Kim, J., Yoo, J., Sangiovanni, G., Di Sante, D., Park, B., Jozwiak, C., Bostwick, A., Rotenberg, E., Kaxiras, E., Wilson, S., Park, J. & Comin, R. Twofold van Hove singularity and origin of charge order in topological kagome superconductor CsV3Sb5. Nature Physics. 18, 301-308 (2022).
* (12) Lou, R., Fedorov, A., Yin, Q., Kuibarov, A., Tu, Z., Gong, C., Schwier, E., Büchner, B., Lei, H. & Borisenko, S. Charge-Density-Wave-Induced Peak-Dip-Hump Structure and the Multiband Superconductivity in a Kagome Superconductor CsV3Sb5. Phys. Rev. Lett.. 128, 036402 (2022).
* (13) Zhao, H., Li, H., Ortiz, B., Teicher, S., Park, T., Ye, M., Wang, Z., Balents, L., Wilson, S. & Zeljkovic, I. Cascade of correlated electron states in the kagome superconductor CsV3Sb5. Nature. 599, 216-221 (2021).
* (14) Zhao, C., Wang, L., Xia, W., Yin, Q., Ni, J., Huang, Y., Tu, C., Tao, Z., Tu, Z., Gong, C., Lei, H., Guo, Y., Yang, X. & Li, S. Nodal superconductivity and superconducting domes in the topological Kagome metal CsV3Sb5. ArXiv:2102.08356v2.
* (15) Wu, X., Schwemmer, T., Müller, T., Consiglio, A., Sangiovanni, G., Di Sante, D., Iqbal, Y., Hanke, W., Schnyder, A., Denner, M., Fischer, M., Neupert, T. & Thomale, R. Nature of Unconventional Pairing in the Kagome Superconductors AV3Sb5 (A=K, Rb, Cs). Phys. Rev. Lett.. 127, 177001 (2021).
* (16) Tan, H., Liu, Y., Wang, Z. & Yan, B. Charge Density Waves and Electronic Properties of Superconducting Kagome Metals. Phys. Rev. Lett.. 127, 046401 (2021).
* (17) Subedi, Alaska. Hexagonal-to-base-centered-orthorhombic $4Q$ charge density wave order in kagome metals KV3Sb5,RbV3Sb5, and CsV3Sb5, Phys. Rev. Materials. 16, 015001 (2022).
* (18) Kiesel, M., Platt, C. & Thomale, R. Unconventional Fermi Surface Instabilities in the Kagome Hubbard Model. Phys. Rev. Lett.. 110, 126405 (2013).
* (19) Wang, W., Li, Z., Xiang, Y. & Wang, Q. Competing electronic orders on kagome lattices at van Hove filling. Phys. Rev. B. 87, 115135 (2013).
* (20) Jiang, Y., Yin, J., Denner, M., Shumiya, N., Ortiz, B., Xu, G., Guguchia, Z., He, J., Hossain, M., Liu, X., Ruff, J., Kautzsch, L., Zhang, S., Chang, G., Belopolski, I., Zhang, Q., Cochran, T., Multer, D., Litskevich, M., Cheng, Z., Yang, X., Wang, Z., Thomale, R., Neupert, T., Wilson, S. & Hasan, M. Unconventional chiral charge order in kagome superconductor KV3Sb5. Nature Materials. 20, 1353-1357 (2021).
* (21) Miao, H., Li, H., Meier, W., Huon, A., Lee, H., Said, A., Lei, H., Ortiz, B., Wilson, S., Yin, J., Hasan, M., Wang, Z., Tan, H. & Yan, B. Geometry of the charge density wave in the kagome metal AV3Sb5. Phys. Rev. B. 104, 195132 (2021).
* (22) Mielke, C., Das, D., Yin, J., Liu, H., Gupta, R., Jiang, Y., Medarde, M., Wu, X., Lei, H., Chang, J., Dai, P., Si, Q., Miao, H., Thomale, R., Neupert, T., Shi, Y., Khasanov, R., Hasan, M., Luetkens, H. & Guguchia, Z. Time-reversal symmetry-breaking charge order in a kagome superconductor. Nature. 602, 245-250 (2022).
* (23) Li, H., Zhao, H., Ortiz, B., Park, T., Ye, M., Balents, L., Wang, Z., Wilson, S. & Zeljkovic, I. Rotation symmetry breaking in the normal state of a kagome superconductor KV3Sb5. Nature Physics. 18, 265-270 (2022).
* (24) Nie, L., Sun, K., Ma, W., Song, D., Zheng, L., Liang, Z., Wu, P., Yu, F., Li, J., Shan, M., Zhao, D., Li, S., Kang, B., Wu, Z., Zhou, Y., Liu, K., Xiang, Z., Ying, J., Wang, Z., Wu, T. & Chen, X. Charge-density-wave-driven electronic nematicity in a kagome superconductor. Nature. 604, 59-64 (2022).
* (25) Chen, H., Yang, H., Hu, B., Zhao, Z., Yuan, J., Xing, Y., Qian, G., Huang, Z., Li, G., Ye, Y., Ma, S., Ni, S., Zhang, H., Yin, Q., Gong, C., Tu, Z., Lei, H., Tan, H., Zhou, S., Shen, C., Dong, X., Yan, B., Wang, Z. & Gao, H. Roton pair density wave in a strong-coupling kagome superconductor. Nature. 599, 222-228 (2021).
* (26) Denner, M., Thomale, R. & Neupert, T. Analysis of Charge Order in the Kagome Metal AV3Sb5 (A=K, Rb, Cs). Phys. Rev. Lett.. 127, 217601 (2021).
* (27) Feng, X., Jiang, K., Wang, Z. & Hu, J. Chiral flux phase in the Kagome superconductor AV3Sb5. Science Bulletin. 66, 1384-1388 (2021).
* (28) Li, H., Zhang, T., Yilmaz, T., Pai, Y., Marvinney, C., Said, A., Yin, Q., Gong, C., Tu, Z., Vescovo, E., Nelson, C., Moore, R., Murakami, S., Lei, H., Lee, H., Lawrie, B. & Miao, H. Observation of Unconventional Charge Density Wave without Acoustic Phonon Anomaly in Kagome Superconductors AV3Sb5 (A=Rb, Cs). Phys. Rev. X. 11, 031050 (20219).
* (29) Luo, H., Gao, Q., Liu, H., Gu, Y., Wu, D., Yi, C., Jia, J., Wu, S., Luo, X., Xu, Y., Zhao, L., Wang, Q., Mao, H., Liu, G., Zhu, Z., Shi, Y., Jiang, K., Hu, J., Xu, Z. & Zhou, X. Electronic nature of charge density wave and electron-phonon coupling in kagome superconductor KV3Sb5. Nature Communications. 13, 273 (2022).
* (30) Zhou, X., Li, Y., Fan, X., Hao, J., Dai, Y., Wang, Z., Yao, Y. & Wen, H. Origin of charge density wave in the kagome metal $\mathrm{CsV_{3}Sb_{5}}$ as revealed by optical spectroscopy. Phys. Rev. B. 104, L041101 (2021).
* (31) Uykur, E., Ortiz, B., Wilson, S., Dressel, M. & Tsirlin, A. Optical detection of the density-wave instability in the kagome metal KV3Sb5. Npj Quantum Materials. 7, 16 (2022).
* (32) Nakayama, K., Li, Y., Kato, T., Liu, M., Wang, Z., Takahashi, T., Yao, Y. & Sato, T. Multiple energy scales and anisotropic energy gap in the charge-density-wave phase of the kagome superconductor CsV3Sb5. Phys. Rev. B. 104, L161112 (2021).
* (33) Riste, T., Samuelsen, E., Otnes, K. & Feder, J. Critical behaviour of SrTiO3 near the 105°K phase transition. Solid State Communications. 9, 1455-1458 (1971).
* (34) Shapiro, S., Axe, J., Shirane, G. & Riste, T. Critical Neutron Scattering in SrTiO3 and KMnF3. Phys. Rev. B. 6, 4332-4341 (1972).
* (35) Moncton, D., Axe, J. & DiSalvo, F. Neutron scattering study of the charge-density wave transitions in 2H-TaSe2 and 2H-Nb$\mathrm{Se_{2}}$. Phys. Rev. B. 16, 801-819 (1977).
* (36) Rouxel, J. Crystal Chemistry and Properties of Materials with Quasi-One-Dimensional Structures. (D. Riedel, Dordrecht, 1986).
* (37) We tentatively associate T* as the temperature where the CDW domains start to percolate.
* (38) Uykur, E., Ortiz, B., Iakutkina, O., Wenzel, M., Wilson, S., Dressel, M. & Tsirlin, A. Low-energy optical properties of the nonmagnetic kagome metal CsV3Sb5. Phys. Rev. B. 104, 045130 (2021).
* (39) Yu, C. & Anderson, P. Local-phonon model of strong electron-phonon interactions in A15 compounds and other strong-coupling superconductors. Phys. Rev. B. 29, 6165-6186 (1984).
* (40) Gor’kov, L. Strong electron-lattice coupling as the mechanism behind charge density wave transformations in transition-metal dichalcogenides. Phys. Rev. B. 85, 165142 (2012).
* (41) Bruce, A. & R.A. Cowley Structural phase transitions III. Critical dynamics and quasi-elastic scattering. Advances In Physics. 29, 219-321 (1980).
* (42) Ilakovac, V., Girard, A., Balédent, V., Foury-Leylekian, P., Winkler, B., Kup, I., Berger, H., Bosak, A. & Pouget, J. Order-disorder type of Peierls instability in BaVS3. Phys. Rev. B. 103, 014306 (2021).
* (43) Girard, A., Ilakovac, V., Stekiel, M., Morgenroth, W., Berger, H., Yoshikazu, T., Hasegawa, T., Bosak, A. & Winkler, B. Soft-phonon driven hexagonal-orthorhombic phase transition in BaVS3. Phys. Rev. B. 99, 144104 (2019).
* (44) Barker, A., Ditzenberger, J. & DiSalvo, F. Infrared study of the electronic instabilities in tantalum disulfide and tantalum diselenide. Phys. Rev. B. 12, 2049-2054 (1975).
* (45) Blinc, R. & B. Žekš Dynamics of order-disorder-type ferroelectrics and anti-ferroelectrics. Advances In Physics. 21, 693-757 (1972).
* (46) Le Tacon, M., Bosak, A., Souliou, S., Dellea, G., Loew, T., Heid, R., Bohnen, K., Ghiringhelli, G., Krisch, M. & Keimer, B. Inelastic X-ray scattering in YBa2Cu3O6.6 reveals giant phonon anomalies and elastic central peak due to charge-density-wave formation. Nature Physics. 10, 52-58 (2014).
* (47) Bernu, S., Fertey, P., Itié, J., Berger, H., Foury-Leylekian, P. & Pouget, J. Vanishing of the metal-insulator Peierls transition in pressurized BaVS3. Phys. Rev. B. 86, 235105 (2012).
* (48) Yu, F., Ma, D., Zhuo, W., Liu, S., Wen, X., Lei, B., Ying, J. & Chen, X. Unusual competition of superconductivity and charge-density-wave state in a compressed topological kagome metal. Nature Communications. 12, 3645 (2021).
* (49) Zhu, C., Yang, X., Xia, W., Yin, Q., Wang, L., Zhao, C., Dai, D., Tu, C., Song, B., Tao, Z., Tu, Z., Gong, C., Lei, H., Guo, Y. & Li, S. Double-dome superconductivity under pressure in the V-based kagome metals AV3Sb5 (A=Rb and K). Phys. Rev. B. 105, 094507 (2022).
* (50) Yu, F., Zhu, X., Wen, X., Gui, Z., Li, Z., Han, Y., Wu, T., Wang, Z., Xiang, Z., Qiao, Z., Ying, J. & Chen, X. Pressure-Induced Dimensional Crossover in a Kagome Superconductor. Phys. Rev. Lett.. 128, 077001 (2022).
* (51) Wang, N., Chen, K., Yin, Q., Ma, Y., Pan, B., Yang, X., Ji, X., Wu, S., Shan, P., Xu, S., Tu, Z., Gong, C., Liu, G., Li, G., Uwatoko, Y., Dong, X., Lei, H., Sun, J. & Cheng, J. Competition between charge-density-wave and superconductivity in the kagome metal RbV3Sb5. Phys. Rev. Research. 3, 043018 (2021).
* (52) Morten H. Christensen, Turan Birol, Brian M. Andersen and Rafael M. Fernandes, Theory of the charge density wave in $A{\mathrm{V}}_{3}{\mathrm{Sb}}_{5}$ kagome metals Phys. Rev. B. 104, 214513, (2021).
* (53) The condition for the susceptibility at the M point to diverge is $\xi^{2}$ $\sim$ $\gamma$ML, where $\xi$ is the correlation length of the L order. If $\xi$ never becomes too large, one would not expect an enhancement of the M susceptibility above TCDW. R. M. Fernandes, Private Communication.
* (54) Mingu Kang, Shiang Fang, Jonggyu Yoo, Brenden R. Ortiz, Yuzki M. Oey, Jonghyeok Choi, Sae Hee Ryu, Jimin Kim, Chris Jozwiak, Aaron Bostwick, Eli Rotenberg, Efthimios Kaxiras, Joseph G. Checkelsky, Stephen D. Wilson, Jae-Hoon Park and Riccardo Comin. Charge order landscape and competition with superconductivity in kagome metals, Nature Materials, (2022). https://doi.org/10.1038/s41563-022-01375-2.
* (55) Luo, H., Gao, Q., Liu, H., Gu, Y., Wu, D., Yi, C., Jia, J., Wu, S., Luo, X., Xu, Y., Zhao, L., Wang, Q., Mao, H., Liu, G., Zhu, Z., Shi, Y., Jiang, K., Hu, J., Xu, Z. & Zhou, X. Electronic nature of charge density wave and electron-phonon coupling in kagome superconductor KV3Sb5. npj Quantum Mater 7, 30 (2022).
* (56) Yaofeng Xie, Yongkai Li, Philippe Bourges, Alexandre Ivanov, Zijin Ye, Jia-Xin Yin, M. Zahid Hasan, Aiyun Luo, Yugui Yao, Zhiwei Wang, Gang Xu, and Pengcheng Dai. Electron-phonon coupling in the charge density wave state of CsV3Sb5. Phys. Rev. B. 105, L140501 (2022).
* (57) Noah Ratcliff, Lily Hallett, Brenden R. Ortiz, Stephen D. Wilson, and John W. Harter. Coherent phonon spectroscopy and interlayer modulation of charge density wave order in the kagome metal CsV3Sb5. Phys. Rev. Materials. 5, L111801 (2021).
* (58) Upreti, P., Krogstad, M., Haley, C., Anitescu, M., Rao, V., Poudel, L., Chmaissem, O., Rosenkranz, S. & Osborn, R. Order-Disorder Transitions in (CaxSr1-x)3Rh4Sn13. Phys. Rev. Lett.. 128, 095701 (2022).

##  VIII Acknowledgments

We acknowledge E. Efremov, A. Schnyder, A. Bernevig, A. Subedi, J. Diego, L.
Elcoro, I. Etxebarria, D. Chernyshov and R. Fernandes for fruitful discussions
and critical reading of the manuscript. S.B-C, D.S and L.S thank the MINECO of
Spain through the project PGC2018-101334-A-C22. S.D.W and B.R.O gratefully
acknowledge support via UC Santa Barbara NSF Quantum Foundry funded via the
Q-AMASE-i program under award DMR-1906325. This research used resources of the
Advanced Photon Source, a U.S. Department of Energy (DOE) Office of Science
user facility operated for the DOE Office of Science by Argonne National
Laboratory under Contract No. DE-AC02-06CH11357.

##  IX Author Contribution

S.B-C conceived and managed the project. B.O and S.D.W. synthesized and
characterized the samples. A.K., A.B. and S.B-C performed the DS experiments.
D.S., A.K., A.H.S., L.S., A.B. and S.B-C carried out the IXS measurements.
S.B-C performed the DFPT calculations, analyzed the data and wrote the
manuscript with input from all coauthors.

##  X Competing Interests

The authors declare no competing interests.