††thanks: These authors contributed equally to this work.††thanks: These
authors contributed equally to this work.††thanks: These authors contributed
equally to this work.

# Revealing kinetically tuned atomic pathways for interfacial strain
relaxation

Sophia B. Betzler  Materials Sciences Division, Lawrence Berkeley National
Laboratory; Berkeley, CA, USA Daewon Lee  Materials Sciences Division,
Lawrence Berkeley National Laboratory; Berkeley, CA, USA Department of
Materials Science and Engineering, University of California, Berkeley;
Berkeley, CA, USA Sam Oaks-Leaf  Department of Chemistry, University of
California, Berkeley; Berkeley, CA, USA Colin Ophus  National Center for
Electron Microscopy, Molecular Foundry, Lawrence Berkeley National Laboratory;
Berkeley, CA, USA  Lin-Wang Wang  State Key Laboratory of Superlattices and
Microstructures, Institute of Semiconductors, Chinese Academy of Science;
Beijing, China  Mark Asta  Materials Sciences Division, Lawrence Berkeley
National Laboratory; Berkeley, CA, USA Department of Materials Science and
Engineering, University of California, Berkeley; Berkeley, CA, USA David T.
Limmer*<EMAIL_ADDRESS>Materials Sciences Division, Lawrence Berkeley
National Laboratory; Berkeley, CA, USA Department of Chemistry, University of
California, Berkeley; Berkeley, CA, USA Chemical Sciences Division, Lawrence
Berkeley National Laboratory; Berkeley, CA, USA Kavli Energy Nanoscience
Institute; Berkeley, CA, USA Haimei Zheng*<EMAIL_ADDRESS>Materials Sciences
Division, Lawrence Berkeley National Laboratory; Berkeley, CA, USA Department
of Materials Science and Engineering, University of California, Berkeley;
Berkeley, CA, USA

###### Abstract

Strain at interfaces may profoundly impact the microstructure and properties
of materials; thus, it is a major consideration when designing and engineering
materials[1, 2, 3, 4, 5, 6, 7, 8]. Dislocation formation is a commonly known
mechanism to release mismatch strain at solid-solid interfaces. However, it is
still unclear about how materials accommodate interfacial strain under
drastically accelerated structural transformation kinetics, since it is
extremely challenging to directly observe the atomic structure evolution of
fast-propagating interfaces. Utilizing liquid phase transmission electron
microscopy (TEM), we have achieved atomic-scale imaging of hydrogen-induced
phase transformations of palladium nanocrystals with different transformation
speeds. Our observation reveals that the fast phase transformation occurs with
an expanded interface of mixed $\alpha$\- and $\beta$-$\mathrm{PdH}_{x}$
phases, and tilting of (020) planes to accommodate mismatch strain. In
contrast, slow phase transformations lead to sharp interfaces with slipping
misfit dislocations. Our kinetic Monte Carlo simulations show that fast phase
transformation pushes the system far-from-equilibrium, generically roughening
the interface; however, a smooth boundary minimizes strain near-equilibrium.
Unveiling the atomic pathways of transformations from near-equilibrium to far-
from-equilibrium, which was previously possible only computationally, this
work holds significant implications for engineering microstructure of
materials through modulating solid-solid transformations in a wide range of
kinetics.

## Introduction

Understanding structural modifications at interfaces by mismatch strain during
growth or structural transformations of materials has been a long-standing
significant topic of study[9, 10, 11, 12]. Compared to the spectroscopic or
computational approaches, in situ transmission electron microscopy (TEM)
allows us to directly observe spatially resolved atomic structure evolution
coupled with the strain relaxation at interfaces during dynamic
transformations. For example, the recent work[2] on imaging of dislocation
movements along $\mathrm{Cu}_{2}\mathrm{O}$/Cu interfaces during the reduction
reaction of $\mathrm{Cu}_{2}\mathrm{O}$ to Cu in $\mathrm{H}_{2}$ marked the
technical advances in observing dynamic dislocation motion under complex
environmental stimuli beyond deformation[13]. However, due to technical
challenges in varying the kinetics of structural transformations over a wide
range while observing the respective atomic structure evolution, it has not
been possible to extend our understanding of the mismatch strain relaxation
mechanisms beyond dislocation formation. Specifically, we have not been able
to uncover alternative actual routes to accommodate mismatch strain at solid-
solid interfaces in the nonequilibrium regime that are relevant to the
practical engineering of materials.

Here, utilizing high-resolution liquid phase TEM, we track the atomic
structural evolution of hydrogen-induced phase transformations in palladium
(Pd) nanocrystals, which has enabled breakthroughs in revealing mismatch
strain relaxation mechanisms at interfaces tuned by phase transformation
kinetics. The electron beam radiolysis of water generates hydrogen gases
inside a liquid cell. By controlling the amount of hydrogen in the cell, we
can modulate the phase transformation kinetics. Two immiscible phases can be
obtained when hydrogen is absorbed into Pd, which are
$\alpha$-$\mathrm{PdH}_{x}$ $(x<0.03)$ and $\beta$-$\mathrm{PdH}_{x}$
$(x>0.6)$, with a large lattice mismatch by 3.5% between them[14, 15]. Thus,
hydrogen absorption and desorption in a palladium nanocrystal generate lattice
mismatch strain at the interfaces between $\alpha$\- and
$\beta$-$\mathrm{PdH}_{x}$. Hydrogen desorption is a reverse process of
hydrogen absorption, providing an additional degree of freedom in tuning the
phase transformation kinetics. Despite significant prior efforts in
investigating hydrogen sorption behaviors of Pd nanoparticles[4, 6, 16, 17,
18, 19, 20, 21, 22, 23, 24, 25, 26, 27] in situ atomic-scale imaging of
hydrogen-induced phase transformations in Pd nanocrystals has not been
achieved. Our unprecedented capability allows us to capture the atomic
structure evolution at interfaces and reveal distinct mismatch strain
relaxation mechanisms under drastically different phase transformation
kinetics.

## Disparate phase transformation kinetics

Fig. 1a illustrates the experimental setup for liquid phase TEM study of
hydrogen-induced phase transformations in Pd nanocrystals. The Pd nanocubes
(40-50 nm) are synthesized following the published synthesis protocol[25] (see
Methods). Representative high-resolution TEM (HRTEM), bright-field TEM, and
scanning TEM images of a Pd nanocube are displayed in Fig. 1b and
Supplementary Figs. 1 and 2. A droplet of Pd nanocubes in 0.1 M KOH aqueous
solution is loaded into a liquid cell, and small liquid pockets are formed
inside the liquid cell[28, 29] (see Methods and Supplementary Fig. 3).
Radiolysis of water yields many products, including the predominant hydrogen
and oxygen gases (molar ratio of hydrogen:oxygen estimated to be 60:40)[30,
31, 32]. By controlling the dose rate used during TEM experiments, we can
control the amount of hydrogen in the liquid cell and thus the transformation
kinetics including the phase propagation pathways and propagation rate (see
Methods and Supplementary Figs. 4 and 5).

Figure 1: In situ liquid phase TEM study of hydrogen absorption and desorption
in Pd nanocubes. a, Schematic illustration of the experimental setup. The
electron beam radiolysis of water generates hydrogen gases in the liquid cell,
and the hydrogen sorption behaviors in Pd nanocubes are observed in real-time
with atomic resolution. b, An HRTEM image of Pd nanocube showing single-
crystal structure. c,d, Sequential HRTEM images and corresponding $d$-spacing
color maps of Pd nanocubes during hydrogen absorption (c) and desorption (d),
leading to the $\beta$-$\mathrm{PdH}_{x}$ and $\alpha$-$\mathrm{PdH}_{x}$
phases, respectively. e, The average propagation speed of the interface along
the [100] direction as a function of time during hydrogen absorption (left
$y$-axis) and desorption (right $y$-axis). The propagation speed of hydrogen
absorption is about one order magnitude faster than the desorption. Details
about the interface position measurements and methods to derive the average
propagation speed are described in Supplementary Figs. 10 and 11.

During hydrogen absorption in a Pd nanocube, the $\beta$-$\mathrm{PdH}_{x}$
phase nucleates at a corner, and it subsequently propagates along the [100]
direction (see Methods, Supplementary Note 1, Supplementary Figs. 6 and 7,
Supplementary Video 1, and Extended Data Fig. 1). This is consistent with the
previously reported phase transformation pathway[22, 23, 24, 25]. To
understand the basic thermodynamic driving forces underpinning the observed
transformation pathway, we establish a minimal model to simulate the nanoscale
dynamics of solid-solid phase transformations. We first calculate both the
energetic and the free energetic cost of forming interfacial structures within
our model (see Supplementary Notes 2 and 3, Supplementary Figs. 8 and 9, and
Supplementary Tables 1 and 2). Our simulations using an elastic Ising model
Hamiltonian[33, 34, 35] indicate that the pathway of corner nucleation
followed by propagation along the [100] direction during hydrogen absorption
can be understood predominantly in terms of the energetics of elastic
interactions (see Extended Data Figs. 2 and 3, Supplementary Note 4, and
Supplementary Video 2). The reversible transformations of hydrogen absorption
and desorption in Pd nanocrystals can be understood by the same set of elastic
interaction terms.

To understand the effects of phase propagation rates on the interface
structures, we focus on the simple configurations of phase transformations
with the interfaces propagating along the [100] directions. We measure the
rates of phase propagation along the [100] directions during both hydrogen
absorption and desorption in Pd nanocrystals. Dramatically different phase
propagation rates are found in the hydrogen absorption and desorption.
Representative sequential HRTEM images and the corresponding $d$-spacing phase
maps show the phase propagation in both scenarios (Fig. 1c,d and Supplementary
Videos 3 and 4). The phase propagation rates, estimated by the interface
position relative to the starting position as a function of time, are
illustrated in Fig. 1e and Supplementary Figs. 10 and 11.

During hydrogen absorption (Fig. 1e, left y-axis and Supplementary Figs. 10a
and 11a), the interface propagates along the [100] direction with an average
speed of 1.3 nm/s at the early stage. When roughly 50% of the nanocube is
transformed, the average propagation speed increases to 3.0 nm/s. The
interphase boundary propagates rapidly in the later stages with the average
propagation speed of 5.9 nm/s. And, the movement of the interface is not
perfectly parallel to the edge of the nanocrystal, with the propagating
interface on the left side slower than the right side, which might be
attributed to the variations of strain along different crystallographic
directions[23]. During the desorption (Fig. 1e, right $y$-axis and
Supplementary Figs. 10b and 11b), a significantly slower propagation speed of
0.2 nm/s is obtained. Similar to the hydrogen absorption, phase propagation
occurs rapidly at the last stages of the desorption, with the measured average
speed of 0.8 nm/s. Our kinetic model predicts that both absorption and
desorption are fastest at the corners of the nanocubes, which could lead to an
increase in propagation speed at the end as the interface nears the cube
surface. We confirm the completion of phase transformations into
$\alpha$-$\mathrm{PdH}_{x}$ after the desorption or $\beta$-$\mathrm{PdH}_{x}$
after the absorption by measuring lattice $d$-spacings (Supplementary Figs. 12
and 13). The final nanocubes are perfect single crystals without noticeable
defects. Moreover, we have observed the swelling or shrinking of the nanocubes
due to the uptake or release of hydrogen (Supplementary Figs. 14 and 15). The
differences of one order of magnitude in phase propagation rates may arise
from the additional kinetic barrier to form hydrogen bonds when hydrogen comes
out from Pd nanocrystals during hydrogen desorption[26].

## A unique expanded interface in fast kinetics

