# Unveiling the Electronic, Transport and Migration properties
of the Te-defect lattice in DyTe1.8

Jinwoong Kim<EMAIL_ADDRESS>Nicholas Kioussis<EMAIL_ADDRESS>Department of Physics and Astronomy, California State University, Northridge,
California 91330, USA

###### Abstract

The rare-earth ditellurides are known to form two-dimensional square lattice
where the strong Fermi surface nesting leads to structural modulation. In
contrast to charge density waves, the supercell modulation is accompanied by
the formation of the periodic Te vacancy network, where the Te deficiency
affects the nesting vector (i.e. the supercell size) via tuning the chemical
potential. In this work, first principles electronic structure calculations
for the $\sqrt{5}\times\sqrt{5}$ supercell, that commonly appears in this
family of tellurides, unveil interesting electronic, transport, and migration
properties of the Te defect lattice in DyTe1.8. The reconstruction of the Te-
deficient square lattice, consisting of a single Te-dimer and a pair Te-
trimers per unit cell, gives rise to an out-of-plane polarization, whose
direction depends on the position of the dimer. This results in various close-
in-energy parallel and antiparallel polarization configurations of successive
Te layers depending on the dimer positions. We predict that the orientation of
the Te dimers, and hence the corresponding structural motifs, can be
reversibly switched between two in-plane perpendicular directions under
tensile epitaxial strain via a piezoelectric substrate, resulting in a
colossal conductivity switching. Furthermore, the Te-dimer orientations result
in asymmetric Fermi surface which can be confirmed by quantum oscillations
measurements. Finally, we present numerical results for the migration paths
and energy landscape through various divacancy configurations in the presence
or absence of epitaxial strain.

##  I Introduction

The rare-earth ditellurides $R$Te2 ($R$ represents a lanthanide element) are a
family of layered materials renowned for a wide range of intriguing properties
at low-temperature. These include charge density wave (CDW)[1, 2]
magnetism,[3] superconductivity under high pressure,[4] as well as the
simultaneous presence of CDW, superconductivity, and antiferromagnetism.[1]
The unit cell of the crystal structure comprises of double-corrugated $R$-Te
quasi-ionic layers separated by single planar square Te sheets stacked along
the $[001]$ direction, as illustrated in Fig. 1(a). From the electronic
structure point of view, each ($R$-Te)1+ pair on the corrugated layer donates
one electron to the Te atoms (Tes) in the square plane.[5, 6]

In contrast to the rare earth tritellurides, $R$Te3, which form stoichiometric
compounds, the $R$Te2-δ ($0.0\lesssim\delta\lesssim 0.2$) compounds host Tes
vacancies on the Tes square layer, giving rise to a wide range of structural
motifs which depends sensitively on the Te deficiency $\delta$ and the rare
earth cation.[7, 8] The Tes vacancies effectively contribute extra charge to
the layer which along with the isolated Te2- anions form differently ordered
patterns within the Te sheets. This results in superstructure modulations
which can be commensurate with the parent structure, such as in CeSe1.9[9] and
Gd8Se15[10] structures, and incommensurate modulated structures such as
DySe1.84[11] and $R$Se1.84 ($R$=La-Nd, Sm)[12].

Recently, we have reported[7] experimental and first principles electronic
structure studies of the supercell formation in epitaxial DyTe2-δ thin films,
where for $\delta$=0.2 the modulation is a ($\sqrt{5}\times\sqrt{5}$) R26.6∘
$\times$ 2 superlattice. The superlattice emerges due to a periodic Te-defect
lattice, which also acts to open a gap in the electronic spectrum and induces
semiconducting transport behavior. First principles calculations point towards
nesting conditions of the Fermi surface at a $\bm{q}$-vector that corresponds
to the $\sqrt{5}\times\sqrt{5}$ condition, suggesting that the formation of
the defect lattice results from a similar driving force to the conventional
picture of CDW formation, where sections of the Fermi surface are gapped out
by the formation of supercells with periodicity corresponding to the nesting
condition[13]. Our total energy calculations[7] of the
$\sqrt{5}\times\sqrt{5}$ supercell with Tes mono- and divacancies in different
configurations for $\delta$ = 0.1 and 0.2, showed that the second nearest-
neighbor (SNN) vacancy configuration (denoted as the A-C configuration in Fig.
5c of Ref. [7]) has the lowest formation energy among the various considered
structures, in agreement with previous X-ray diffraction measurements in
SmTe1.8.[14] The SNN di-vacancy configuration shown in Fig. 2(c), where the
di-vacancies nucleate on sites 0 and 4 $\\{0,4\\}$ in Fig. 2(a), consists of
an ordered lattice of Tes dimers and trimers.

Furthermore, first principles calculations of the
($\sqrt{5}\times\sqrt{5}$)$\times$ 2 modulation along the $c$ axis
corresponding to the nucleation of a pair of A-C (0,4) divacancies on two Te
square layers, reveals that the SNN A-C divacancy appears at a laterally
displaced position $(0.5,0.5)$, relative to the divacancy on the adjacent Te
square net layer. [7]

The objective of this work is to unveil the emergence of novel interesting
properties of the periodic Te-defect lattice associated with the onset of
ferroelectric polarization which depends sensitively on the position of the
Te-dimer. In Sec. II we outline the methodology employed. In Sec. III we
present results of (a) the position and orientation of the various Te-dimer
configurations, (b) their polarization properties, (c) the electronic
structure, (d) the Fermi surface anisotropy, (e) and (f) the transport and
migration properties in the absence or presence of epitaxial strain.
Conclusions are summarized in Sev. IV.

##  II Methodology

The density functional theory calculations are performed using the Vienna ab
intio simulation package (VASP)[15, 16] with the projector augmented wave
method.[17, 18] The Dy$–4f$ states are treated as core with valence
configurations $5p^{6}5d^{1}6s^{2}$ and $5s^{2}5p^{4}$ for Te. All structures
are optimized using the PBEsol exchange correlation functional,[19] and the
Bloch states are calculated and Wannierized[20, 21] with the modified Becke-
Johnson potential[22, 23] which provides accurate band gaps, effective masses,
and frontier-band ordering. The momentum space is sampled at a $80\times
80\times 80$ k-point mesh with the Wannier interpolation scheme. The
conductivity is calculated using the semi-classical Boltzmann transport theory
assuming a constant relaxation time. The tetrahedron method is employed to
obtain the Fermi surface and its cross section area. In the polarization
calculation, the topmost 51 valence bands are Wannierized for each spin
channel where three and two $p$ orbitals are initially projected on the Tec
and Tes atoms, respectively, and five $s$ orbitals are projected at the center
of the Te-Te bonds in the square lattice. The integer ionic charges $Z_{n}$ in
Eq. 3 are $+3$ for Dy and $+4$ for Te atoms, respectively.