Figure 2: Analysis of an expanded interface between $\alpha$\- and
$\beta$-$\mathrm{PdH}_{x}$ during hydrogen absorption with fast phase
propagation kinetics. a, An overview image of the Pd nanocube with the visible
interface. The respective FFT pattern (inset) demonstrates the co-existing two
crystal phases by the splitting reflections. b, The (200) $d$-spacing map of
the Pd nanocube. c, An enlarged view of the selected region in (a)
demonstrating the atomic structure at the interface. The last unaltered atomic
columns at the interphase boundary are highlighted in red in both $\alpha$\-
and $\beta$-$\mathrm{PdH}_{x}$ regions. d, GPA performed on the [100]
($\mathrm{E}_{yy}$) direction based on the high-resolution image displayed in
(c). Compared with the smooth contrast in the regions of
$\alpha$-$\mathrm{PdH}_{x}$ (selected as GPA reference area) and
$\beta$-$\mathrm{PdH}_{x}$, drastic contrast changes within the expanded
interface suggest significant $d$-spacing variations. e,f, The atomic
arrangement analysis corresponding to the enlarged region in (c). (e) The
$d$-spacing variations of (200) planes. We can identify mixed $\beta$\- and
$\alpha$-$\mathrm{PdH}_{x}$ phases within the expanded interface. (f) The
tilting angle of (020) planes (obtained by linear fitting of four atomic
columns).

To examine the atomic structure evolution of the interfaces during hydrogen
absorption, we analyze a selected representative HRTEM image (from
Supplementary Video 3) to elucidate the atomic arrangements (Fig. 2a). The
presence of two crystal phases is reflected in the splitting of {200} and
{220} diffraction spots in the fast Fourier transformation (FFT) pattern (Fig.
2a inset and Supplementary Fig. 16), and it can also be visualized by tracking
the $d$-spacing within the nanocube (Fig. 2b and Supplementary Fig. 16). A
spatially expanded interface is found as illustrated in Fig. 2c, where the
last atomic columns without obvious displacements in both
$\alpha$-$\mathrm{PdH}_{x}$ and the newly formed $\beta$-$\mathrm{PdH}_{x}$
are highlighted. We employ geometric phase analysis (GPA)[36, 37] to analyze
the lattice mismatch strain. We analyze the strain in both [010]
($\mathrm{E}_{xx}$) and [100] ($\mathrm{E}_{yy}$) directions. An overlay of
the $\mathrm{E}_{yy}$ strain map and the corresponding high-resolution image
is displayed in Fig. 2d. It is well known that it is challenging to obtain the
quantitative GPA results at interfaces[38], which is reflected by the
significantly wide range of strain values in Fig. 2d. More careful
investigation of the $d$-spacing values in the real space is carried out in
Fig. 2e. The $\mathrm{E}_{yy}$ GPA strain map qualitatively demonstrates large
variations of $d$-spacing at the expanded interface. Relatively smooth strain
variations are observed at the interface in the $\mathrm{E}_{xx}$ direction
(Supplementary Fig. 17). We examine the highlighted area shown in Fig. 2c, and
no edge dislocations can be identified. Dislocations are found within the
expanded interface outside the field of view. These illustrate more complex
interphase boundary structures during hydrogen absorption with the rapid phase
propagation kinetics.

We further calculate the $d$-spacing of the (200) planes as well as the
tilting of (020) lattice planes. The color map of the $d$-spacing values
corresponding to the image in Fig. 2c is displayed in Fig. 2e. The region with
the missing data within the interphase boundary is due to the atomic blurring,
likely from the fast phase propagation. Interestingly, we can identify mixed
$\beta$-$\mathrm{PdH}_{x}$ and $\alpha$-$\mathrm{PdH}_{x}$ phases within the
expanded interface and the strain variations. Fig. 2f illustrates the local
tilt angle of the (020) planes, obtained by fitting a line to the position of
four atoms in the [100] direction. We find the lattice tilts around +5°
relative to the vertical axis close to the upper phase propagation front and
around -5° near the lower part of the interface. Remarkably, the observed
lattice imperfection within the expanded interface (lattice tilting and
mixture of $\alpha$\- and $\beta$-$\mathrm{PdH}_{x}$ phases) can be fully
recovered to form a single crystal as the phase propagation proceeds. This
further confirms that the formation of a unique expanded interface results
from the lattice mismatch at the interface.

We conclude that the expanded interface with mixed $\alpha$\- and
$\beta$-$\mathrm{PdH}_{x}$ phases accompanied by lattice tilting allows the
crystal to compensate for the lattice mismatch strain during hydrogen
absorption due to the fast phase propagation kinetics.

## Dislocation motion at interfaces in slow kinetics

Figure 3: Analysis of the atomic structure, strain, and misfit dislocation
evolution at interfaces between $\alpha$\- and $\beta$-$\mathrm{PdH}_{x}$
during hydrogen desorption with slow phase propagation kinetics. a, High-
resolution TEM image of a Pd nanocrystal during hydrogen desorption. The
corresponding FFT pattern (inset) shows the presence of two phases. b, The
calculated (200) $d$-spacing map of the Pd nanocrystal. c, GPA strain map
overlaid on the magnified high-resolution image, which corresponds to the
selected region in (a), demonstrates the sharp interface boundary between two
smooth contrast regions of $\alpha$-$\mathrm{PdH}_{x}$ and
$\beta$-$\mathrm{PdH}_{x}$ in the [100] direction ($\mathrm{E}_{yy}$). We set
the bottom $\beta$-$\mathrm{PdH}_{x}$ area as reference for this GPA for the
desorption case. d, Bragg filtering of the (020) reflection reveals
dislocations at the interface in three consecutive frames during the hydrogen
desorption reaction. The dislocation moves in both [100] and [010] directions.
e, Schematic representation of the atomic structure evolution of the crystal
lattice across the interface based on the experimental observation as
displayed in (d). The propagation of the edge dislocation along with the
interface is highlighted.

A representative high-resolution image of a Pd nanocube during hydrogen
desorption is shown in Fig. 3a (from Supplementary Video 4). The co-existence
of two phases is presented in the splitting {200} and {220} spots in the FFT
pattern (Fig. 3a inset and Supplementary Fig. 18). The calculated $d$-spacing
map of the (200) lattice planes clearly shows a sharp interface between the
$\beta$\- and $\alpha$-$\mathrm{PdH}_{x}$ phases (Fig. 3b and Supplementary
Fig. 18). The GPA ($\mathrm{E}_{yy}$) strain map (Fig. 3c), superimposed on
the interface region of the corresponding HRTEM image, reveals the sharp
interface separating the $\beta$\- and $\alpha$-$\mathrm{PdH}_{x}$. In the
$\mathrm{E}_{xx}$ direction, the GPA data shows no significant strain
variations at the sharp interface (Supplementary Fig. 19). A closer
examination of the atomic arrangement at the interface reveals misfit
dislocations (Fig. 3d,e). The inverse FFT images using the (020) spot (Fig.
3d) and the corresponding schematic representation (Fig. 3e) illustrate the
positions of an edge dislocation in three consecutive frames. The analysis
demonstrates the dynamic motion of misfit dislocation in both [100] and [010]
directions at the propagating sharp interface between $\beta$\- and
$\alpha$-$\mathrm{PdH}_{x}$ during hydrogen desorption.

## Simulations of kinetically tuned interface configurations

Figure 4: Simulations of the dynamic phase transformations of
$\alpha$/$\beta$-$\mathrm{PdH}_{x}$ system under different phase propagation
speeds using an elastic Ising model. a, Representative configurations from
interfacial propagation trajectories under $\Delta\mu$ of 0.1 and 1.0
$k_{\mathrm{B}}T$, in rods (16 × 16 × 50 unit cells): all unit cells and only
those with the local value of $c_{\beta}>0.5$. $\delta h$ represents the
instantaneous interfacial width in units of ‘$a$’, the model’s lattice
constant. Propagation videos are Supplementary Videos 5 and 6. b, The
interfacial regions (dotted boxes in (a)). Color indicates the location along
the $z$-axis: zero set at the topmost (001) plane, with all unit cells having
the local value of $c_{\beta}>0.5$. c, Averaged time-series of the interface
height for $\Delta\mu$ of 0.1 and 1.0 $k_{\mathrm{B}}T$, with one standard
deviation shaded. d,e, The average velocity of the interfacial surface in the
$z$-direction (d) and its characteristic width (e) under increasing
$\Delta\mu$. Error bars display two standard errors. (c-e) Angled brackets
represent an ensemble average over at least 115 independent trajectories. The
subscript zero denotes the average under the ensemble with $\Delta\mu=0.1$
$k_{\mathrm{B}}T$. When multiplying $\Delta\mu$, $\beta$ is the inverse
temperature.

We employ an elastic Ising model to understand the observed differences in
interface configurations under different phase propagation speeds. Our
simulations are conducted in a rod-like geometry with one end of the rod held
in $\alpha$-$\mathrm{PdH}_{x}$ phase and the other in
$\beta$-$\mathrm{PdH}_{x}$ phase[39]. In every trajectory, a rod of
$\alpha$-$\mathrm{PdH}_{x}$ is transformed into a rod of
$\beta$-$\mathrm{PdH}_{x}$ with an interface that propagates in the positive
$z$-direction. Kinetic Monte Carlo dynamics are used to evolve the system,
with rates that include diffusion of $\beta$-phase within the rod, and the
absorption and desorption of $\beta$-phase only at the interface. Within this
framework of our constrained dynamic model, this rod can be interpreted as the
central part of a nanocube that possesses a greater spatial extent than the
rod geometry (see Supplementary Note 5 and Extended Data Fig. 4).

Within the model, the chemical potential difference, $\Delta\mu$, between the
$\mathrm{H}_{2}$ molecules in solution and the H atoms in the Pd lattice sets
the driving force for interface propagation (see Supplementary Note 5). The
configurations shown in Fig. 4a,b represent typical interfacial structures
from trajectories under a chemical driving force of $\Delta\mu$ = 0.1
$k_{\mathrm{B}}T$ and $\Delta\mu$ = 1.0 $k_{\mathrm{B}}T$. When the chemical
driving force is small compared to the thermal energy, diffusion occurs more
quickly than absorption, and a narrower boundary between the two phases is
maintained during propagation. However, when the chemical driving force is
comparable to the thermal energy, absorption occurs too quickly for the
surface to relax, and the boundary between the two phases becomes expanded.
The instantaneous interface between the two phases can be defined more
quantitatively as the isosurface where the concentration of
$\beta$-$\mathrm{PdH}_{x}$ at coexistence ($c_{\beta}$) is 0.5 (see
Supplementary Note 5). Fig. 4c shows a time series of the average position of
this instantaneous interface in the rod, averaged over an ensemble of
independent trajectories for both $\Delta\mu$ = 0.1 $k_{\mathrm{B}}T$ and
$\Delta\mu$ = 1.0 $k_{\mathrm{B}}T$. Under each chemical driving force, we
measure the average speed (see Supplementary Note 5) of the interface
propagation, which increases monotonically with the driving force as expected
(Fig. 4d). The increase in the average speed of propagation across the range
of $\Delta\mu$ = 0.1-1.0 $k_{\mathrm{B}}T$ is approximately one order of
magnitude, akin to the observed difference in propagation speeds between
absorption and desorption in our liquid phase TEM experiments. The width of
the interface is computed from the mean squared deviations of the
instantaneous interface along the direction of propagation, which is
representative of the average size of the interphase region (see Supplementary
Note 5). Fig. 4e shows that as the interface propagates more quickly it also
becomes expanded, increasing by approximately 30% for the largest value of
$\Delta\mu$ shown. In Fig. 4b, it is apparent that a typical configuration in
the expanded regime ($\Delta\mu$ = 1.0 $k_{\mathrm{B}}T$) contains more (001)
planes, which have both $\alpha$\- and $\beta$-phase unit cells. While the
fact that our model contains only two components prevents the direct
observation of mixed phases, it may be reasonable to assume that such
incomplete layers would represent mixed phases within a more flexible model.
The characterization of the degree of mixing between the $\alpha$\- and
$\beta$-phases warrants further study with a finer-grained modeling approach.
This model also includes only coherent interfaces. Thus, the experimentally
observed expanded interface upon hydrogen absorption, relative to the slower
desorption, can be attributed to a kinetic roughening[40, 41] and associated
purely with faster propagation.

Our simulations of the dynamic phase transformations with a minimal model
provide critical insights into the different interfacial structures observed
in our liquid phase TEM experiments. The fast phase transformation pushes the
system far-from-equilibrium, and the interface expansion and imperfection,
analogous to interface roughening, is a generic route to accommodate mismatch
strain under fast transformation kinetics. Under slow transformations, smooth
interfaces are observed, suggesting the system is near-equilibrium.

## Conclusion

By utilizing liquid phase TEM, we have achieved unprecedented atomic
resolution imaging of hydrogen-induced phase transformations in Pd
nanocrystals, which allows us to reveal distinct atomic pathways for mismatch
strain relaxation at interfaces under phase transformation kinetics differing
by about an order of magnitude. We have discovered a generic route to
compensate for mismatch strain by expansion and imperfection ($e.g.$, mixing
of phases, lattice tilting, or other defects) of the interfaces under fast
phase transformations in a far-from-equilibrium regime. This extends our
conventional understanding of dislocation formation as the standard mechanism
to release mismatch strain at interfaces. The expanded interface is configured
as the lattice tilting up to ±5° near the phase propagation front, two-phase
mixture, and other defects. This is distinctly different from the near-
equilibrium phase transformation, where the interface is sharp, accompanied by
the misfit dislocation movements at the interfaces. Our findings provide
valuable insights into engineering materials with desired structures and
functionalities. For example, through far-from-equilibrium material processing
techniques, we can synthesize materials rich in lattice imperfections, thus
providing abundant catalytic active sites. Conversely, near-equilibrium
processing may create defect-free materials more effective for hydrogen
storage and purification, or electrochemical ion insertion devices. Moreover,
our capability to directly observe atomic-level dynamics in systems far-from-
equilibrium opens future opportunities for identifying metastable structures
and leveraging intermediate phases to tune kinetic pathways across a broad
range of significant physicochemical processes of materials.

###### Acknowledgements.

This work was supported by the U.S. Department of Energy, Office of Science,
Office of Basic Energy Sciences (BES), Materials Sciences and Engineering
Division under Contract No. DE-AC02-05-CH11231 within the in-situ TEM program
(KC22ZH). Work at the Molecular Foundry of Lawrence Berkeley National
Laboratory (LBNL) was supported by the U.S. Department of Energy under
Contract No. DE-AC02-05CH11231. S.B.B. acknowledges financial support from the
Alexander-von-Humboldt Association. D.L. acknowledges the Kwanjeong Study
Abroad Scholarship from the KEF (Kwanjeong Educational Foundation) (KEF-2019).
S.B.B. thanks Prof. A. Paul Alivisatos for discussions. We thank Dr. Karen C.
Bustillo for her help with the electron microscope setup.

### Author Contributions

S.B.B. did the liquid phase TEM experiments. S.B.B., D.L., and C.O. performed
the experimental data analysis. S.O.L. carried out the theoretical modeling
and simulations and analyzed results under the supervision of D.T.L. This
paper was written by S.B.B., D.L., S.O.L., D.T.L., and H.Z. with input from
all co-authors. H.Z. conceived and supervised the project.

### Competing Interests

Authors declare that they have no competing interests.

### Correspondence

Correspondance and requests for materials should be addressed to H. Zheng
<EMAIL_ADDRESS>or D. T. Limmer (dlimmer@berkeley.edu).

## Methods

### Synthesis of palladium nanocubes

Palladium nanocubes with sizes between 40 and 50 nm are synthesized using a
reported two-step synthesis strategy[25]. In the first step, the seed solution
is synthesized: 91.1 mg cetrimonium bromide (CTAB from Sigma-Aldrich) is
dissolved in 20 mL $\mathrm{H}_{2}\mathrm{O}$ by stirring (400 rpm) at 95 °C
for 5 minutes. A 10 mM $\mathrm{H}_{2}\mathrm{PdCl}_{4}$ solution is prepared
by dissolving 10.7 mg anhydrous $\mathrm{PdCl}_{2}$ (Sigma-Aldrich) in 5 mL of
20 mM HCl (Sigma-Aldrich). 1 mL of this 10 mM
$\mathrm{H}_{2}\mathrm{PdCl}_{4}$ solution is added to the CTAB solution,
which subsequently turns red. The mixture is stirred for 5 minutes before
adding 160 $\mu$L of a 100 mM L-ascorbic acid solution (created by dissolving
17.6 mg L-ascorbic acid from Sigma-Aldrich in 1 mL
$\mathrm{H}_{2}\mathrm{O}$). The final mixture is stirred for 5 min at 95 °C
before being placed in an ice bath to stop the reaction.

In the second step, the nanocube size is increased to 40 and 50 nm from the
seeds; for this, a 0.1 M CTAB solution is prepared by dissolving 364.5 mg CTAB
in 10 mL $\mathrm{H}_{2}\mathrm{O}$. This CTAB solution is heated to 60 °C
before adding 125 $\mu$L of the 10 mM $\mathrm{H}_{2}\mathrm{PdCl}_{4}$
solution, 100 $\mu$L of the seed solution, and 50 $\mu$L of the 100 mM
L-ascorbic acid solution. The mixture is slightly shaken and allowed to sit
for 2 hours at 60 °C. The final product is collected by centrifugation at 6000
rpm for 5 minutes and washed twice with 5 mL $\mathrm{H}_{2}\mathrm{O}$ using
a centrifugation rotation speed of 10000 rpm. Overview images, including
bright-field transmission electron microscopy (TEM) and energy-dispersive
X-ray spectroscopy (EDS) map images, of synthesized nanocubes are displayed in
Supplementary Fig. 1. Detailed crystallographic properties of a synthesized Pd
nanocube are demonstrated in Supplementary Fig. 2.

### High-resolution liquid cell fabrication

The liquid cells used in this study are based on Formvar-coated TEM grids. The
Formvar films are fabricated according to the following protocol. First, a
microscope glass slide is dipped into a 0.25% Formvar solution in ethylene
dichloride (Electron Microscopy Sciences), resulting in a thin film coating
the glass slide (Supplementary Fig. 3a). The film is dried in the air for 5
minutes. Subsequently, the surface tension of water is exploited to separate
the Formvar film from the glass slide (Supplementary Fig. 3b). The Formvar
film floats on the water’s surface, and Cu TEM grids (Electron Microscopy
Sciences) are dropped onto the film’s surface (Supplementary Fig. 3c).
Finally, the Formvar-coated grids are collected with a second glass slide and
dried at ambient conditions (Supplementary Fig. 3d). Electron energy-loss
spectroscopy demonstrates a reproducible thickness of 17 $\pm$ 3 nm for the
synthesized Formvar films.

Then, 1 $\mu$L of the aqueous solution containing palladium nanocubes is drop-
casted onto the Formvar-coated TEM grids. The TEM grids are subsequently
plasma-cleaned in an Argon-Oxygen atmosphere (Fischione Model 1020 plasma
cleaner) for 8 seconds. The plasma cleaning makes the films hydrophilic,
reduces the film thickness by about 40%, and removes the residual organic
ligands from the nanoparticle surface. Then, 0.5 $\mu$L of 0.1 M KOH (Sigma-
Aldrich) is dropped onto one TEM grid. A filter paper is used to remove the
excess liquid so that only a thin liquid layer covers the grid. Subsequently,
a second TEM grid is placed, and the liquid cell is sealed by the van-der-
Waals attraction between the films. Liquids are trapped within small pockets
formed between the membranes of the two grids.

### Transmission electron microscopy

The FEI ThemIS 60-300 (Thermo Fisher Scientific) transmission electron
microscope with an X-FEG gun operated at 200 and 300 kV was used for in situ
TEM experiments. The microscope is equipped with an image aberration corrector
(spatial resolution limit of 70 pm), Bruker SuperX energy-dispersive X-ray
spectroscopy (EDS) detector (solid angle of 0.7 steradians), and a Ceta2 CMOS
camera. The FEI F20 Tecnai microscope was used to record electron energy loss
spectroscopy data.

### Determination of the $\mathrm{\textbf{H}}_{\mathbf{2}}$ concentration as a
function of the electron dose rate

A quasi-stationary model[42] is used to describe the growth/dissolution of gas
bubbles in supersaturated/undersaturated solutions. The model describes the
changes in the bubble radius $r(t)$ due to the mass transport of the gaseous
species into/out of the bubble, neglecting the surface tension. Bubbles grow
if the concentration of the gaseous species dissolved in a liquid ($c_{i}$)
exceeds its saturation concentration in the liquid ($c_{s}$)[43]. By depleting
the supersaturated $c_{i}$, the bubbles grow until the gas concentration
inside the liquid reaches its saturation concentration $c_{s}$. We track the
changes in the bubble radius with time as a function of the electron dose rate
to establish an understanding of the $\mathrm{H}_{2}$ and $\mathrm{O}_{2}$
concentrations during the in situ experiments. An exemplary TEM image is shown
in Supplementary Fig. 4a highlighting the identified gas bubbles. Ellipses are
fitted to the bubbles in the image to account for the non-spherical shape of
some bubbles. The spherical radius $r_{s}$ of the bubble is related to the two
radii $r_{1}$ and $r_{2}$ characterizing an ellipse (see Supplementary Fig.
4b): $r_{s}=\sqrt{r_{1}r_{2}}$. Dynamic changes of individual bubbles are
reflected by a large size distribution within one TEM image. Supplementary
Fig. 4c illustrates the radii detected at distinct times during an in situ
liquid phase experiment. Considering that the supersaturation of the solution
is reflected by the size of the largest bubbles, we track the average of the
25% largest bubbles as a function of the dose rate or illumination time
(highlighted in red in Supplementary Fig. 4c). The analysis of the radius
change as a function of the reaction time at a constant electron dose rate
demonstrates that the bubble radius approaches a steady-state size
(Supplementary Fig. 4d), which depends on the electron dose rate used during
the experiment (Supplementary Fig. 5). The time for establishing the steady-
state concentration varies randomly and does not follow the expected trend;
for example, the higher the dose rate, the shorter the time. This is likely
due to the different sizes of the investigated liquid pockets in the liquid
cell. Since only a small subsection of the liquid cell is illuminated by the
beam during liquid phase TEM experiments, the mass transport outside the
illuminated area impacts the steady-state concentrations[28, 31].

### High-resolution TEM (HRTEM) image analysis – a rolling window algorithm

As shown in Supplementary Fig. 6, palladium’s lattice expansion due to
hydrogen absorption is represented in the changes in radial distances of
reflections in the diffractogram. In other words, we can differentiate
$\alpha$\- and $\beta$-$\mathrm{PdH}_{x}$ phases during hydrogen-induced phase
transformations based on the differences in $d$-spacing values. To trace local
variations of the $d$-spacing within the palladium nanocubes during the phase
transformations, we developed a rolling window algorithm inspired by previous
work[44] (Supplementary Fig. 7). This algorithm utilizes dynamically moving
windows over HRTEM images, which generate fast Fourier transformation (FFT)
patterns corresponding to local regions of Pd nanocubes, determining the
$d$-spacing variations ($i.e.$, the spatial distributions of $\alpha$\- and
$\beta$-$\mathrm{PdH}_{x}$ phases) based on the FFT patterns. Fitting 2D
Gaussian functions to the FFT pattern allows the determination of peak
positions with sub-pixel accuracy. The determined $d$-spacing values are used
to generate maps visualizing the variations in the lattice arrangements within
the nanoparticles. The fitted patterns with variations more significant than
the angle variation of ±5° and the $d$-spacing variation of $\pm$0.1 Å are
discarded. The window size of the FFT pattern determines the spatial
resolution, and the values of 10×10 pixels are averaged to generate the maps.
Additionally, a mask segmenting the nanoparticle from the amorphous background
is used to remove artifacts created by the amorphous structure.