##  III RESULTS AND DISCUSSION

###  III.1 Crystal structure

Figure 1: (a) Unit cell of pristine DyTe2 with two Te sites (Tes) on the
square net layer. (b) Calculated nesting function of DyTe2. The back ground
color map illustrates Gaussian-smeared nesting function for an electron doping
level corresponding to Te deficiency, $\delta=0.2$. The cross point indicates
$\sqrt{5}\times\sqrt{5}$ wave vector. The two local nesting function peaks are
marked with dots upon the doping level, that moves away from the origin with
increasing electron doping.

The crystal structure of pristine DyTe2 is shown in Fig. 1(a), where the
buckled DyTe slabs are separated by the planar Te square lattice that tends to
undergo a supercell modulation whose periodicity interestingly depends on the
Te vacancy concentration.[7, 8] The nesting function, $N(\bm{q})$, of the
Fermi surface for an electron doping corresponding to the Te deficiency
$\delta=0.2$, reported in our recent calculations [Fig. S7(h) in Supporting
information of Ref. [7]], is replotted in Fig. 1(b). It is zoomed-in near the
$\bm{q}$ point corresponding to the $\sqrt{5}\times\sqrt{5}$ modulation,
marked as a cross point at $(0.2,0.4)$. Two strong peaks near the cross point,
illustrated as dark orange colors, are traced as a function of doping which
are displayed with the blue and green dots, where their separation increases
with electron doping. The green peak is found to pass the cross point at
higher Te deficient level ($\delta>0.2$) that does not necessarily imply the
need for additional Te vacancies, since the nesting function calculations
involve several assumptions such as the rigid band approximation, constant
matrix elements, and the exclusion of strain effect on the chemical potential.
Our calculations demonstrate that the nesting function peaks associated with
the formation of various supercell modulations change with chemical doping.
Especially, the nesting function peaks (denoted by blue and green dots)
approach the $\sqrt{5}\times\sqrt{5}$ modulation $\bm{q}$ vectors with
increasing Te deficiency in a good agreement with experiments.

Figure 2: Structural modulation of DyTe2. (a) Projection of the Te-square
network on the (001) plane, where the numerical indices of the
$\sqrt{5}\times\sqrt{5}$-modulated cell denote the plausible positions of the
Tes divacancies. The $\bm{a}^{p}_{i}$ are the primitive lattice vectors and
$\bm{a}^{s\pm}_{i}$ are two sets of supercell lattice vectors with $\pm$
chirality. (b) Side view of $\sqrt{5}\times\sqrt{5}\times 2$ modulated
DyTe1.8. Arrows depict the polarization configuration of individual Te square
lattice layers in the ground state. (c-e) Reconstructed ordered Te-square
network of DyTe1.8 for several Te divacancy positions and orientations at the
(c) {0,4}, (d) {1,5}, and (e) {1,3} sites. Small squares in (c) and (d) denote
the pristine cell enclosing the Te-dimer. (f) Mixed phase of the reconstructed
Te-square network comprised of the structural motifs of (c) and (e),
respectively, where the large square denotes the $\sqrt{5}\times\sqrt{5}$
cell.

Consequently, in this work we focus on the $\sqrt{5}\times\sqrt{5}$ modulation
that commonly occurs for various rare-earth elements with two Tes vacancies
per modulated cell ($\delta=0.2$)[7]. The in-plane lattice vectors of the
supercell exhibiting chirality, shown in Fig. 2(a), are