### Geometric phase analysis (GPA)

GPA was used to investigate strain at interphase boundaries within palladium
nanocubes during hydrogen absorption and desorption reactions. The GPA method
measures the displacement of lattice fringes shown in HRTEM images of
palladium nanocubes during the hydrogen-induced phase transformations. We
select the parts of nanocrystals below their interphase boundaries as
references. In other words, we choose the $\alpha$-$\mathrm{PdH}_{x}$ and
$\beta$-$\mathrm{PdH}_{x}$ as references for the absorption and desorption
cases, respectively. Then, compute the differences in $d$-spacing values
throughout an HRTEM image to generate a corresponding strain map. Finally, the
strain distributions are revealed in Figs. 2d and 3c and Supplementary Figs.
17 and 19.

### Data and Code Availability

All data and codes that support the findings of this study are available from
the corresponding authors upon reasonable request.

## Extended Data

Extended Data Fig. 1. The nucleation and propagation of the
$\beta$-$\mathrm{PdH}_{x}$ phase in Pd nanocubes during hydrogen absorption.
a, HRTEM image of a Pd nanocube at the initial state of the hydrogen
absorption reaction. The respective FFT pattern shows the presence of two
crystal phases ($\alpha$ and $\beta$ phases) as indicated by a splitting of
the (-220) diffraction spots. b, The analysis of the spatial variation of the
$d$-spacing of (-220) reveals the local lattice expansion at the upper right
corner of the nanocube. The local lattice expansion induces strain and a local
thickness variation, resulting in the dark fringes in the respective TEM
image. c, Phase transformation dynamics determined based on the relative area
of the nanocube covered by the dark fringes ($i.e.$, the relative area of the
$\beta$-$\mathrm{PdH}_{x}$) and the bubble radius, which is proportional to
hydrogen concentration (Supplementary Video 1). The reaction happens in three
stages, highlighted in purple, yellow, and red. The displayed TEM images are
characteristic of the nanocube during the three different stages (designated
by colored arrows). The experiment time is given relative to the initial data
acquisition. Extended Data Fig. 2. Free energy curves for hydrogen absorption
and configurations along a reversible absorption pathway in the model. a, The
Helmholtz free energy as a function of $c_{\beta}$ computed from standard
umbrella sampling techniques and the weighted histogram analysis method. All
simulations are run at temperature (T) = 300 K (complete simulation details in
Supplementary Note 3). b, The same free energy curves as in (a), but with
linear shifts due to chemical potential differences, $\Delta\mu$, between
$\mathrm{H}_{2}$ molecules in the surrounding solution and H atoms in the
palladium lattice. One sees the corner nucleated configuration, (2), is a
representative critical nucleus when chemical driving for absorption is on the
order of one-tenth the thermal energy. Note that, when multiplying an energy,
$\beta$ represents the inverse temperature. c,d, Two perspectives on
representative configurations taken from umbrella sampling trajectories, with
$c_{\beta}$ increasing in each configuration from left to right. e, The same
configurations, (2) and (3) from (c and d), where only unit cells with a local
value of $c_{\beta}$ greater than 0.5 are shown to visualize the full
interphase boundaries. Extended Data Fig. 3. An energetic analysis of
idealized configurations along hypothetical absorption pathways. a, The
elastic energies (Supplementary Equation 2 in Supplementary Notes 2.2) of
idealized (100), (110), and (111) interfacial configurations as a function of
$c_{\beta}$, at three different sizes of nanocube. $N$ is the total number of
sites, and $\beta$ represents the inverse temperature when multiplying an
energy. b,c,d, Selected configurations with labels corresponding to the
markers in the third panel of (a). Extended Data Fig. 4. Analysis of the size
dependence of the interfacial propagation in rod geometries. a, A schematic of
how a rod can be viewed as the core of a larger nanocube. b,c,d, $L_{s}$
denotes the width of the rod in the $x$-$y$ dimensions in units of ‘$a$’, the
model’s lattice constant. The length in the $z$-direction is 50 unit cells for
all data. b Average interfacial widths under increasing chemical driving
force. c Average velocities under increasing chemical driving force. (b and c)
Averages recorded over the range $\bar{h}/a=25$ to $\bar{h}/a=31$. d Time
series of the width under three chemical driving forces. $\Delta t$ is
recorded from the time the interface reaches an average height of
$\bar{h}/a=10$. The zero subscript denotes the average under the ensemble with
$\Delta\mu=0.1k_{\mathrm{B}}T$.

## Supplementary Notes

##  Hydrogen absorption pathway in palladium nanocubes

Extended Data Fig. 1a shows a high-resolution TEM image of a Pd nanocube
recorded during the initial period of the hydrogen absorption reaction. The
FFT pattern shows the co-existence of $\alpha$\- and
$\beta$-$\mathrm{PdH}_{x}$ phases indicated by the splitting of the (-220)
diffraction spots (as described in Supplementary Fig. 6c). High-resolution
image analysis with the rolling window algorithm tracking the spatial
distribution of the two crystal phases in the Pd nanocube (see Methods and
Supplementary Fig. 7) demonstrates that the lattice expansion
($\beta$-$\mathrm{PdH}_{x}$ phase) is localized at the upper right corner of
the nanocube (Extended Data Fig. 1b), coinciding with the dark fringes in the
image. Similar fringes were previously observed in the dark-field scanning TEM
(STEM) experiments on hydrogen absorption in a gas atmosphere[22]. These
fringes result from either the strain induced by the lattice mismatch between
$\alpha$ and $\beta$-$\mathrm{PdH}_{x}$ or the thickness increase induced by
the local swelling of the crystal lattice.

We utilize the fringes to identify the hydrogen uptake into palladium
nanocubes at lower magnification (Supplementary Video 1). Extended Data Fig.
1c shows the kinetics of the phase transformation induced by the hydrogen
absorption in a Pd nanocube, quantifying the area of the dark fringes. During
this in situ experiment, the electron flux was stepwise increased to
manipulate the $\mathrm{H}_{2}$ concentration. The hydrogen concentration can
be correlated with the average bubble radius (see Methods and Supplementary
Figs. 4 and 5). Coupled with the time-series data of the average bubble
radius, the analysis of the proportion of $\beta$-$\mathrm{PdH}_{x}$
(calculated by the projected area of dark fringes) with time demonstrates
three stages of hydrogen absorption in the palladium nanocrystal (Extended
Data Fig. 1c).

From the nucleation of $\beta$-$\mathrm{PdH}_{x}$ at the corner of the
palladium nanocube to 25% of the nanocube being transformed into
$\beta$-$\mathrm{PdH}_{x}$ (Stage 1), the $\beta$-$\mathrm{PdH}_{x}$ is
localized at the corner. Increasing the electron dose rate promotes the
hydrogen absorption reactions, since more hydrogen can be produced, as
indicated by the increased bubble radius. After a short reaction stagnation,
the phase propagation proceeds rapidly along the [010] direction until about
50% of the nanocube is transformed to $\beta$-$\mathrm{PdH}_{x}$ (Stage 2).
The phase transformation slows again until the rapid phase propagation along
the [100] direction completes the phase transformation within a few seconds
(Stage 3). Although our analysis is limited to the 2D projections of three-
dimensional nanocubes, it is evident that the phase transformation from Stage
1 to Stage 2 and Stage 3 requires a change of the phase propagation direction.
Our experimental data demonstrate the corner nucleation and [100] propagation
of $\beta$-$\mathrm{PdH}_{x}$ during hydrogen absorption in palladium
nanocubes, which aligns with the results reported in preceding studies[22,
25].

Supplementary Fig. 1. S/TEM characterizations of as-synthesized palladium
nanocubes. a, Bright-field TEM image of several palladium nanocubes on an
amorphous carbon film synthesized following the described procedure. b, STEM
image of palladium nanocubes overlaid with EDS maps (left) and corresponding
Pd-$L$ and C-$K$ edges EDS maps (right), which confirm the purity of the
palladium nanocubes and the absence of a carbon ligand shell. Supplementary
Fig. 2. Crystallographic characteristics of a synthesized palladium nanocube.
HRTEM image (left) of a palladium nanocube with a 40 to 50 nm size range,
exposing the {100} surface facets. The corresponding FFT pattern and magnified
HRTEM image (right) exhibit the single-crystalline nature of the nanocube.
Supplementary Fig. 3. Schematic representation of the fabrication of Formvar-
coated TEM grids. a, A glass slide is dipped into a 0.25% Formvar solution,
resulting in a thin Formvar film coating the glass slide. b, The surface
tension of water is used to separate the Formvar film from the glass slide. c,
The Formvar film floats on water, and the TEM grids are placed on top of it.
d, A second glass slide is used to catch the coated TEM grids. Supplementary
Fig. 4. The $\mathrm{H}_{2}$ concentration inside the liquid cell is
indirectly estimated based on the size of the bubbles formed inside the
liquid. a, A representative TEM image used for the bubble size analysis. The
identified bubbles are highlighted in red. A dose rate of 40
$e^{-}\text{\AA}^{-2}s^{-1}$ was used to record this TEM image. b, Schematic
presentation of the liquid cell conditions. Ellipses are fitted to the
identified bubbles to account for their non-perfect spherical shapes. c,
Scatter plot illustrating the identified bubble radii at different stages of
the liquid phase experiment. The 25% largest radii are highlighted in red.
These radii are used to determine the average radius as a function of reaction
time, which is represented by the red line in the graph. The dose rate was
increased from 10 to 40 $e^{-}\text{\AA}^{-2}s^{-1}$ during the experiment. d,
Change of the average bubble radius as a function of time during an experiment
performed at a constant dose rate of 90 $e^{-}\text{\AA}^{-2}s^{-1}$. The
radii change rapidly during the first 60 s of the illumination and
subsequently approach a steady state. Supplementary Fig. 5. Electron-dose-
rate-dependent gas bubble sizes. The steady-state average radius of the gas
bubbles depends on the electron dose rate. Supplementary Fig. 6. Crystal
structures of palladium and $\beta$-$\mathrm{PdH}_{x}$ and their corresponding
electron diffractograms. a,b, Schematic representations of the crystal
lattices of Pd (a) and $\beta$-$\mathrm{PdH}_{x}$ (b). c, The incorporation of
hydrogen into the lattice causes the lattice expansion of palladium by about
3.5% in the lattice parameter[15, 16]. This expansion is visible in the
simulated diffraction patterns, in terms of radial distances of reflections
for Pd (blue) and $\beta$-$\mathrm{PdH}_{x}$ (red). Please note that given
that the difference in such radial distances between Pd and
$\alpha$-$\mathrm{PdH}_{x}$ is immeasurably small, under our in situ liquid
phase experimental conditions, we can regard Pd as
$\alpha$-$\mathrm{PdH}_{x}$. Supplementary Fig. 7. Overview of the steps
involved in the determination of the local variations of the $d$-spacing
within individual Pd nanoparticles. a, A rolling window algorithm uses
dynamically moving windows (a red box showing an example window) over the
nanocube to determine the position of reflections in generated FFT patterns
from those windows. b, An exemplary FFT pattern generated based on a moving
window, which shows an allowed angle variation and fitting 2D Gaussian
functions to a reflection. c, The initially obtained spatial $d$-spacing
variation within the individual nanoparticle. d, Falsely identified values are
removed by allowing an angle variance of $\pm 5$° and a $d$-spacing variance
of $\pm 0.1$ Å. e, Since the resolution is limited by the window size used
during the analysis, the $d$-spacing values are averaged for several pixel
values. f, Artifacts created by the amorphous background are removed using a
mask tool separating the nanoparticle from the background.

##  Details of elastic Ising model for the $\mathrm{PdH}_{x}$ system

###  Parameterization of the Hamiltonian

The elastic Ising model was used to model the interconversion of the
$\alpha$-$\mathrm{PdH}_{x}$ and $\beta$-$\mathrm{PdH}_{x}$ phases. This simple
model includes harmonic interactions between both nearest and next-nearest
neighbors on a simple cubic lattice. The Hamiltonian is given by a sum over
lattice positions, $\mathbf{R}$,

$\displaystyle
H=\sum_{\mathbf{R}}\Bigg{[}\frac{K_{1}}{4}\sum_{\hat{\mathbf{\alpha}}(\mathbf{R})}[|a\hat{\mathbf{\alpha}}+\mathbf{u}_{\mathbf{R}}-\mathbf{u}_{\mathbf{R}+a\hat{\mathbf{\alpha}}}|-l(\sigma_{\mathbf{R}},\sigma_{\mathbf{R}+a\hat{\mathbf{\alpha}}})]^{2}$
$\displaystyle+\frac{K_{2}}{4}\sum_{\hat{\mathbf{\gamma}}(\mathbf{R})}[|\sqrt{2}a\hat{\mathbf{\gamma}}+\mathbf{u}_{\mathbf{R}}-\mathbf{u}_{\mathbf{R}+\sqrt{2}a\hat{\mathbf{\gamma}}}|-l(\sigma_{\mathbf{R}},\sigma_{\mathbf{R}+\sqrt{2}a\hat{\mathbf{\gamma}}})]^{2}\Bigg{]}$
(1)

The set of unit vectors $\\{\hat{\mathbf{\alpha}}(\mathbf{R})\\}$ points from
site $\mathbf{R}$ to the nearest neighbors of $\mathbf{R}$ on the simple cubic
lattice. For a bulk lattice site, there are six vectors, while generally there
are fewer for sites at the edge of the cube. The set of unit vectors
$\\{\mathbf{\hat{\gamma}}(\mathbf{R})\\}$ point to the next-nearest neighbors
of site $\mathbf{R}$, totaling 12 for a bulk lattice site. The variables
$\\{\mathbf{u}_{\mathbf{R}}\\}$ represent displacements of site $\mathbf{R}$
from its equilibrium position. When site $\mathbf{R}$ is viewed as a unit cell
of the Pd nanocube, $\\{\mathbf{u}_{\mathbf{R}}\\}$ represents a compression
or expansion of the cell. The variables $\\{\sigma_{\mathbf{R}}\\}$ are the
spins variables and take values $\pm 1$, with $+1$ corresponding to a unit
cell that is locally $\beta$-$\mathrm{PdH}_{x}$ and $-1$ corresponding to a
unit cell that is locally $\alpha$-$\mathrm{PdH}_{x}$.
$l(\sigma_{\mathbf{R}},\sigma_{\mathbf{R}^{\prime}})$ is the preferred
distance between two sites, which, due to the expansion of the lattice in the
$\beta$ phase, depends on the values of the spin variables. We take
$l(\sigma_{\mathbf{R}},\sigma_{\mathbf{R}^{\prime}})=a_{0}+\Delta/2+\Delta/4(\sigma_{\mathbf{R}}+\sigma_{\mathbf{R}^{\prime}})$,
where $a_{0}$ is the equilibrium lattice constant of bulk Pd and $\Delta$ is
approximately 3.5% of this value (Supplementary Table 1).

The values for the spring constants, $K_{1}$ and $K_{2}$, were chosen such
that the model could reproduce the response of bulk Pd to a longitudinal
stretch and a pure shear. To do so, we considered a 2 × 2 × 2 cube of sites in
the model. The corners of our cube are connected by springs with stiffnesses
given by $K_{1}$ and $K_{2}$ (Supplementary Fig. 8). We assume our system is
tiled by such cubes, and we account for the fact that the springs of one cube
are shared by the neighboring cubes. This means that the springs with
stiffness $K_{1}$ have a weight of one-fourth, and that the springs of
stiffness $K_{2}$ have a weight of one-half.

If all sites are the same spin, then a configuration with no strain will make
all side lengths of the cube equal to $l$, the lattice constant for that
phase. We apply a stretch and a shear strain to the cube and calculate the
energetic cost according to Supplementary Equation 1. This gives us $\Delta
E_{1}=(K_{1}+2K_{2})\Delta y^{2}/2$ and $\Delta E_{2}=K_{2}\Delta x^{2}/2$,
while elastic theory gives $\Delta E_{1}/l^{3}=c_{11}(\Delta y/l)^{2}/2$ and
$\Delta E_{2}/l^{3}=c_{44}(\Delta x/l)^{2}/2$. Therefore, a consistent choice
for the spring constants are, $K_{1}=l(c_{11}-2c_{44})$, and $K_{2}=lc_{44}$.
Hsu and Leisure[45] have calculated the low-temperature limits of the elastic
moduli of $\mathrm{PdH}_{x}$ in both the $\alpha$ and $\beta$ phases. Since we
have created our mapping considering only the energetic cost rather than the
free energy, we take the 0 K extrapolations for $c_{11}$ and $c_{44}$. Since
we wish to model the regime of coexistence, we take averages of the $\alpha$
and $\beta$ phase values for both the value of $l$ and the elastic moduli.
This procedure provides the parameters listed in Supplementary Table 1.

Supplementary Fig. 8. Depiction of the two applied strains used to parameterize the model Hamiltonian as described in Supplementary Note 2.1. Parameter | Value
---|---
$\Delta$ | 0.1387 Å
$a_{0}$ | 0.3890 Å
$K_{1}$ | 2.2303 eV/$\text{\AA}^{2}$
$K_{2}$ | 1.7307 eV/$\text{\AA}^{2}$
Table 1: Parameters used for the Hamiltonian (Supplementary Equation 1) based
on experimental data from the reference.[45]

###  Small mismatch approximation and effective interaction potential

By linearizing the absolute value terms in Supplementary Equation 1, one
obtains a Hamiltonian in the limit of a small lattice mismatch between the two
phases[46] that is quadratic in the displacement variables. We can make the
further assumption that vibrations of the lattice relax quickly compared to
typical timescales of diffusion of hydrogen between lattice sites. This allows
us to analytically integrate out the displacement variables, yielding an
effective interaction potential among only the spin variables[35], of the
form,

$\displaystyle
H_{\mathrm{eff}}=\frac{1}{2}\sum_{\mathbf{R},\mathbf{R}^{\prime}}\sigma_{\mathbf{R}}\mathbf{V}_{\mathbf{R},\mathbf{R}^{\prime}}\sigma_{\mathbf{R}^{\prime}}$
(2)

where spins are coupled by an effective long-range potential. The spin-spin
potential, $V_{\mathbf{R},\mathbf{R}^{\prime}}$, is given by a sum of two
terms,

$\displaystyle\mathbf{V}_{\mathbf{R},\mathbf{R}^{\prime}}=\Bigg{(}\mathbf{S}-\frac{1}{4}\mathbf{C}^{T}\cdot\mathbf{D}^{-1}\cdot\mathbf{C}\Bigg{)}_{\mathbf{R},\mathbf{R}^{\prime}}$

a term,

$\displaystyle\mathbf{S}_{\mathbf{R},\mathbf{R}^{\prime}}=\frac{K_{1}\Delta^{2}}{16}\sum_{\hat{\mathbf{\alpha}}(\mathbf{R})}(\delta_{\textbf{
R
},\mathbf{R}^{\prime}}+\delta_{\mathbf{R}+a\hat{\mathbf{\alpha}},\mathbf{R}^{\prime}})+\frac{K_{2}\Delta^{2}}{8}\sum_{\hat{\mathbf{\gamma}}(\mathbf{R})}(\delta_{\textbf{
R
},\mathbf{R}^{\prime}}+\delta_{\mathbf{R}+\sqrt{2}a\hat{\mathbf{\gamma}},\mathbf{R}^{\prime}})$

and a product term,

$\displaystyle\mathbf{C}_{\mathbf{R},\mathbf{R}^{\prime}}=\frac{K_{1}\Delta}{2}\sum_{\hat{\mathbf{\alpha}}(\mathbf{R})}\hat{\mathbf{\alpha}}(\delta_{\textbf{
R
},\mathbf{R}^{\prime}}+\delta_{\mathbf{R}+a\hat{\mathbf{\alpha}},\mathbf{R}^{\prime}})+\frac{\sqrt{2}K_{2}\Delta}{2}\sum_{\hat{\mathbf{\gamma}}(\mathbf{R})}\hat{\mathbf{\gamma}}(\delta_{\textbf{
R
},\mathbf{R}^{\prime}}+\delta_{\mathbf{R}+\sqrt{2}a\hat{\mathbf{\gamma}},\mathbf{R}^{\prime}})$
$\displaystyle\mathbf{D}_{\mathbf{R},\mathbf{R}^{\prime}}=K_{1}\sum_{\hat{\mathbf{\alpha}}(\mathbf{R})}\hat{\mathbf{\alpha}}\hat{\mathbf{\alpha}}(\delta_{\textbf{
R
},\mathbf{R}^{\prime}}-\delta_{\mathbf{R}+a\hat{\mathbf{\alpha}},\mathbf{R}^{\prime}})+K_{2}\sum_{\hat{\mathbf{\gamma}}(\mathbf{R})}\hat{\mathbf{\gamma}}\hat{\mathbf{\gamma}}(\delta_{\textbf{
R
},\mathbf{R}^{\prime}}-\delta_{\mathbf{R}+\sqrt{2}a\hat{\mathbf{\gamma}},\mathbf{R}^{\prime}})$

requiring the inversion of the matrix
$\mathbf{D}_{\mathbf{R},\mathbf{R}^{\prime}}$.

###  Kinetic Monte Carlo algorithm

For all simulations, the standard Gillespie algorithm[47] was used to
propagate the dynamics. The rate to go from configuration $\mathbf{C}$ to
$\mathbf{C}^{\prime}$ took the form,

$\displaystyle
k_{\mathbf{C},\mathbf{C}^{\prime}}=k_{0}\exp\Bigg{\\{}-\frac{\beta}{2}\bigg{(}E(\mathbf{C}^{\prime})-E(\mathbf{C})\bigg{)}+\frac{\beta\Delta\mu}{2}\bigg{(}\big{(}N(\mathbf{C}^{\prime})-N(\mathbf{C})\big{)}\bigg{)}\Bigg{\\}}$
(3)

$E(\mathbf{C})$ and $N(\mathbf{C})$ respectively denote the energy and the
number of spin up sites in configuration $\mathbf{C}$. The pre-factor,
$k_{0}$, was constant for all processes, and it sets the units of time for the
simulation, $k_{0}^{-1}=\tau_{0}$. The allowed processes were different for
the umbrella sampling and interfacial propagation simulations, which will be
described in Supplementary Notes 2.4, 3.1 and 5.2 .

###  Phase diagrams of the model

Supplementary Fig. 9 shows four approximate phase diagrams for nanocubes with
side lengths of 10, 15, 20, and 30 sites. At each temperature, at least four
independent trajectories were initialized in random configurations with
average values of
$c_{\beta}=\sum_{\mathbf{R}}\frac{\sigma_{\mathbf{R}}+1}{2N}$ of 0.75 and
0.25. The trajectories were propagated for at least 500 rejection-free kinetic
Monte Carlo sweeps, and up to $10^{4}$ sweeps in the case of the smaller
cubes. The allowed moves were spin flips and nearest neighbor exchanges. The
average value of $c_{\beta}$ was computed along each trajectory after
discarding 200 sweeps for equilibration. The two markers for each temperature
in Supplementary Fig. 9 are the averages across the independent trajectories
beginning near $c_{\beta}=0.25$ and $c_{\beta}=0.75$, respectively. The error
bars are three standard errors in the mean across the independent
trajectories. The dashed lines are quartic fits to the markers, where the
fourth-order coefficient is enforced to be negative.