$\left(\begin{aligned}
\bm{a}^{p}_{1}&=&\tfrac{a_{0}}{\sqrt{2}}(\hat{x}-\hat{y})\\\
\bm{a}^{p}_{2}&=&\tfrac{a_{0}}{\sqrt{2}}(\hat{x}+\hat{y})\\\
\end{aligned}\right.\qquad\left(\begin{aligned}
\bm{a}^{s\pm}_{1}&=&2\bm{a}^{p}_{1}\pm\bm{a}^{p}_{2}\\\
\bm{a}^{s\pm}_{2}&=&\mp\bm{a}^{p}_{1}+2\bm{a}^{p}_{2},\end{aligned}\right.$
(1)

where $\bm{a}^{p}_{i}$ are the primitive lattice vectors and
$\bm{a}^{s\pm}_{i}$ are two sets of supercell lattice vectors with $\pm$
chirality. In general, a synthesized sample may have both types of domains
leading to domain walls across which the chirality gets reversed. The two
supercell structures with different chirality and their physical properties
can be transformed to each other by a mirror operator,
$\mathcal{M}_{y}:(x,y,z)\rightarrow(x,-y,z)$. Unless otherwise stated
throughout the remainder of this manuscript we focus only the $+$ chirality
structure and ignore for convenience the superscript s+, namely,
$\bm{a}_{i}\equiv\bm{a}_{i}^{s+}$.

The numerical indices in Fig. 2(a) represent the ten Tes sites of the
$\sqrt{5}\times\sqrt{5}$ supercell where the divacancy may be placed. The SNN
divacancy (such as the $\\{0,4\\}$, denoted as the A-C configuration in Ref.
[7]) was found to have the lowest formation energy for the
$\sqrt{5}\times\sqrt{5}$ modulation. The reconstructed square lattice consists
of one Te-dimer and two Te-trimers per unit cell, resulting in several iso-
energetic and ordered structural motifs. Figures 2(c) and (d) show the first
and second structural motifs, associated with the $\\{0,4\\}$ and $\\{1,5\\}$
SNN divacancies, respectively, that look almost identical except for their
positions relative to the DyTe slab. The thick squares in each panel emphasize
the position of the Tes dimers whose centers lie atop of the Dy (Tec) atoms
for the first (second) structural motif. The two structural motifs transform
into each other through the mirror operation, $\mathcal{M}_{z}$, followed by a
proper in-plane translation indicating a finite vertical polarization (see
Fig. 3) whose sign depends on the position of the Te dimer, as will be
discussed in the next section. The divacancy also reduces the
$\mathcal{C}_{4}$ rotation symmetry of the pristine supercell to the
$\mathcal{C}_{2}$ symmetry, whose rotational axis can be placed at the center
of the Te dimer. The third structural motif associated with the $\\{1,3\\}$
divacancy can be transformed into the first motif of the $\\{0,4\\}$ divacancy
through a $\mathcal{C}_{4}$ rotation, as illustrated in Figs. 2(c) and (e).
The actual DyTe1.8 samples may consist of a mixed phase of the $\\{0,4\\}$\-
and $\\{1,3\\}$-based structural motifs as shown in Fig. 2(f). The mixed phase
can then preserve the tetragonal symmetry in macroscopic scale while the
microscopic local domains are under uniaxial in-plane strain depending on the
Te dimer orientation. This in turn raises an intriguing question whether an
in-plane strain can align the orientation of the Te dimers, which is discussed
in detail in Sec. III.5.

In general, CDWs with a commensurate cell modulation of $n\times m$ can
introduce $nm$ distinct types of domains, associated with the relative
translations of the modulated supercell. Therefore, there can be five
different types of domains for the $\sqrt{5}\times\sqrt{5}$ modulation of
DyTe1.8. It is noteworthy that the three structural motifs discussed above do
not transform into one another via lateral translation alone. The total number
of distinct domains amounts to $2\times 2\times 5=20$, where the first
(second) factor of 2 arises from the two possible positions (orientations) of
the Te dimer in addition to the factor of 5 due to the lateral translations.

###  III.2 Polarization

Figure 3: Polarization of the Te square network. (a) Center of Te-dimer is atop of the Dy atom corresponding to the {0,4} divacancy configuration. (b) Center of Te-dimer is atop of the Te atom corresponding to the {1,5} divacancy configuration. The interlayer distances of the corrugated Dy-Te layers sandwiching the Te-dimer sheet are different. The Te-trimers are omitted for clarity. (c) Wannier charge centers lying on the covalent bonds, denoted by the purple dots on the reconstructed Te-square lattice for the {0,4} divacancy configuration. The position of Wannier charge centers affects the polarization as well as the dipole moments.  Table 1:  Dipole moment and polarization of DyTe1.8 for the {0,4} divacancy along the supercell lattice vectors, $\bm{a}_{i}$, of the $\sqrt{5}\times\sqrt{5}\times 1$ cell.  | Dipole moment ($Q_{e}\bm{a}_{i}$) | Polarization (C/m2)
---|---|---
$\bm{a}_{1}$ | 0.500 | 0.095
$\bm{a}_{2}$ | 0.500 | 0.095
$\bm{a}_{3}$ | 0.893 | 0.160

The emergence of the ($\sqrt{5}\times\sqrt{5}$) superlattice induced by the
periodic Te-vacant sites renders the DyTe1.8 insulator.[7] The formation of
the divacancy breaks the inversion symmetry which in turn results in finite
polarization. Figures 3(a) and (b) clearly show the Dy-Te bonds of the Tes
dimers for the first and second structural motifs, which depend on their
lateral position. Namely, in the first motif, where the center of the Te dimer
is atop of the Dy atom, the dimer forms four (two) bonds with the Dy atoms on
the upper (lower) layer, resulting in its proximity to the upper Dy-Te
corrugated layer. Conversely, the second structural motif exhibits the
opposite behavior.

The polarization $\bm{p}$ and dipole moment $\bm{d}$ of an insulator can be
computed by considering the position of charged ions and Wannier charge
centers (WCC) of the occupied Bloch states below the insulating band gap.[21,
24] As described in the methodology section, the topmost 51 valence bands are
Wannierized which are well separated from both the semi-core states and the
conduction bands. The Wannier functions are found exponentially localized near
the initial position of projectors. Three and two Wannier functions are thus
attached to Tec and Tes atoms, respectively. The positions of five
interstitial Wannier functions, illustrated in Fig. 3(c), lie on the Tes-Tes
bonds. The Wannier functions that deviate from the atomic sites render the
system an obstructed atomic insulator.[25, 26, 27, 28] The polarization and
dipole moment are given by,

$\displaystyle\bm{d}$ $\displaystyle=$ $\displaystyle\bm{p}V$ (2)
$\displaystyle=$ $\displaystyle
Q_{e}\left(\sum_{n}^{\mathrm{ions}}Z_{n}\left(\bm{R}_{n}-\bm{r}_{0}\right)-\sum_{n}^{\mathrm{WCC}}D\left(\bm{r}_{n}-\bm{r}_{0}\right)\right)$
(3) $\displaystyle=$ $\displaystyle
Q_{e}\left(\sum_{n}^{\mathrm{ions}}Z_{n}\bm{R}_{n}-\sum_{n}^{\mathrm{WCC}}D\bm{r}_{n}-\bm{r}_{0}\left\\{N^{+}-N^{-}\right\\}\right)$
(4) $\displaystyle=$ $\displaystyle
Q_{e}\left(\sum_{n}^{\mathrm{ions}}Z_{n}\bm{R}_{n}-\sum_{n}^{\mathrm{WCC}}D\bm{r}_{n}\right),$
(5)

where $V$ is the unit cell volume, $Q_{e}$ is the electron charge, $Z_{n}$ and
$\bm{R}_{n}$ are the ionic charge number and atomic position of the $n$-th
atom, $\bm{r}_{n}$ is the position of the $n$th WCC, and $D=2$ is the spin
degeneracy of the Wannier functions. The origin $\bm{r}_{0}$ does not affect
the dipole moment as it is charge neutral, $N^{+}=N^{-}$, where the net ionic
charge, $N^{+}=\sum_{n}^{\mathrm{ions}}Z_{n}$ is equal to the number of WCCs,
$N^{-}=D\sum_{n}^{\mathrm{WCC}}$. Note that the WCC position, $\bm{r}_{n}$,
can be shifted by integer multiples of lattice vectors,
$\bm{r}_{n}\rightarrow\bm{r}_{n}+\sum_{i}m_{i}\bm{a}_{i}$ without affecting
the physical observables that is associated with the gauge freedom leading to
the modulo definition of the modern theory of polarization. [29, 30, 31] The
dipole moment is thus defined under modulo $DQ_{e}\bm{a}_{i}$ as long as the
time reversal symmetry is preserved on the surface as well as in the bulk.

The calculated dipole moment and polarization of the first structural motif of
the $\\{0,4\\}$ divacancy, are listed in Table 1. The $\mathcal{C}_{2}$
symmetry of the crystal structure confines the in-plane components of the
dipole moments to be quantized[32] either to 0 or $1/2$ in units of
$Q_{e}\bm{a}_{1,2}$. The calculated in-plane dipole moments values of, $1/2$,
indicate nontrivial Zak phase[32] of DyTe1.8 in the
$\sqrt{5}\times\sqrt{5}\times 1$ supercell that induces metallic (side)
surfaces due to the emergent in-gap states. We note that the physical samples,
however, exhibit additional $\times 2$ structural modulation along the
$\hat{z}$ direction[7] resulting in $\sqrt{5}\times\sqrt{5}\times 2$ supercell
which cancels out the in-plane dipole moment. The concept of quantized
polarization under crystal symmetry was first elucidated by Zak[32] and has
recently been integrated into the framework of topological insulators.[33, 30,
34, 35]

This additional modulation also controls the net dipole moment along the
$\hat{z}$ direction as each layer may have upward or downward dipole moment
depending on the position of the Tes dimers of each layer. In our previous
study,[7] the total energy of the $\sqrt{5}\times\sqrt{5}\times 2$ supercell
was calculated (see Table 2) for various pairs of divacancy configurations,
$\\{0,4\\}-\\{d,d+4\\}$, $d=\\{0,1,2,3,4,5\\}$ where the first divacancy
nucleates at sites $\\{0,4\\}$ on the first Tes layer while the second pair
nucleates on sites $\\{d,d+4\\}$ on the second Tes layer. The supercell
structure with $d=5$ is found to exhibit the lowest energy in which the
polarization direction of the upper layer is opposite to that of the lower
layer as illustrated in Fig. 2(b). The next stable structure with $d=2$ has
2.4 meV/Dy higher total energy and the polarization directions of the two
layers are parallel leading to finite net polarization, suggesting that an
out-of-plane external electric field, $E\gtrsim 15$ MV/m may stabilize the
polarized structure.

Table 2: Calculated total energy of the $\sqrt{5}\times\sqrt{5}\times 2$ supercell for the lateral position, $d$, of the divacancy on the second Te layer relative to the first. The relative polarization directions of the two layers are labeled as P and AP for parallel and antiparallel configurations, respectively.  | Displacement of 2nd layer, d
---|---
| 0 | 1 | 2 | 3 | 4 | 5
$\Delta E$ (meV/Dy) | 0.0 | -6.0 | -8.4 | -5.8 | -7.3 | -10.8
Polarization | P | AP | P | AP | P | AP

###  III.3 Electronic structure - Van Hove singularities

Figure 4: Calculated electronic structure of DyTe1.8 and Van Hove
singularities near the Fermi level. (a) Band structure of the $\\{0,4\\}$
divacancy configuration along high symmetry lines. (b) BZ of the
$\sqrt{5}\times\sqrt{5}$-modulated supercell cell and the high symmetry
points. Van Hove singularities near the Fermi level. (c) Calculated DOS near
the Fermi level showing the emergence of Van Hove singularities at the
energies, $e_{1}$ and $e_{2}$, of the saddle points in both the hole and
electron doped regimes. (d) Two dimensional band structure on the
$\hat{b}_{1}-\hat{b}_{3}$ plane, where the orange and blue colors denote the
Tes-$p$ and Dy-$d$ orbital characters, respectively. Valence (conduction) band
exhibits saddle point(s) at the $\Gamma$ (near the $Z$) point marked with red
curves. The energy levels, $e_{1}$ and $e_{2}$, of the saddle points are
illustrated with transparent planes.

Figure 4 shows the calculated electronic structure of DyTe1.8 with the
$\\{0,4\\}$ divacancy configuration. In Fig. 4(b) we plot the Brillouin zone
(BZ) of the $\sqrt{5}\times\sqrt{5}$-modulated cell of DyTe1.8 with the high
symmetry directions used in the band structure plot. The blue, red and orange
colors denote projections on the Dy-$d$, Tes-$p$, and Tec-p orbitals. The
indirect (direct) band gap is 0.14 (0.33) eV in qualitative agreement with the
activation gap measured in transport.[7]

In our previous study[7] of the electronic structure of pristine DyTe2, we
found that the Tes-$p$ bands, originating from the square lattice, intersect
the Fermi level and give rise to the nesting peaks near the $\bm{q}$ vectors
corresponding to the $\sqrt{5}\times\sqrt{5}$ cell modulation. These nesting
peaks explain the cell modulation and resulting insulating band gap in
DyTe1.8, which is further enhanced by the formation of the vacancy network,
which increases the inter-nested-states coupling strengths. The topmost
valence band of DyTe1.8 retains the 2D nature of the square lattice, with
highly dispersive bands along the in-plane direction and weakly dispersive
bands along the out-of-plane direction such as the $\Gamma-Z$ and $R-S$
symmetry lines. The dimensional reduction in the band structure likely
introduces saddle points, leading to Van Hove singularities (VHS). Figure 4(d)
illustrates the saddle points on the $\bm{b}_{1}-\bm{b}_{3}$ plane and the
corresponding divergent peaks of the density of states (DOS). The conduction
band, dominated by the Dy-$d$ orbital character, also exhibits saddle points
on this plane which slightly deviate from the $Z$ point by $\pm
0.1\bm{b}_{3}$. These VHS in both the valence and conduction bands are
expected to give rise to significant optical responses due to the high joint
DOS, where the energy difference, $e_{2}-e_{1}=0.44$ eV, corresponds to the
infrared wavelength, $\lambda=2.8\,\mu$m. However, these saddle points are
separated in momentum space, necessitating a phonon coupling to complement the
momentum difference, $\Delta k_{z}\sim 0.4|\bm{b}_{3}|$. Fortunately, the
system undergoes another double cell modulation along the $\hat{z}$ direction,
which in turn folds the $Z$ point to the $\Gamma$ point. Thus, in the long-
wavelength limit, optical absorption between the two saddle points can occur
without phonon assistance. Additionally, the vertical polarization of each
square lattice layer induces a Rashba spin splitting in the presence of spin-
orbit coupling. The alignment of these layers’ polarizations determines the
spin texture chirality of the spin-split bands, influencing the band alignment
and consequent optical responses. Further theoretical and experimental
investigations are needed to fully understand these effects.

###  III.4 Electronic structure - Fermi surface

Figure 5: Asymmetric Fermi surface (FS) sheets (left columns) and their (001)
projections (right columns), of DyTe1.8 with the $\\{0,4\\}$ divacancy
configuration for chemical potential shifts, $\mu=E-E_{F}$, ranging from (a)
$+0.3$ eV to (e) $-0.5$ eV in steps of 0.2 eV. Orange and blue colors denote
the Tes-$p$ and Dy-$d$ orbital characters, respectively. (f) FS sheets (left
panel) and their (001) projections (right panel) for the $\\{1,3\\}$ divacancy
configuration for $\mu=-0.5$ eV. (g) Fermi surface at $\mu=-0.1$ eV. Green
circular area shows the $\hat{k}_{x}-\hat{k}_{y}$ plane on which the external
magnetic field is swept with an azimuthal angle $\phi$ from the $\hat{k}_{x}$
axis. The corresponding maximal cross-section paths under the external field
along $B_{1}$ and $B_{2}$ are illustrated with blue and red stripes,
respectively, on the Fermi surface. (h) Calculated maximal cross-section area
versus $\phi$.

Figures 5(a)-(e) display the calculated 3D Fermi surface (FS) sheets and their
(001) projections, of DyTe1.8 with the $\\{0,4\\}$ divacancy configuration for
chemical potential shifts, $\mu=E-E_{F}$, ranging from (a) $+0.3$ eV to (e)
$-0.5$ eV in steps of 0.2 eV. The orange and blue colors denote the Tes-$p$
and Dy-$d$ orbital characters, respectively. For $\mu=-0.5$ eV the chemical
potential lies in the Tes-$p$-derived valence bands [see Fig. 4(a)], and the
FS consists of highly-nested parallel sheets stacked normal to the $\Gamma$-S
symmetry direction of the BZ [Fig. 4(b)]. As the chemical potential shifts
upward in energy, the FS area gradually gets reduced and isolated pockets
appear at $Z$ point. Above the Fermi level, two elongated electron pockets
emerge on the side surface of the BZ as shown in Fig. 5(b) and eventually a
Dy-$d$ derived electron pocket emerges at $\Gamma$ point [Fig. 5(a)]. In Fig.
5(f) we display the FS sheets and their (001) projections for the $\\{1,3\\}$
divacancy configuration for $\mu=-0.5$ eV. Comparison of Figs. 5(e) and (f)
clearly shows the 90∘ rotation of the FS associated with the Tes dimer
rotation by 90∘ about [001], shown in Fig. 2(e). Since the group velocity is
always normal to the FS, the parallel FS sheets indicate 1D conducting
channels where their direction can be switched by the Te dimer orientations,
that will be discussed in the next section.

As alluded above, the valence bands are dominated by the Tes atoms of the 2D
square lattice. The Te-dimer orientations switch the hopping directions within
the square lattice that further refine the 2D like band structure into an
effective 1D channel resulting in asymmetric Fermi surface as shown in Fig. 5.
The FS asymmetry can be confirmed by measuring the quantum oscillations under
external magnetic field whose frequency provides information of the maximal
cross-section area of the FS. Figure 5(g) shows the disk shape Fermi surface
centered at the $Z$ point under hole doping for $\mu=-0.1$ eV. The maximal
cross-section area shown in Fig. 5(h) is expected to strongly depend on the
magnetic field direction, where a peak appears at the azimuthal angle
(measured from $\hat{k}_{x}$), $\phi\sim 103^{\circ}$, corresponding to
$B_{2}$. The red stripe on the Fermi surface shows the cyclotron motion path
under the external field $B_{2}$ while the blue stripe is under $B_{1}$ whose
cross-section area is smaller than that of $B_{2}$ by a factor of five.

###  III.5 Quasi 1D conducting channel

Figure 6: Effect of uniaxial in-plane strain on Te-dimer orientation. (a)
Schematic figure of DyTe1.8 epitaxially grown on a piezoelectric substrate.
Small (large) square shows the in-plane unit cell of primitive DyTe2
($\sqrt{5}\times\sqrt{5}$-modulated DyTe1.8) structure. (b) Landscape of total
energy difference, $\Delta E=E(\\{0,4\\})-E(\\{1,3\\})$, between the
$\\{0,4\\}$\- and $\\{1,3\\}$ divacancy structural motifs of the
$\sqrt{5}\times\sqrt{5}$ supercell on two-dimensional epitaxial strain,
$r_{x}\equiv\frac{\alpha_{x}}{\alpha_{0}}$ and
$r_{y}\equiv\frac{\alpha_{y}}{\alpha_{0}}$. Yellow and blue dashed curves
indicate energy contours of the $\\{0,4\\}$ and $\\{1,3\\}$ phases,
respectively. (c) and (d) Relatively stable Te-dimer orientations under
uniaxial strain. The $\\{0,4\\}$ ($\\{1,3\\}$) phase is preferred under
uniaxial strain along the $\hat{y}$ ($\hat{x}$) direction. Gold (black) leads
indicate higher (lower) conductivity channel. (e) Calculated conductivity of
DyTe1.8 with the two vacancies at $\\{0,4\\}$ versus the chemical potential
$\mu$. Conductivity along $\hat{y}$ [see (a) or Fig. 2(a) for the direction]
is higher more than one order of magnitude compared to the other directions
under hole doping. The dashed line denotes the in-plane conductivity ratio,
$\sigma_{yy}/\sigma_{xx}$. The band gap near $\mu=0$ is marked with the gray
box.

The results of the various reconstructed iso-energetic divacancy
configurations invite the intriguing question whether epitaxial strain, which
is ubiquitous in DyTe1.8 samples grown on atomically flat piezoelectric
substrates, shown schematically in Fig. 6(a), can stabilize a specific ordered
structural motif, and whether the divacancy configuration network can undergo
a transition from one motif to another. Fig. 6(b) shows the landscape of the
total energy difference, $\Delta E=E(\\{0,4\\})-E(\\{1,3\\})$, between the
$\\{0,4\\}$\- and $\\{1,3\\}$ divacancy structural motifs of the
$\sqrt{5}\times\sqrt{5}$ reconstructed DyTe1.8 on the two-dimensional
epitaxial strain, $r_{x}\equiv\frac{\alpha_{x}}{\alpha_{0}}$ and
$r_{y}\equiv\frac{\alpha_{y}}{\alpha_{0}}$. Here, $\alpha_{0}$ is the in-plane
equilibrium lattice constant of the square lattice, and the x and y refer to
the [110] and [$\bar{1}$10] directions, respectively. The calculations reveal
that tensile epitaxial strain, $\epsilon_{yy}\approx$ 0.7%, renders the
$\\{0,4\\}$ divacancy-based configuration as the ground state [Fig. 6(c)]. On
the other hand, the $\\{1,3\\}$ divacancy-based structural motif becomes the
ground state under tensile strain, $\epsilon_{xx}\approx$ 0.7% along [110]
[Fig. 6(d)]. Interestingly, these two distinct Tes-vacancy ordered networks
have substantial different transport properties. Figure 6(e) shows the matrix
elements of the longitudinal conductivity, $\sigma_{ii}$ ($i$=x,y,z), of
DyTe1.8 with the $\\{0,4\\}$ divacancy configuration versus the chemical
potential, $\mu$. The conductivity ratio, $\frac{\sigma_{yy}}{\sigma_{xx}}$
versus $\mu$ is also plotted as a dashed line whose scale appears on the right
vertical axis. Interestingly, $\sigma_{yy}$, is higher by more than one order
of magnitude compared to $\sigma_{xx}$ under hole doping.

This notable conductivity switching between the two phases upon the Te dimer
orientation, suggests potential application as a memory device. Control of the
dimer orientation can be achieved by applying epitaxial strain via a
piezoelectric substrate, which can in turn align all or most of the sample’s
Te dimers along the $\hat{x}$ or $\hat{y}$ directions, resulting in a colossal
conductivity contrast, ranging from 1000% to 2500% between the ON and OFF
phases.

Consequently, the Te dimer in the supercell acts as a molecular switch not
only in a macroscopic scale but also in a microscopic level by controlling the
electron hopping directions within the $\sqrt{5}\times\sqrt{5}$ supercell. The
overall conductivity of the two-dimensional percolation network[36, 37] will
be determined by the microscopic configurations of the Te dimer orientations,
which can be tuned by epitaxial strain. We acknowledge the significance of
prior nanotube percolation experiments[38, 39], particularly those where the
nanotube orientation was controllable, which garnered considerable attention
within the research community. Given the potential for strain-controlled dimer
orientations, DyTe1.8 emerges as a promising alternative platform for
conducting percolation experiments.

###  III.6 Migration Properties

Figure 7: (a-h) Te-square net with various di-vacancy configurations and their
migrations. Background blue (orange) dots denote the positions of Dy (Tec)
atoms below the Te-square net layer. Arrows illustrate the considered paths of
ionic migrations. Four different types of migration paths are colored by red,
green, cyan, and violet. (i,j) Calculated migration energy landscape of the
various transitions for (i) symmetric, $r_{1}=r_{2}=1.0$ and (j) asymmetric,
$r_{1}=0.988$, $r_{2}=1.012$, in-plane strain.

In order to understand the diffusion properties of the divacancy
configurations in the presence of absence of epitaxial strain, we have carried
first principles calculations using the nudged elastic band (NEB) method[40,
41, 42, 43] to determine the energy barriers and the minimum energy paths
associated with the atomic migrations within the Te square network. The energy
saddle point identified by this method represents the barrier between two
stable configurations and the migration is assumed to occur simultaneously
across all periodic (both in-plane and out-of-plane) cells, rendering the
calculated energy barrier an upper bound. Consequently, the energy required to
incrementally expand a domain area could be less than the barrier determined
by this method.

There are numerous atomic migration paths, some of which are symmetric
equivalents with identical energy barriers. Therefore, only four irreducible
migration paths require computation which collectively represent all possible
migrations. The eight dotted arrows depicted in the Fig. 7(a), (b), (e), and
(f) represent these four migration paths (colored red, green, cyan, and
purple) adjacent to the vacant sites. Paths of the same color are reverse
equivalents of each other. A pairwise combination of these migrations can
describe both polarization reversal and dimer rotation.

For instance, the migration from (a) $\\{5,9\\}$ to (b) $\\{5,6\\}$ to (c)
$\\{2,6\\}$ changes the position of the dimer and thus the polarization
direction. In the first step, the right side Te of a dimer moves to the vacant
site, denoted by the red dotted arrow in Fig. 7(a) where the intermediate
positions of the Te atom are illustrated as transparent red spheres in Fig.
7(b). Next, the middle Te of a trimer moves upwards, marked as red dotted
arrow in Fig. 7(b) that completes the migration path to the opposite
polarization [Fig. 7(c)]. If the second step migration is replaced by the
purple dotted arrow [Fig. 7(b)] where the end Te atom moves instead of the
middle one of the trimer, it arrives to the (d) $\\{3,5\\}$ configuration
exhibiting opposite polarization and rotated dimer orientation compared to the
initial state (a) $\\{5,9\\}$.

The calculated energies during each migration step are shown in Fig. 7(i)
where the red (purple) colored line corresponds to the migration marked with
the red (purple) arrows with an energy barrier of 0.75 eV (0.37 eV). The
horizontal axis is for images in the NEB method between stable configurations
listed below the tics. The possible migration paths between the configurations
are denoted as colored lines connecting them. For an instance, migrations from
(e) $\\{6,8\\}$ to (f) $\\{5,8\\}$ to (a) $\\{5,9\\}$ are connected by cyan
and green lines below Fig. 7(i) corresponding to the energy curves colored
cyan and green, respectively.

We note that the previously considered polarization reversal migration,
$\\{5,9\\}-\\{5,6\\}-\\{2,6\\}$ takes the red energy curves back and forth
passing the 0.75 eV energy barrier twice. The Te-dimer then moves towards the
$-\hat{y}$ direction. The polarization can be reversed by another migration of
$\\{5,9\\}-\\{5,8\\}-\\{4,8\\}$ towards the $-\hat{x}$ direction that needs to
overcome the green energy curves passing the 0.37 eV barrier twice instead.
The latter one experiences a lower energy barrier indicating that there is a
preferred direction in the dimer migration which is determined by the
supercell chirality and the dimer orientation.

The migration energies are considered within square symmetry, where pairs of
paths exhibit the same energy barrier. Figure 7 (j) shows the change of the
migration energy landscape under uniaxial strain $r_{2}>r_{1}$ (see the
caption). The previously degenerate total energies of the ground-state
divacancy configurations are lifted, and the $\\{5,9\\}$ configuration becomes
more stable than the $\\{6,8\\}$ due to distinct Te-dimer orientations.
Additionally, the metastable $\\{5,8\\}$ configuration becomes more stable
than its $\\{5,6\\}$ counterpart.

##  IV Summary

The structural and electronic properties of DyTe1.8 are investigated by
employing first-principles calculations. In addition to the
$\sqrt{5}\times\sqrt{5}$ modulation, the detailed position and orientation of
the Te dimers are found to play a crucial role in determining the polarization
and the electronic structure near the Fermi level. The predicted small
insulating gap with VHS in both the valence and conduction bands imply
characteristic optical responses with fine structures related with the
polarization direction of each Te square lattice layer. The Te-dimer-
orientation-induced asymmetric Fermi surface may be directly confirmed by
quantum oscillation measurements and paves a way to memory device application
with mechanically switchable quasi one-dimensional conductivity. Finally, the
ionic migration paths and energy landscape results reveal a plethora of
structural transitions between iso-energetic motifs. Our predictions of these
interesting and diverse properties of DyTe1.8 will hopefully motivate further
experimental and theoretical investigations.

###### Acknowledgements.

We appreciate insightful discussions with Joseph Falson and Adrian Llanos. The
work is supported by the NSF-PREP CSUN/Caltech-IQIM Partnership (grant number
2216774) and the NSF-Partnership in Research and Education in Materials (PREM)
(grant number DMR-1828019).

## References

  * Shin  _et al._ [2005] K. Y. Shin, V. Brouet, N. Ru, Z. X. Shen, and I. R. Fisher, Electronic structure and charge-density wave formation in $\mathrm{La}\mathrm{Te}_{1.95}$ and $\mathrm{Ce}\mathrm{Te}_{2.00}$, Phys. Rev. B 72, 085132 (2005).
  * Garcia  _et al._ [2007] D. R. Garcia, G.-H. Gweon, S. Y. Zhou, J. Graf, C. M. Jozwiak, M. H. Jung, Y. S. Kwon, and A. Lanzara, Revealing charge density wave formation in the $\mathrm{LaTe}_{2}$ system by angle resolved photoemission spectroscopy, Phys. Rev. Lett. 98, 166403 (2007).
  * Shin  _et al._ [2000] Y. Shin, C. Han, B. Min, H. Lee, C. Choi, Y. Kim, D. Kim, and Y. Kwon, Anisotropic magnetization in $\it{R}\mathrm{Te}_{2}$ ($\it{R}$: $\mathrm{Ce}$, $\mathrm{Pr}$, $\mathrm{Gd}$ and $\mathrm{Sm}$), Physica B: Condensed Matter 291, 225 (2000).
  * Jung  _et al._ [2003] M. H. Jung, A. Alsmadi, H. C. Kim, Y. Bang, K. H. Ahn, K. Umeo, A. H. Lacerda, H. Nakotte, H. C. Ri, and T. Takabatake, Superconductivity in magnetically ordered $\mathrm{CeTe}_{1.82}$, Phys. Rev. B 67, 212504 (2003).
  * Tremel and Hoffmann [1987] W. Tremel and R. Hoffmann, Square nets of main-group elements in solid-state materials, Journal of the American Chemical Society 109, 124 (1987).
  * Kikuchi [1998] A. Kikuchi, Electronic structure of lanthan ditellurides, Journal of the Physical Society of Japan 67, 1308 (1998), https://doi.org/10.1143/JPSJ.67.1308 .
  * Llanos  _et al._ [2024] A. Llanos, S. Salmani-Rezaie, J. Kim, N. Kioussis, D. A. Muller, and J. Falson, Supercell formation in epitaxial rare-earth ditelluride thin films, Crystal Growth & Design 24, 115 (2024).
  * Poddig  _et al._ [2021] H. Poddig, P. Gebauer, K. Finzel, K. Stöwe, and T. Doert, Structural variations and bonding analysis of the rare-earth metal tellurides $\it{RE}\mathrm{Te}_{1.875\pm\delta}$ ($\it{RE}$ = $\mathrm{Ce}$, $\mathrm{Pr}$, $\mathrm{Sm}$, $\mathrm{Gd}$; $0.004\leq\delta\leq 0.025$), Inorganic Chemistry 60, 11231 (2021).
  * Plambeck-Fischer  _et al._ [1989] P. Plambeck-Fischer, W. Abriel, and W. Urland, Preparation and crystal structure of $\it{RE}\mathrm{Se}_{1.9}$ ($\it{RE}$ = $\mathrm{Ce}$, $\mathrm{Pr}$), Journal of Solid State Chemistry 78, 164 (1989).
  * Dashjav  _et al._ [2000] E. Dashjav, O. Oeckler, T. Doert, H. Mattausch, and P. Böttcher, $\mathrm{Gd}_{8}\mathrm{Se}_{15}$—$\mathrm{A}$ 24-fold superstructure of the $\mathrm{ZrSSi}$ type, Angewandte Chemie International Edition 39, 1987 (2000), https://onlinelibrary.wiley.com/doi/pdf/10.1002/1521-3773 .
  * van der Lee  _et al._ [1997] A. van der Lee, L. M. Hoistad, M. Evain, B. J. Foran, and S. Lee, Resolution of the 66-fold superstructure of $\mathrm{DySe}_{1.84}$ by x-ray diffraction and second-moment scaled hückel calculations, Chemistry of Materials 9, 218 (1997).
  * Graf and Doert [2009] C. Graf and T. Doert, $\mathrm{LaSe}_{1.85}$, $\mathrm{CeSe}_{1.83}$, $\mathrm{NdSe}_{1.83}$ and $\mathrm{SmSe}_{1.84}$ – four new rare earth metal polyselenides with incommensurate site occupancy and displacive modulation, Zeitschrift für Kristallographie 224, 568 (2009).
  * Lee and Foran [1996] S. Lee and B. J. Foran, Defective lattice charge density waves in $\mathrm{La}_{10}\mathrm{Se}_{19}$, $\mathrm{Cs}_{3}\mathrm{Te}_{22}$, $\mathrm{RbDy}_{3}\mathrm{Se}_{8}$, and $\mathrm{Dy}_{65.33}\mathrm{Se}_{120}$, Journal of the American Chemical Society 118, 9139 (1996).
  * Ijjaali and Ibers [2006] I. Ijjaali and J. A. Ibers, Two new binary lanthanide polytellurides: Syntheses and crystal structures of $\mathrm{CeTe}_{1.90}$ and $\mathrm{SmTe}_{1.80}$, Journal of Solid State Chemistry 179, 3456 (2006).
  * Kresse and Furthmüller [1996a] G. Kresse and J. Furthmüller, Efficient iterative schemes for $\it{abinitio}$ total-energy calculations using a plane-wave basis set, Phys. Rev. B 54, 11169 (1996a).
  * Kresse and Furthmüller [1996b] G. Kresse and J. Furthmüller, Efficiency of ab-initio total energy calculations for metals and semiconductors using a plane-wave basis set, Computational Materials Science 6, 15 (1996b).
  * Blöchl [1994] P. E. Blöchl, Projector augmented-wave method, Phys. Rev. B 50, 17953 (1994).
  * Kresse and Joubert [1999] G. Kresse and D. Joubert, From ultrasoft pseudopotentials to the projector augmented-wave method, Phys. Rev. B 59, 1758 (1999).
  * Perdew  _et al._ [2008] J. P. Perdew, A. Ruzsinszky, G. I. Csonka, O. A. Vydrov, G. E. Scuseria, L. A. Constantin, X. Zhou, and K. Burke, Restoring the density-gradient expansion for exchange in solids and surfaces, Phys. Rev. Lett. 100, 136406 (2008).
  * Mostofi  _et al._ [2014] A. A. Mostofi, J. R. Yates, G. Pizzi, Y.-S. Lee, I. Souza, D. Vanderbilt, and N. Marzari, An updated version of wannier90: A tool for obtaining maximally-localised wannier functions, Comput. Phys. Commun. 185, 2309 (2014).
  * Marzari  _et al._ [2012] N. Marzari, A. A. Mostofi, J. R. Yates, I. Souza, and D. Vanderbilt, Maximally localized Wannier functions: Theory and applications, Rev. Mod. Phys. 84, 1419 (2012).
  * Becke and Johnson [2006] A. D. Becke and E. R. Johnson, A simple effective potential for exchange, The Journal of Chemical Physics 124, 10.1063/1.2213970 (2006), 221101, https://pubs.aip.org/aip/jcp/article-pdf/doi/10.1063/1.2213970/15385734/221101_1_online.pdf .
  * Tran and Blaha [2009] F. Tran and P. Blaha, Accurate band gaps of semiconductors and insulators with a semilocal exchange-correlation potential, Phys. Rev. Lett. 102, 226401 (2009).
  * Kim  _et al._ [2019] J. Kim, K. M. Rabe, and D. Vanderbilt, Negative piezoelectric response of van der waals layered bismuth tellurohalides, Phys. Rev. B 100, 104115 (2019).
  * Po  _et al._ [2018] H. C. Po, H. Watanabe, and A. Vishwanath, Fragile topology and wannier obstructions, Phys. Rev. Lett. 121, 126402 (2018).
  * Radha and Lambrecht [2021] S. K. Radha and W. R. L. Lambrecht, Topological obstructed atomic limit insulators by annihilating dirac fermions, Phys. Rev. B 103, 075435 (2021).
  * Li  _et al._ [2022] X.-P. Li, D.-S. Ma, C.-C. Liu, Z.-M. Yu, and Y. Yao, From atomic semimetal to topological nontrivial insulator, Phys. Rev. B 105, 165135 (2022).
  * Xu  _et al._ [2024] Y. Xu, L. Elcoro, Z.-D. Song, M. G. Vergniory, C. Felser, S. S. P. Parkin, N. Regnault, J. L. Mañes, and B. A. Bernevig, Filling-enforced obstructed atomic insulators, Phys. Rev. B 109, 165139 (2024).
  * King-Smith and Vanderbilt [1993] R. D. King-Smith and D. Vanderbilt, Theory of polarization of crystalline solids, Phys. Rev. B 47, 1651 (1993).
  * Vanderbilt and King-Smith [1993] D. Vanderbilt and R. D. King-Smith, Electric polarization as a bulk quantity and its relation to surface charge, Phys. Rev. B 48, 4442 (1993).
  * Resta [1994] R. Resta, Macroscopic polarization in crystalline dielectrics: the geometric phase approach, Rev. Mod. Phys. 66, 899 (1994).
  * Zak [1989] J. Zak, Berry’s phase for energy bands in solids, Phys. Rev. Lett. 62, 2747 (1989).
  * Aihara  _et al._ [2020] Y. Aihara, M. Hirayama, and S. Murakami, Anomalous dielectric response in insulators with the $\pi$ zak phase, Phys. Rev. Research 2, 033224 (2020).
  * Smeu  _et al._ [2012] M. Smeu, H. Guo, W. Ji, and R. A. Wolkow, Electronic properties of si(111)-$7\times 7$ and related reconstructions: Density functional theory calculations, Phys. Rev. B 85, 195315 (2012).
  * Kim  _et al._ [2023] J. Kim, C.-Y. Huang, H. Lin, D. Vanderbilt, and N. Kioussis, Bismuth antiphase domain wall: A three-dimensional manifestation of the su-schrieffer-heeger model, Phys. Rev. B 107, 045135 (2023).
  * Sahimi [2023] M. Sahimi,  _Application of Percolation Theory_, Vol. 213 (SpringerCham, 2023) pp. XXI,680.
  * Cherkasova  _et al._ [2010] V. A. Cherkasova, Y. Y. Tarasevich, N. I. Lebovka, and N. V. Vygornitskii, Percolation of aligned dimers on a square lattice, The European Physical Journal B 74, 205 (2010).
  * Du  _et al._ [2005] F. Du, J. E. Fischer, and K. I. Winey, Effect of nanotube alignment on percolation conductivity in carbon nanotube/polymer composites, Phys. Rev. B 72, 121404(R) (2005).
  * Park  _et al._ [2006] C. Park, J. Wilkinson, S. Banda, Z. Ounaies, K. E. Wise, G. Sauti, P. T. Lillehei, and J. S. Harrison, Aligned single-wall carbon nanotube polymer composites using an electric field, Journal of Polymer Science Part B: Polymer Physics 44, 1751 (2006).
  * Mills  _et al._ [1995] G. Mills, H. Jónsson, and G. K. Schenter, Reversible work transition state theory: application to dissociative adsorption of hydrogen, Surface Science 324, 305 (1995).
  * [41] H. Jónsson, G. Mills, and K. W. Jacobsen, Nudged elastic band method for finding minimum energy paths of transitions, in  _Classical and Quantum Dynamics in Condensed Phase Simulations_, pp. 385–404.
  * Henkelman  _et al._ [2000] G. Henkelman, B. P. Uberuaga, and H. Jónsson, A climbing image nudged elastic band method for finding saddle points and minimum energy paths, The Journal of Chemical Physics 113, 9901 (2000).
  * Henkelman and Jónsson [2000] G. Henkelman and H. Jónsson, Improved tangent estimate in the nudged elastic band method for finding minimum energy paths and saddle points, The Journal of Chemical Physics 113, 9978 (2000).