For a system of infinite extent, one can obtain an analytic form for the
interaction potential of Supplementary Equation 2, in a Fourier basis[48].
Using a numerical Fourier transform, we can construct a Hamiltonian of the
form,

$\displaystyle
H_{\mathrm{eff}}=\frac{1}{2}\sum_{\mathbf{R},\mathbf{R}^{\prime}\neq\mathbf{R}}\sigma_{\mathbf{R}}V^{\mathrm{bulk}}(\mathbf{R}-\mathbf{R}^{\prime})\sigma_{\mathbf{R}^{\prime}}$
(4)

The simplest mean-field theory for this Hamiltonian yields the self-consistent
equation,

$\displaystyle\tanh(\beta\bar{V}m)=m$ (5)
$\displaystyle\bar{V}=\sum_{\mathbf{R}^{\prime}\neq\mathbf{R}}V^{\mathrm{bulk}}(\mathbf{R}-\mathbf{R}^{\prime}),\hskip
5.69054ptm=\frac{1}{N}\sum_{\mathbf{R}}\sigma_{\mathbf{R}}$

Calculating $\bar{V}$ out to 256 lattice sites in all three Cartesian
directions, we obtain a bulk, mean-field critical point,
$T_{c}^{\mathrm{MF}}=644.7\mathrm{K}$. The solid line in Supplementary Fig. 9
is the numerical solution to Supplementary Equation 5 with this value of
$\bar{V}$. All simulations are run at 300 K.

Supplementary Fig. 9. Phase diagrams for the elastic Ising model, for
different size nanocubes (with side length, $L$). Markers are averages over at
least 4 independent trajectories of lengths of at least 500 rejection-free KMC
sweeps. Dashed lines are quartic fits to the markers. The solid line is the
numerical solution to an equation of state under a simple mean-field
approximation in a bulk system.

##  Free energy calculations

###  Computational details

Umbrella sampling calculations were carried out using $c_{\beta}$ as the order
parameter on cubes of side length, 15 sites or approximately 6 nm. 415
umbrellas were used with centers between $c_{\beta}^{*}=-0.05$ and
$c_{\beta}^{*}=1.05$. The biasing potential was of the form
$U=KN(c_{\beta}-c_{\beta}^{*})^{2}/2$. Three spring constants were used in
different ranges of $c_{\beta}^{*}$, which are summarized in Supplementary
Table 2. Trajectories were run with kinetic Monte Carlo (KMC) dynamics, with
spin flips allowed at the boundary of the nanocube and exchanges allowed
between any nearest-neighbor sites. 200 rejection-free sweeps were discarded
for equilibration in each trajectory. The correlation time for the
concentration was assumed to be less than five rejection-free sweeps in each
trajectory. Trajectories were between $5\times 10^{3}$ and $5\times 10^{4}$
rejection-free sweeps in length. The weighted histogram analysis method
(WHAM)[49] was used to construct the free energy profile with a tolerance of
$10^{-3}$. To estimate the statistical error, trajectories were split into
five segments of equal length, and separate WHAM calculations were performed
on each data set. Three standard errors in the mean of the five independent
calculations are smaller than the line widths of Extended Data Fig. 2a,b.

$\mathbf{c_{\beta}^{*}}$ range | Spring constant (eV)
---|---
$-0.05\text{ to }0.09$ | 0.8402
$0.10\text{ to }0.40$ | 0.6463
$0.40\text{ to }0.60$ | 0.4524
$0.60\text{ to }0.90$ | 0.6463
$0.90\text{ to }1.05$ | 0.8402
Table 2: Spring constants used for umbrella sampling for umbrellas with
centers in different ranges.

###  The energetics of idealized configurations

To further understand the observed absorption pathway, we studied the elastic
energies of Supplementary Equation 2 from idealized configurations along
potential pathways of absorption. We consider perfect interfaces in the
$\langle 100\rangle$, $\langle 110\rangle$, and $\langle 111\rangle$
directions moving across the cube. In Extended Data Fig. 3, we compute the
energies for different size nanocubes for each of the interfaces. We observe
that the energetic cost to form any interface scales extensively, and that for
each of the sizes the $\langle 100\rangle$ interface is significantly
energetically preferred near $c_{\beta}=0.5$. However, for $c_{\beta}$ below
$0.25$, the most energetically preferred interface is along the $\langle
110\rangle$ direction, corresponding to the corner-nucleus configurations
observed in both the model and in the experiment. Thus, it seems the observed
absorption pathway can be understood mainly in terms of minimizing the cost in
elastic energy due to the interface between the two phases at each value of
$c_{\beta}$. In addition, the energetic cost to form any interface increases
in proportion to system size, while the differences in the energetics among
the $\langle 100\rangle$, $\langle 110\rangle$, and $\langle 111\rangle$
directions remain constant. Therefore, we can expect that the reversible
pathway observed in Extended Data Fig. 2 will be consistent as the nanocube
approaches the size used in the liquid phase TEM experiments.

##  Theoretical simulation of hydrogen absorption pathway in palladium
nanocubes

Simulations of hydrogen absorption and desorption are conducted with kinetic
Monte Carlo of the elastic Ising model Hamiltonian[33, 34, 35] (see
Supplementary Note 2.1-2.3). The elastic Ising model accounts for the
energetic cost of the strain that is induced due to the difference in lattice
constants between regions of the nanocube rich in $\alpha$-$\mathrm{PdH}_{x}$
and those rich in $\beta$-$\mathrm{PdH}_{x}$. In the model, a nanocube is
coarse-grained such that each unit cell is defined to be either
$\alpha$-$\mathrm{PdH}_{x}$ or $\beta$-$\mathrm{PdH}_{x}$ phase, and the
Hamiltonian is parameterized with experimental values of the elastic moduli of
$\mathrm{PdH}_{x}$[45] (see Supplementary Note 2.1) in order to describe the
elastic mediated effective interaction potential between each site.
Supplementary Fig. 9 shows that with this parameterization, the simple model
exhibits a phase diagram that is quantitatively similar to the experimental
diagram in the bulk[50] (see Supplementary Note 2.4). It also captures the
expected effects of nanoscale confinement, in a suppressed critical
temperature and enhanced equilibrium concentration of hydrogen at the
interface in the $\alpha$-$\mathrm{PdH}_{x}$ phase for smaller cubes[17].

To investigate the pathways of reversible absorption in the model, free energy
calculations are carried out as a function of the concentration of
$\beta$-$\mathrm{PdH}_{x}$ ($c_{\beta}$) at coexistence (Extended Data Fig.
2a) and away from it (Extended Data Fig. 2b), using umbrella sampling (see
Supplementary Note 3.1). Representative configurations from four points along
the absorption pathway (Extended Data Fig. 2c,d and Supplementary Video 2)
recapitulate the observed experimental pathway of nucleation at one corner
followed by propagation of the interface along a $\langle 100\rangle$
direction. Extended Data Fig. 2e shows that the corner nucleated configuration
contains an interface, approximately in a $\langle 110\rangle$ direction, that
spans the length of the cube. This interface then moves in either of the two
equivalent $\langle 100\rangle$ directions, in order to join with an opposite
face of the cube and form a stable {100} interface that spans the cube.
Extended Data Fig. 2b depicts the shift in the free energy curve Extended Data
Fig. 2a in response to a chemical potential difference, $\Delta\mu$, between
the $\mathrm{H}_{2}$ molecules in solution and the H atoms in the palladium
lattice, which would act as a chemical driving force for absorption. The
curves show that the single-corner configuration represents a critical nucleus
when the chemical driving force is significantly smaller than the thermal
energy. While free energy curves are only presented for an approximately 6 nm
cube due to computational limitations, we expect that the reversible
absorption pathway would be recapitulated in larger cubes, as is observed in
the liquid phase TEM (see Supplementary Note 3.2). Additionally, the model
includes neither effects from diffusion of $\mathrm{H}_{2}$ in solution nor
defects in the nanocube. This indicates that the corner nucleation and [100]
propagation pathway can be understood purely in terms of the energetics of
elastic interactions (see Supplementary Note 3.2). The free energy curves have
only a single barrier, meaning that the stagnation in the reaction observed
experimentally after the corner nucleation (see Supplementary Note 1 and
Extended Data Fig. 1c) is likely due to either a depletion of $\mathrm{H}_{2}$
concentration near the cube, or the need to overcome an energetic barrier
caused by a defect.

##  Interfacial propagation simulations

###  Construction of density field and interface

A configuration in the elastic Ising model is represented by a discrete spin-
field, where each site has a value of $\pm 1$. However, in the physical
system, hydrogen density should vary continuously over space. To produce the
density maps of Extended Data Fig. 2c-e and Fig. 4a,b, we used the procedure
of Willard and Chandler[51], convoluting the discrete density field with
Gaussian distributions, such that the density for a given site is equal to,

$\displaystyle\tilde{c}_{\beta}(\mathbf{R})=\frac{1}{\mathcal{N}}\sum_{\mathbf{R^{\prime}}}\sigma_{\mathbf{R^{\prime}}}e^{\frac{-|\mathbf{R^{\prime}}-\mathbf{R}|^{2}}{2a^{2}}}$
(6)

In Supplementary Equation 6, $a$, is the lattice constant of the model. In
practice, the sum was cut off at a distance of five unit cells, and
$\mathcal{N}$ is the normalization of the discretized Gaussian over this
range. Periodic boundary conditions were used in the $x$ and $y$ directions,
to compute the density field for the rod geometries. For the cube geometries
of Extended Data Fig. 2c-e, and for the $z$-direction of the rod geometries,
hard-wall boundaries are used. This means that, in Supplementary Equation 6,
$\sigma_{\mathbf{R}}=0$ if $\mathbf{R}$ is outside the nanocrystal.

This convoluted density field was also used to compute the field, $h(x,y)$,
that defines the interface in the simulations of interfacial propagation (Fig.
4). For each pair of $x$ and $y$ values defining a column of sites along the
$z$-axis, the $z$-value, $z^{*}$, was identified such that
$z^{*}(x,y)=\mathrm{min}\\{z\hskip 2.84526pt|\hskip
2.84526pt\tilde{c}_{\beta}(x,y,z)<0.5\\}$. Then, via linear interpolation
between $\tilde{c}_{\beta}(x,y,z^{*}-a)$ and $\tilde{c}_{\beta}(x,y,z^{*})$,
the value of $h(x,y)$ was determined such that
$\tilde{c}_{\beta}(x,y,h(x,y))=0.5$. We then define the position of the
interface at a time, $t$, in Supplementary Equation 7.

$\displaystyle\frac{\bar{h}(t)}{a}=\frac{1}{N_{xy}}\sum_{x,y}h(x,y,t)$ (7)

Here, $N_{xy}$ is the number of total unit cells in the $x$ and $y$
dimensions, or number of columns.

We also define the instantaneous width of the interface in Supplementary
Equation 8,

$\displaystyle(\delta
h)^{2}/a^{2}=\frac{1}{N_{xy}}\sum_{x,y}\bigg{(}h(x,y,t)-\bar{h}(t)\bigg{)}^{2}$
(8)

as the variance of the height field. And, lastly, we define the average
velocity of the interface over a given time period in Supplementary Equation
9.

$\displaystyle\bar{v}=\frac{\bar{h}(t_{b})-\bar{h}(t_{a})}{t_{b}-t_{a}}$ (9)

To compute the velocities in Fig. 4 and Extended Data Fig. 4, $t_{b}$ and
$t_{a}$ were taken to be the time when $\bar{h}$ reached $25a$ and $31a$,
respectively, for each independent trajectory.

###  Computational details of rod trajectories

To decouple the nucleation process from the interfacial dynamics, the
simulations of interfacial propagation of Fig. 4 were carried out with
constrained dynamics in a rod geometry of 16 × 16 × 50 unit cells. Propagation
of $\beta$-$\mathrm{PdH}_{x}$ was directed in the positive $z$-direction by
holding the spins at the $z/a=1$ layer of the rod fixed at $\pm 1$ and the
spins at the $z/a=50$ layer of the rod fixed as $-1$. Due to computational
constraints, the size of the rod in the $x$ and $y$ dimensions is
significantly smaller than the nanocube used in the experiment. However, with
our constrained dynamics, the rod can be seen as the core of a larger nanocube
(Extended Data Fig. 4a). We allow the exchange of spins between any two
nearest neighbors, and we allow any single spin to be flipped if it is a site
at an interface, which is defined as any site that has both spin up and spin
down nearest neighbors. This is representative of the dynamics in the larger
nanocube in which hydrogen is absorbed at the surface of the nanocube and then
preferentially diffuses along the interface, towards the core of the cube. It
is consistent with the well-separated phase domains observed in the liquid
phase TEM experiments.

To verify that the observed behavior of interfacial propagation remains
qualitatively consistent across different sizes of the rod geometry, we ran
simulations with rods of sizes 12, 14, and 16 unit cells in the $x$-$y$
dimensions, keeping the length in the $z$-dimension fixed at 50 unit cells. At
low chemical driving force, $\Delta\mu<0.5\text{ }k_{\mathrm{B}}T$, there is
no statistically significant dependence of the average interfacial width or
average velocity on the size in the $x$-$y$ dimensions, for the sizes we have
tested (Extended Data Fig. 4b-d). Under a large chemical driving force,
$\Delta\mu>0.5\text{ }k_{\mathrm{B}}T$, the interfacial width does acquire a
non-trivial dependence on the size of the $x$-$y$ dimensions. Under these
conditions, the interfacial width diverges over time, and the rate of this
divergence appears to depend on the width of the rod, as would be expected
from a general, theoretical description of kinetic roughening[40]. In all such
cases, there is a clear distinction between a smooth regime, in which the
interfacial width is constant over time, and a kinetically rough regime, in
which the interfacial width increases over time (Extended Data Fig. 4d).
Additionally, the relative increase in the speed of propagation under a large
chemical driving force, $\Delta\mu>0.5\text{ }k_{\mathrm{B}}T$, decreases
slightly with the increasing width of the rod. However, extrapolating the
trend from 12 to 14 to 16 unit cells in rod-width (Extended Data Fig. 4c), it
appears that the speed will converge as the size increases and that the
relative increase in propagation speed will remain of the same order of
magnitude as the difference in propagation speed between absorption and
desorption observed in the liquid phase TEM. Therefore, the qualitative
conclusion that kinetic roughening contributes to the expanding of the
interface observed during faster propagation in the liquid phase TEM
experiments remains valid.

To compute average quantities like the velocity and the interfacial width,
data was taken from at least 115 independent trajectories under each value of
$\Delta\mu$. Every trajectory starts from the same initial condition, which is
an entirely spin down ($-1$) rod, except for the fixed layer of spin up ($+1$)
at one end. For the plots in Fig. 4d,e and Extended Data Fig. 4b,c, an average
is first taken for each individual trajectory over the range of time from
which the trajectory first reaches $\bar{h}/a=25$ to the time when it reaches
$\bar{h}/a=31$. This range was chosen such that it was large enough to
decrease the sampling error, but small enough so that the interface would not
feel effects from the boundaries at the ends of the rods. The time-averages
from each of the trajectories are themselves, then further averaged across
multiple independent trajectories. Error bars show two standard errors in the
mean from the independent trajectories. To compute the time series plots in
Fig. 4c and Extended Data Fig. 4d, a histogram in time is made for each
trajectory, and the average value of the observable is computed in each time
window. The elapsed time, $\Delta t$, is recorded from the time the interface
reaches an average height of $\bar{h}/a=10$. The histograms for each
trajectory are then averaged to produce the solid lines. The shaded region of
Extended Data Fig. 4d represents two standard errors.

Supplementary Fig. 10. Interface position as a function of time during the
hydrogen absorption and desorption reactions. a,b, Plots illustrating the
interface position along the [100] direction with reaction time during the
absorption (a) and desorption (b) reactions. As illustrated in the inset
images, four and three different interface regions (‘Regs’) were analyzed for
the absorption and desorption cases, respectively. Please note that the speed
of interface movement notably differs at the different Regs for the absorption
case. Supplementary Fig. 11. Comprehensive analysis of phase propagation
speeds during the absorption and desorption. a,b, Graphs depicting the average
speed of interface movement ($i.e.$, phase propagation) along the [100]
direction as a function of reaction time observed during the absorption (a)
and desorption (b) processes. The propagation speed as a function of time was
calculated using a linear fit applied to three consecutive data points from
Supplementary Fig. 10. Subsequently, these calculated speeds were averaged
across all regions to provide a comprehensive view of the phase propagation
dynamics for both absorption and desorption. Error bars display standard
deviations. Supplementary Fig. 12. Co-existence of two $\mathrm{PdH}_{x}$
phases in the Pd nanocubes during the hydrogen absorption and desorption. a,b,
HRTEM images of the nanocubes in the middle of the hydrogen absorption (a) and
desorption (b) reactions. These images are from Supplementary Videos 3 and 4.
The corresponding FFT patterns from the selected regions (colored rectangles)
within the nanocubes are also displayed along with the assigned $d$-spacing
values, demonstrating the co-existence of both $\alpha$\- and
$\beta$-$\mathrm{PdH}_{x}$ phases. Supplementary Fig. 13. Complete
transformation into single-crystalline $\alpha$\- or
$\beta$-$\mathrm{PdH}_{x}$ phases after the absorption and desorption
reactions. a,b, HRTEM images of the nanocubes after the hydrogen absorption
(a) and desorption (b). These images are from Supplementary Videos 3 and 4.
The corresponding FFT patterns from the entire nanocubes are presented with
the respective $d$-spacing values, showing the complete phase transformations
into $\alpha$\- or $\beta$-$\mathrm{PdH}_{x}$ phases. The nanocubes are
single-crystalline after the hydrogen absorption and desorption reactions
Supplementary Fig. 14. Swelling of the palladium nanocube due to hydrogen
absorption. a,b,c, This figure compares the 2D projected area of the palladium
nanocube at the start of the recorded absorption process (a) with its area
after complete phase transformation (b). The transformation results in a 4.7%
increase in the projected area (c).  Supplementary Fig. 15. Shrinking of the
palladium nanocube due to hydrogen desorption. a,b,c, This figure illustrates
the change in the 2D projected area of the palladium nanocube at the onset of
the imaged desorption process (a) and after complete phase transformation (b).
Overall, there is a 1.9% reduction in the nanocube’s projected area. However,
since only 50% of the nanocube is $\beta$-$\mathrm{PdH}_{x}$ at the start of
the HRTEM imaging (a), a more precise comparison involves only the
$\beta$-$\mathrm{PdH}_{x}$ portion (marked by white rectangles in (a) and (b),
and a black rectangle in (c)), showing a 5.4% decrease in this specific area
(c).  Supplementary Fig. 16. Local variations in $d$-spacing values in terms
of various Bragg reflections from the nanocube during the absorption. The
$d$-spacing color maps from the (200), (220), and (020) reflections
(represented in the FFT pattern) of the Pd nanocube during the absorption
reaction. Please note that the (200) color map is the same as Fig. 2b.
Supplementary Fig. 17. GPA $\mathrm{E}_{xx}$ strain map overlaid on the HRTEM
image. GPA analysis performed on the [010] ($\mathrm{E}_{xx}$) direction
related to the $\mathrm{E}_{yy}$ strain map in Fig. 2d. No significant strain
variations at the interphase boundary are visible in the [010] direction.
Supplementary Fig. 18. Local variations in $d$-spacing values in terms of
various Bragg reflections from the nanocube during the desorption. The
$d$-spacing color maps from the (200), (220), and (020) reflections
(represented in the FFT pattern) of the Pd nanocube during the desorption
reaction. Please note that the (200) color map is the same as Fig. 3b.
Supplementary Fig. 19 GPA $\mathrm{E}_{xx}$ strain map overlaid on the HRTEM
image. GPA analysis at the interface, which corresponds to the
$\mathrm{E}_{yy}$ strain map in Fig. 3c, demonstrates a relatively smooth
strain variation in the [010] ($\mathrm{E}_{xx}$) direction.

## References

## References

  * Zhang  _et al._ [2022a] H. Zhang, N. Pryds, D.-S. Park, N. Gauquelin, S. Santucci, D. V. Christensen, D. Jannis, D. Chezganov, D. A. Rata, A. R. Insinga, I. E. Castelli, J. Verbeeck, I. Lubomirsky, P. Muralt, D. Damjanovic, and V. Esposito, “Atomically engineered interfaces yield extraordinary electrostriction,” Nature 609 (2022a).
  * Sun  _et al._ [2022] X. Sun, D. Wu, L. Zou, S. D. House, X. Chen, M. Li, D. N. Zakharov, J. C. Yang, and G. Zhou, “Dislocation-induced stop-and-go kinetics of interfacial transformations,” Nature 607 (2022).
  * Zheng  _et al._ [2004] H. Zheng, J. Wang, S. E. Lofland, Z. Ma, L. Mohaddes-Ardabili, T. Zhao, L. Salamanca-Riba, S. R. Shinde, S. B. Ogale, F. Bai, D. Viehland, Y. Jia, D. G. Schlom, M. Wuttig, A. Roytburd, and R. Ramesh, “Multiferroic BaTiO3-CoFe2O4 Nanostructures,” Science 303 (2004).
  * Baldi  _et al._ [2014] A. Baldi, T. C. Narayan, A. L. Koh, and J. A. Dionne, “In situ detection of hydrogen-induced phase transitions in individual palladium nanocrystals,” Nature Mater 13 (2014).
  * Zou  _et al._ [2018] L. Zou, C. Yang, Y. Lei, D. Zakharov, J. M. K. Wiezorek, D. Su, Q. Yin, J. Li, Z. Liu, E. A. Stach, J. C. Yang, L. Qi, G. Wang, and G. Zhou, “Dislocation nucleation facilitated by atomic segregation,” Nature Mater 17, 56–63 (2018).
  * Syrenova  _et al._ [2015] S. Syrenova, C. Wadell, F. A. A. Nugroho, T. A. Gschneidtner, Y. A. Diaz Fernandez, G. Nalin, D. Świtlik, F. Westerlund, T. J. Antosiewicz, V. P. Zhdanov, K. Moth-Poulsen, and C. Langhammer, “Hydride formation thermodynamics and hysteresis in individual Pd nanocrystals with different size and shape,” Nature Mater 14, 1236–1244 (2015).
  * Steele  _et al._ [2019] J. A. Steele, H. Jin, I. Dovgaliuk, R. F. Berger, T. Braeckevelt, H. Yuan, C. Martin, E. Solano, K. Lejaeghere, S. M. J. Rogge, C. NOOPnotebaert, W. Vandezande, K. P. F. Janssen, B. Goderis, E. Debroye, Y.-K. Wang, Y. Dong, D. Ma, M. Saidaminov, H. Tan, Z. Lu, V. Dyadkin, D. Chernyshov, V. Van Speybroeck, E. H. Sargent, J. Hofkens, and M. B. J. Roeffaers, “Thermal unequilibrium of strained black CsPbI3 thin films,” Science 365, 679–684 (2019).
  * Chen  _et al._ [2020] Y. Chen, Y. Lei, Y. Li, Y. Yu, J. Cai, M.-H. Chiu, R. Rao, Y. Gu, C. Wang, W. Choi, H. Hu, C. Wang, Y. Li, J. Song, J. Zhang, B. Qi, M. Lin, Z. Zhang, A. E. Islam, B. Maruyama, S. Dayeh, L.-J. Li, K. Yang, Y.-H. Lo, and S. Xu, “Strain engineering and epitaxial stabilization of halide perovskites,” Nature 577, 209–215 (2020).
  * Liu  _et al._ [2024] G. Liu, A. J. Shih, H. Deng, K. Ojha, X. Chen, M. Luo, I. T. McCrum, M. T. M. Koper, J. Greeley, and Z. Zeng, “Site-specific reactivity of stepped Pt surfaces driven by stress release,” Nature 626, 1005–1010 (2024).
  * Bai, Misra, and Fan [2022] Z. Bai, A. Misra, and Y. Fan, “Universal Trend in the dynamic relaxations of tilted metastable grain boundaries during ultrafast thermal cycle,” Materials Research Letters 10, 343–351 (2022).
  * Oh  _et al._ [2020] M. H. Oh, M. G. Cho, D. Y. Chung, I. Park, Y. P. Kwon, C. Ophus, D. Kim, M. G. Kim, B. Jeong, X. W. Gu, J. Jo, J. M. Yoo, J. Hong, S. McMains, K. Kang, Y.-E. Sung, A. P. Alivisatos, and T. Hyeon, “Design and synthesis of multigrain nanocrystals via geometric misfit strain,” Nature 577, 359–363 (2020).
  * Xu  _et al._ [2019] Y. Xu, J.-H. Park, Z. Yao, C. Wolverton, M. Razeghi, J. Wu, and V. P. Dravid, “Strain-Induced Metastable Phase Stabilization in Ga2O3 Thin Films | ACS Applied Materials & Interfaces,”  (2019).
  * Kacher, Eftink, and Robertson [2019] J. Kacher, B. P. Eftink, and I. M. Robertson, “In situ transmission electron microscopy investigation of dislocation interactions,” in  _Handbook of Mechanics of Materials_ (Unknown Publisher, 2019) pp. 131–166.
  * Wyckoff [1924] R. W. G. Wyckoff,  _The Structure of Crystals_ (Chemical Catalog Company, Incorporated, 1924).
  * Wicke and Brodowsky [1978] E. Wicke and H. Brodowsky, “Hydrogen in palladium and palladium alloys,” , 73–155 (1978).
  * Ingham  _et al._ [2008] B. Ingham, M. F. Toney, S. C. Hendy, T. Cox, D. D. Fong, J. A. Eastman, P. H. Fuoss, K. J. Stevens, A. Lassesson, S. A. Brown, and M. P. Ryan, “Particle size effect of hydrogen-induced lattice expansion of palladium nanoclusters,” Phys. Rev. B 78, 245408 (2008).
  * Yamauchi  _et al._ [2008] M. Yamauchi, R. Ikeda, H. Kitagawa, and M. Takata, “Nanosize Effects on Hydrogen Storage in Palladium,” J. Phys. Chem. C 112, 3294–3299 (2008).
  * Liu  _et al._ [2011] N. Liu, M. L. Tang, M. Hentschel, H. Giessen, and A. P. Alivisatos, “Nanoantenna-enhanced gas sensing in a single tailored nanofocus,” Nature Mater 10, 631–636 (2011).
  * Bardhan  _et al._ [2013] R. Bardhan, L. O. Hedges, C. L. Pint, A. Javey, S. Whitelam, and J. J. Urban, “Uncovering the intrinsic size dependence of hydriding phase transformations in nanocrystals,” Nature Mater 12, 905–912 (2013).
  * Li  _et al._ [2014] G. Li, H. Kobayashi, J. M. Taylor, R. Ikeda, Y. Kubota, K. Kato, M. Takata, T. Yamamoto, S. Toh, S. Matsumura, and H. Kitagawa, “Hydrogen storage in Pd nanocrystals covered with a metal–organic framework,” Nature Mater 13, 802–806 (2014).
  * Griessen, Strohfeldt, and Giessen [2016] R. Griessen, N. Strohfeldt, and H. Giessen, “Thermodynamics of the hybrid interaction of hydrogen with palladium nanoparticles,” Nature Mater 15, 311–317 (2016).
  * Narayan  _et al._ [2017] T. C. Narayan, F. Hayee, A. Baldi, A. Leen Koh, R. Sinclair, and J. A. Dionne, “Direct visualization of hydrogen absorption dynamics in individual palladium nanoparticles,” Nat Commun 8, 14020 (2017).
  * Ulvestad  _et al._ [2017] A. Ulvestad, M. J. Welland, W. Cha, Y. Liu, J. W. Kim, R. Harder, E. Maxey, J. N. Clark, M. J. Highland, H. You, P. Zapol, S. O. Hruszkewycz, and G. B. Stephenson, “Three-dimensional imaging of dislocation dynamics during the hydriding phase transformation,” Nature Mater 16, 565–571 (2017).
  * Schneemann  _et al._ [2018] A. Schneemann, J. L. White, S. Kang, S. Jeong, L. F. Wan, E. S. Cho, T. W. Heo, D. Prendergast, J. L. Urban, B. C. Wood, M. D. Allendorf, and V. Stavila, “Nanostructured Metal Hydrides for Hydrogen Storage | Chemical Reviews,”  (2018).
  * Sytwu  _et al._ [2018] K. Sytwu, F. Hayee, T. C. Narayan, A. L. Koh, R. Sinclair, and J. A. Dionne, “Visualizing Facet-Dependent Hydrogenation Dynamics in Individual Palladium Nanoparticles,” Nano Lett. 18, 5357–5363 (2018).
  * Johnson  _et al._ [2019] N. J. J. Johnson, B. Lam, B. P. MacLeod, R. S. Sherbo, M. Moreno-Gonzalez, D. K. Fork, and C. P. Berlinguette, “Facets and vertices regulate hydrogen uptake and release in palladium nanocrystals,” Nat. Mater. 18, 454–458 (2019).
  * Koo  _et al._ [2024] K. Koo, Z. Li, Y. Liu, S. M. Ribet, X. Fu, Y. Jia, X. Chen, G. Shekhawat, P. J. M. Smeets, R. dos Reis, J. Park, J. M. Yuk, X. Hu, and V. P. Dravid, “Ultrathin silicon nitride microchip for in situ/operando microscopy with high spatial resolution and spectral visibility,” Science Advances 10, eadj6417 (2024).
  * Yin  _et al._ [2019] Z.-W. Yin, S. B. Betzler, T. Sheng, Q. Zhang, X. Peng, J. Shangguan, K. C. Bustillo, J.-T. Li, S.-G. Sun, and H. Zheng, “Visualization of facet-dependent pseudo-photocatalytic behavior of TiO2 nanorods for water splitting using In situ liquid cell TEM,” Nano Energy 62, 507–512 (2019).
  * Zhang  _et al._ [2022b] Q. Zhang, X. Peng, Y. Nie, Q. Zheng, J. Shangguan, C. Zhu, K. C. Bustillo, P. Ercius, L. Wang, D. T. Limmer, and H. Zheng, “Defect-mediated ripening of core-shell nanostructures,” Nat Commun 13, 2211 (2022b).
  * Hart [1964] E. J. Hart, “The Hydrated Electron,” Science 146, 19–25 (1964).
  * Grogan  _et al._ [2014] J. M. Grogan, N. M. Schneider, F. M. Ross, and H. H. Bau, “Bubble and Pattern Formation in Liquid Induced by an Electron Beam,” Nano Lett. 14, 359–364 (2014).
  * Schneider  _et al._ [2014] N. M. Schneider, M. M. Norton, B. J. Mendel, J. M. Grogan, F. M. Ross, and H. H. Bau, “Electron–Water Interactions and Implications for Liquid Cell Electron Microscopy | The Journal of Physical Chemistry C,”  (2014).
  * Fratzl, Penrose, and Lebowitz [1999] P. Fratzl, O. Penrose, and J. L. Lebowitz, “Modeling of Phase Separation in Alloys with Coherent Elastic Misfit,” Journal of Statistical Physics 95, 1429–1503 (1999).
  * Frechette, Dellago, and Geissler [2020] L. B. Frechette, C. Dellago, and P. L. Geissler, “Origin of mean-field behavior in an elastic Ising model,” Phys. Rev. B 102, 024102 (2020).
  * Frechette, Dellago, and Geissler [2021] L. B. Frechette, C. Dellago, and P. L. Geissler, “Elastic forces drive nonequilibrium pattern formation in a model of nanocrystal ion exchange | PNAS,”  (2021).
  * Hÿtch, Snoeck, and Kilaas [1998] M. J. Hÿtch, E. Snoeck, and R. Kilaas, “Quantitative measurement of displacement and strain fields from HREM micrographs,” Ultramicroscopy 74, 131–146 (1998).
  * Rouvière and Sarigiannidou [2005] J. L. Rouvière and E. Sarigiannidou, “Theoretical discussions on the geometrical phase analysis,” Ultramicroscopy 106, 1–17 (2005).
  * Zhu  _et al._ [2013] Y. Zhu, C. Ophus, J. Ciston, and H. Wang, “Interface lattice displacement measurement to 1 pm by geometric phase analysis on aberration-corrected HAADF STEM images,” Acta Materialia 61, 5646–5663 (2013).
  * Nishino  _et al._ [2013] M. Nishino, T. Nakada, C. Enachescu, K. Boukheddaden, and S. Miyashita, “Crossover of the roughness exponent for interface growth in systems with long-range interactions due to lattice distortion,”  (2013).
  * Tong and Williams [1994] W. M. Tong and R. S. Williams, “Kinetics of Surface Growth: Phenomenology, Scaling, and Mechanisms of Smoothening and Roughening,” Annual Review of Physical Chemistry 45, 401–438 (1994).
  * Mukherji and Bhattacharjee [1997] S. Mukherji and S. M. Bhattacharjee, “Nonlocality in Kinetic Roughening,” Phys. Rev. Lett. 79, 2502–2505 (1997).
  * Epstein and Plesset [1950] P. S. Epstein and M. S. Plesset, “On the Stability of Gas Bubbles in Liquid‐Gas Solutions,” The Journal of Chemical Physics 18, 1505–1509 (1950).
  * Kim  _et al._ [2021] Q. Kim, D. Shin, J. Park, D. A. Weitz, and W. Jhe, “Initial growth dynamics of 10 nm nanobubbles in the graphene liquid cell,” Appl Nanosci 11, 1–7 (2021).
  * Jany, Janas, and Krok [2020] B. R. Jany, A. Janas, and F. Krok, “Automatic microscopic image analysis by moving window local Fourier Transform and Machine Learning,” Micron 130, 102800 (2020).
  * Hsu and Leisure [1979] D. K. Hsu and R. G. Leisure, “Elastic constants of palladium and $\ensuremath{\beta}$-phase palladium hydride between 4 and 300 K,” Phys. Rev. B 20, 1339–1344 (1979).
  * Frechette, Dellago, and Geissler [2019] L. B. Frechette, C. Dellago, and P. L. Geissler, “Consequences of Lattice Mismatch for Phase Equilibrium in Heterostructured Solids,” Phys. Rev. Lett. 123, 135701 (2019).
  * Gillespie [2007] D. T. Gillespie, “Stochastic Simulation of Chemical Kinetics,” Annual Review of Physical Chemistry 58, 35–55 (2007).
  * Frechette [2020] L. B. Frechette,  _Chemical Transformations of Nanocrystals: Theory and Molecular Simulation_ , Ph.D., University of California, Berkeley, United States – California (2020).
  * Grossﬁeld [2011] A. Grossﬁeld, “An implementation of WHAM: the Weighted Histogram Analysis Method Version 2.0.10,” (2011).
  * Flanagan and Oates [1991] T. B. Flanagan and W. A. Oates, “The Palladium-Hydrogen System,” Annual Review of Materials Research 21, 269–304 (1991).
  * Willard and Chandler [2010] A. P. Willard and D. Chandler, “Instantaneous Liquid Interfaces,” J. Phys. Chem. B 114, 1954–1958 (2010).