# Tight-Binding Bandstructure of $\beta-$ and $\alpha-$ phase Ga2O3 and Al2O3

Y. Zhang<EMAIL_ADDRESS>Electrical and Computer Engineering, Cornell
University, Ithaca, NY 14850, USA  M. Liu  Materials Science and Engineering,
Cornell University, Ithaca, NY 14850, USA  G. Khalsa<EMAIL_ADDRESS>Materials Science and Engineering, Cornell University, Ithaca, NY 14850, USA
D. Jena  Electrical and Computer Engineering, Cornell University, Ithaca, NY
14850, USA  Materials Science and Engineering, Cornell University, Ithaca, NY
14850, USA

###### Abstract

Rapid design and development of the emergent ultra-wide bandgap semiconductors
Ga2O3 and Al2O3 requires a compact model of their electronic structures,
accurate over the broad energy range accessed in future high-field, high-
frequency, and high-temperature electronics and visible and ultraviolet
photonics. A minimal tight-binding model is developed to reproduce the first-
principles electronic structures of the $\beta-$ and $\alpha-$ phases of Ga2O3
and Al2O3 throughout their reciprocal spaces. Accurately reproducing the
bandgap, orbital character, and effective mass and high-energy features of the
conduction band, this compact model will assist in the investigation and
design of the electrical and optical properties of bulk materials, devices,
and quantum confined heterostructures.

††preprint: AIP/123-QED

##  I Introduction

The recent integration of Ga2O3 with Al2O3 has the potential to revolutionize
high-power electronics. The availability of large, inexpensive, single-crystal
substrates Higashiwaki  _et al._ (2012); Kuramata  _et al._ (2016), recent
advances in thin film growthGogova  _et al._ (2014); Okumura  _et al._ (2014);
Jinno  _et al._ (2021), and the ability to dope these wide-bandgap
semiconductors have enabled transistors and Schottky diodes based on
$\beta$-Ga2O3 with breakdown fields as large as 5.45 MV/cm (Ref. Roy  _et al._
, 2021) and 5.7 MV/cm (Ref. Kalarickal  _et al._ , 2021) and approaching the
projected theoretical estimate of 8 MV/cm (Ref. Pearton  _et al._ , 2018).
Comparing these breakdown fields with the existing technological
semiconductors Si (0.3 MV/cm), SiC (3.1 MV/cm), and GaN (3.3 MV/cm) Meneghini,
Meneghesso, and Zanoni (2017), $\beta$-Ga2O3 promises new high-frequency,
high-voltage, and high-temperature electronics applications. $\alpha$-Ga2O3
and $\alpha$-Al2O3 further expand the bandgap to 5.2 eV and 8.8 eV, signifying
the potential for oxide semiconductors to expand the future electronics
materials tool-set.

The successful design of these future electronic devices requires an accurate
modeling and understanding of the electronic structure and bonding of Ga2O3
and Al2O3. The tight-binding method provides a flexible, chemically motivated
description of the electronic structure of materialsHarrison (1989). When
compared with modern computational approaches to materials physics like
density functional theory (DFT), tight-binding models are compact, intuitive,
and require less computational resources. As a result, tight-binding models
are ubiquitous in device engineering and development and have successfully
described electronic transportDatta (2005); Lake  _et al._ (1997); Nehari  _et
al._ (2005); Popov (2004); Sols  _et al._ (1989) and optical propertiesGraf
and Vogl (1995); Popov (2004); Pi and Delerue (2013); Jirauschek and Kubis
(2014) of bulk materials, heterostructures, and devices. To aid in the
development of new high-power electronics, we derive semi-empirical tight-
binding models in this work for three technologically relevant oxide
semiconductors: $\beta$-Ga2O3, $\alpha$-Ga2O3, and $\alpha$-Al2O3.

While we are unaware of a tight-binding model describing these three oxide
semiconductors, a recent study reports a tight-binding model of $\beta$-Ga2O3
using atomic orbitals as a basis, with parameters drawn from DFT
calculationsLee  _et al._ (2019). The authors employ the model to study the
surface energy of $\beta$-Ga2O3 and formation energy of Ga and O vacancy
defects. We derive an alternative tight-binding model with the goal of
accurate parameterization of the conduction band and fundamental optical gaps
of $\beta$-Ga2O3, $\alpha$-Ga2O3 and $\alpha$-Al2O3 so that electrical and
optical properties can be faithfully simulated.

We derive tight-binding models using a Wannier functions basis. Wannier
functions are a convenient basis for tight-binding models because they are
derived from the underlying band structure of the material, are formally
orthogonal, can be localized to atomic sites, and preserve the site symmetry
and coordination. This approach of DFT-derived tight-binding has been used
successfully to describe the electronic structure of broad classes of
technologically important materials including siliconSouza, Marzari, and
Vanderbilt (2002), III-V semiconductorsGresch  _et al._ (2018), and 2D
materialsSouza  _et al._ (2012).

##  II Structural and Electrical Properties

The crystal symmetry and bonding environment constrain the tight-binding
description of the electronic structure. When compared to a conventional
semiconductor like Si, Ga2O3 and Al2O3 have relatively low symmetry and
complicated bonding networks. $\beta$-Ga2O3 has a monoclinic structure (space
group C2/m, No. 12). The monoclinic structure contains two pairs of symmetry
inequivalent Ga sites, each coordinated by O, forming two distorted GaO4
tetrahedra and two distorted GaO6 octahedra per unit cell (Fig. 1(a)).
$\alpha$-Al2O3 and $\alpha$-Ga2O3 crystallize in the sapphire structure
(rhombohedral, space group R-3c, No. 167). In the $\alpha$ phase, the Al(Ga)
atoms occupy 4 equivalent sites, each coordinated by six O, forming distorted
AlO6(GaO6) octahedra (Fig. 1(b)). The structural information obtained from DFT
structural optimization and experimental data are given in Table 1 where DFT
is shown to describe the experimental structure surprisingly well. The valence
configurations of O, Al, and Ga are 2s22p4, 3s23p1, and 4s24p1, respectively.
The Al(Ga) is expected to donate its valence electrons in order to fill the O
valence shell, leading to an O-2p derived valence band with Al-3s(Ga-4s) and
Al-3p(Ga-4p) derived conduction band.

Figure 1:  Crystal structure of monoclinic $\beta$-Ga2O3, and rhombohedral
$\alpha$-Ga2O3 and $\alpha$-Al2O3. Coordination octahedra and tetraheda
highlight the different Ga(Al)-O bonding environments in the $\beta$ (a) and
$\alpha$ (b) phases. This bonding difference is evident in the symmetry of
Ga(Al)-site Wannier functions (c,d) for the conduction band where clear
hybridization of the Ga(Al)-s and O-p is seen. In (c,d), the positive lobes of
the Wannier functions are shown in yellow and the negative lobes are shown in
blue.

Structural optimization and electronic band structure are calculated by DFT
using Quantum EspressoGiannozzi  _et al._ (2009). We choose projector-
augmented wave pseudopotentials and generalized-gradient approximation of
exchange correlation using the Perdew-Burke-Ernzerhof functional generalized
for solidsPerdew, Burke, and Ernzerhof (1996); Perdew  _et al._ (2008).
Structural convergence is found for a $3\times 3\times 3$ k-point mesh and a
60 Ry plane-wave cutoff. We find good agreement between the DFT-relaxed
structure and the experimental structure (see Table 1).

Table 1:  $\beta$-Ga2O3, $\alpha$-Ga2O3, and $\alpha$-Al2O3 structural data.
In the $\beta$-phase (left column), Ga and O occupy the $i$ Wyckoff site with
fractional coordinates $(x,0,z)$ and $(-x,0,-z)$. In the $\alpha$ phase (right
column), Ga(Al) occupy the $c$ Wyckoff site with fractional coordinates
$(z,z,z)$, $(-z+\frac{1}{2},-z+\frac{1}{2},-z+\frac{1}{2})$, $(-z,-z,-z)$, and
$(z+\frac{1}{2},z+\frac{1}{2},z+\frac{1}{2})$ and $O$ occupies the $e$ Wyckoff
site with fractional coordinates $(x+\frac{1}{4},-x+\frac{1}{4},\frac{1}{4})$,
$(\frac{1}{4},x+\frac{1}{4},-x+\frac{1}{4})$,
$(-x+\frac{1}{4},\frac{1}{4},x+\frac{1}{4})$,
$(-x+\frac{3}{4},x+\frac{3}{4},\frac{3}{4})$,
$(\frac{3}{4},-x+\frac{3}{4},x+\frac{3}{4})$, and
$(x+\frac{3}{4},\frac{3}{4},-x+\frac{3}{4})$.

$\beta$-Ga2O3 |  $a$ (Å) |  $b$ (Å) |  $c$ (Å) |  $\beta$ (∘)
---|---|---|---|---
DFT | 12.237 | 3.062 | 5.813 | 103.81
ExperimentZade  _et al._ (2019) | 12.214 | 3.037 | 5.798 | 103.83
| Ga(I) (4i) | Ga(II) (4i)
---|---|---
$\beta$-Ga2O3 | $x$ | $z$ | $x$ | $z$
DFT | 0.0904 | 0.7949 | 0.6585 | 0.3137
ExperimentZade  _et al._ (2019) | 0.0895 | 0.7938 | 0.6585 | 0.3100
| O(I) (4i) | O(II) (4i) | O(III) (4i)
---|---|---|---
$\beta$-Ga2O3 | $x$ | $z$ | $x$ | $z$ | $x$ | $z$
DFT | 0.1647 | 0.1094 | 0.1733 | 0.5632 | 0.4959 | 0.2563
ExperimentZade  _et al._ (2019) | 0.1519 | 0.1001 | 0.1722 | 0.5640 | 0.4920 | 0.2645

$\alpha$-Ga2O3 |  $a$ (Å) |  $b$ (Å) | Ga(I) 4c | O(I) 6e
---|---|---|---|---
DFT | 5.3221 | 55.82 |  $z$=0.1446 |  $x$=0.3049
ExperimentMarezio and Remeika (1967) | 5.3221 | 55.82 |  $z$=0.1446 |  $x$=0.3049
$\alpha$-Al2O3 |  $a$ (Å) |  $b$ (Å) | Al(I) 4c | O(I) 6e
---|---|---|---|---
DFT | 5.1779 | 55.28 |  $z$=0.1479 |  $x$=0.3056
ExperimentIshizawa  _et al._ (1980) | 5.126 | 55.25 |  $z$=0.3522 |  $x$=0.3064

The DFT band structure and DOS are shown in Fig. 2(c,f,i) and found to be in
agreement with previous workPeelaers and Van de Walle (2015). Orbital
projected DOS shows the valence band is primarily O-2p and the conduction band
is primarily Ga-4s(Al-3s) and Ga-4p(Al-3p). The bottom of the conduction band
is dominated by Ga-4s(Al-3s) states, transitioning to mainly Ga-4p(Al-3p) in
character at 6 eV above the conduction band minimum. All three materials share
several broad features of their band structure: broad O-2p valence bands with
a flat valence band edges, conduction band edge features near the $\Gamma$
point, and indirect bandgaps. In all three materials, the valence band edge is
populated by flat O-2p bands. DFT predicts the valence band maximum between M2
and D in $\beta$-Ga2O3, and between $\Gamma$ and S0 in both $\alpha$-Ga2O3 and
$\alpha$-Al2O3. The conduction band minimum at the $\Gamma$-point has nearly
isotropic dispersion, suggesting that it is insensitive to symmetry and
chemistry. We find an effective mass of 0.26 $m_{e}$ in $\beta$-Ga2O3, 0.28
$m_{e}$ in $\alpha$-Ga2O3, and 0.42 $m_{e}$ in $\alpha$-Al2O3. We found that
the conduction band effective mass of $\beta$-Ga2O3 is in good agreement with
recent transportZhang  _et al._ (2018) and angle-resolved photoemission
spectroscopyMohamed  _et al._ (2010); Janowitz  _et al._ (2011) (ARPES)
experiments. Away from the band minimum, the conduction band transitions from
parabolic to linear dispersion, a band structure feature used to explain high
field transportGhosh and Singisetti (2017, 2018) and optical absorptionSingh
_et al._ (2020) experiments. The second conduction band at the $\Gamma$ point
is above the conduction band minimum by 3.3 eV in $\beta$-Ga2O3, 3.5 eV in
$\alpha$-Ga2O3, and 3.1 eV in $\alpha$-Al2O3, respectively. The value of
$\beta$-Ga2O3 agrees with the experimental observation of 3.55 eVSingh  _et
al._ (2020). DFT predicts an indirect bandgap with the difference between the
direct and indirect bandgap less than 20 meV. This is confirmed by the ARPES
measurementMohamed  _et al._ (2010); Janowitz  _et al._ (2011). As expected,
the bandgap predicted by DFT underestimates the experimental bandgaps. In
order to match the experimental bandgap, we have conducted a “scissor cut” by
shifting the conduction band states up in energy so that the bandgaps match
the experimental values of 4.9 eVJanowitz  _et al._ (2011), 5.2 eVShinohara
and Fujita (2008), and 8.8 eVFrench (1990) in $\beta$-Ga2O3, $\alpha$-Ga2O3,
and $\alpha$-Al2O3, respectively. In constructing a tight-binding model, we
aim to describe the above key features of the DFT band structure.

##  III Tight-Binding Derivation

We derive Wannier functions from the DFT band structure using Wannier90Pizzi
_et al._ (2020) to construct the tight-binding basis. Wannier functions are
initialized on the Ga(Al) sites with s- and p-orbital symmetry and on the O
sites with p-orbital symmetry. We selectively localize only the Ga-s(Al-s)
Wannier functions onto Ga(Al) sites using a Lagrange multiplierWang  _et al._
(2014). Selectively-localized Wannier functions improve the localization of
Ga-s(Al-s) Wannier functions at the cost of delocalizing Ga-p(Al-p) Wannier
functions which describe the higher conduction bands (> 5 eV above the
conduction band minimum). We find that the O-p Wannier functions have similar
localization in both schemes.

We show isosurface plots of Ga-s Wannier functions in $\beta$ and $\alpha$
phases resulting from selective localization in Fig. 1(c,d). (The Al-s and
Ga-s Wannier functions in the $\alpha$ phase are qualitatively similar.) The
Wannier functions reproduce the distorted tetrahedral and octahedral symmetry
of the coordination polyhedra as expected from the hybridization between the
Ga-4s and O-2p atomic orbitals in Ga2O3.

We construct the tight-binding model from the Wannier function basis by
extracting the Hamiltonian matrix element using PythTBCoh and Vanderbilt
(2013). We parameterize our tight-binding model in the following form:

$\begin{split}\hat{H}&=\sum_{i}\epsilon_{i}\left.|i\right>\left<i|\right.\\\
&+\sum_{i,j}t_{ij}e^{-i\vec{k}\cdot\vec{\Delta
d_{ij}}}\left.|j\right>\left<i|\right.+c.c.\end{split}$ (1)

The right-hand side of the first line describes the contribution of individual
Wannier functions to the total energy, commonly called the onsite energy. The
second line describes the kinetic energy of the electrons, commonly called the
“hopping" energy. In Equation 1, $\left.|i\right>$ represents the $i^{th}$
Wannier function with on-site energy $\epsilon_{i}$.
$t_{ij}=\tilde{t}_{ij}e^{i\phi_{ij}}$ is the hopping term between $i^{th}$ and
$j^{th}$ Wannier functions. Here $\phi_{ij}$ characterizes the phase of
$t_{ij}$ and encompasses the symmetry of the local bonding environment.
$\vec{k}$ is the wavevector, and $\Delta\vec{d}_{ij}$ is the displacement
vector of the two Wannier functions $i$ and $j$ defined as $\vec{\Delta
d_{ij}}=\vec{d_{j}}-\vec{d_{i}}+\vec{R}$. Here $\vec{d}_{i}$ and $\vec{d}_{j}$
are the centers of Wannier function $i$ and $j$ within the same unit cell, and
$\vec{R}$ is the unit cell translation to indicate coupling across adjacent
unit cells.

In our tight-binding model, we include only the s-orbital-derived Wannier
functions from the 4 Ga(Al) atoms and the 3 p-orbital-derived Wannier
functions from the 6 O atoms. This translates to a 22 $\times$ 22 matrix when
the model is implemented numerically. To describe the band structure with a
minimal set of parameters, we include: 1) the onsite energy ($\epsilon_{i}$
terms for Ga-s(Al-s) and O-p, 2) the nearest-neighbor Ga-s(Al-s) to O-p and,
3) the dominant Ga-s(Al-s) to Ga-s(Al-s) terms ($t_{ij}$). This amounts to a
nearest-neighbor tight-binding model augmented by the Ga-s(Al-s) to Ga-s(Al-s)
next-nearest-neighbor hopping. These next-nearest-neighbor terms aid in the
accuracy of the higher Ga-4s(Al-3s) derived conduction bands. With these
simplifications, the model contains 60 parameters depending on the structural
phases (see Table SIII to SV). Including the Ga-p(Al-p) Wannier functions and
O-p to O-p coupling terms gives a satisfactory description of valence band DOS
at the cost of significantly increasing the number of parameters (around 300
terms), but provides little improvement to the description of the lower
conduction bands. Thus, we neglect these terms in our model.

To reproduce the experimental bandgaps and conduction band effective masses,
we adjust the Ga-s(Al-s) onsite energy and tune the Ga-s(Al-s) to O-p coupling
term. (We note that in their model of $\beta$-Ga2O3, Lee, et al.Lee  _et al._
(2019) focus on surface and defect formation energies which depend on the
broad features of the band structure. As a result, their tight-binding model
overestimates the conduction band effective mass.) The lists of parameters are
given in Table SIII, SIV, and SV.

When implementing the tight-binding model on a computer, Eq. (1) can be
written as a 22 $\times$ 22 matrix defined over the basis of 4 Ga-s(Al-s) and
18 O-p Wannier functions. In evaluating the model, the on-site energy
$\epsilon_{i}$ become diagonal terms while the hopping terms $t_{ij}$ populate
the ith row and jth column and must be multiplied by the phase factor
$e^{-i\vec{k}\cdot\vec{\Delta d_{ij}}}$ at each k-point. To illustrate the
model and gain insight into the material physics, we explicitly evaluate the
tight-binding model of $\beta$-Ga2O3 at the $\Gamma$-point
($\vec{k}=\vec{0}$). We find that the tight-binding model at the
$\Gamma$-point can be written in the block-matrix form:

$H=\begin{pmatrix}H_{Ga:s}&H_{Ga:s,O:p_{x}}^{\dagger}&0&H_{Ga:s,O:p_{z}}^{\dagger}\\\
H_{Ga:s,O:p_{x}}&H_{O:p_{x}}&0&0\\\ 0&0&H_{O:p_{y}}&0\\\
H_{Ga:s,O:p_{z}}&0&0&H_{O:p_{z}}\end{pmatrix}$ (2)

where the $H_{Ga:s}$, $H_{O:p_{x}}$, $H_{O:p_{y}}$, and $H_{O:p_{z}}$ blocks
define the coupling within the Ga-s, O-px, O-py, and O-pz Wannier function
sub-spaces, respectively. The $H_{Ga:s,O:p_{x}}$ and $H_{Ga:s,O:p_{z}}$ blocks
describe the coupling between the Ga-s and O-px, and O-pz Wannier functions,
respectively. We construct the tight-binding model so that it reflects the
crystal symmetry. This can be seen clearly in the lack of coupling between the
Ga-s and O-py Wannier functions. The two-fold rotation about the
crystallographic b-axis and the mirror operation through the plane
perpendicular to the b-axis guarantee this coupling is zero at the
$\Gamma$-point. The coupling between different O-p Wannier function blocks
(e.g. O-px and O-py) are zero because we have neglected the coupling within
the valence band. At the $\Gamma$ point, the phase factor
$e^{-i\vec{k}\cdot\vec{\Delta d_{ij}}}\rightarrow 1$, leaving the matrix real
and symmetric.

$H_{Ga:s}$ takes the form:

$\begin{split}H_{Ga:s}&=\begin{pmatrix}\epsilon_{Ga1:s}+2t_{Ga1:s,Ga1:s}&0&0&0\\\
0&\epsilon_{Ga1:s}+2t_{Ga1:s,Ga1:s}&0&0\\\
0&0&\epsilon_{Ga3:s}+2t_{Ga3:s,Ga3:s}&t_{Ga3:s,Ga4:s}\\\
0&0&t_{Ga3:s,Ga4:s}&\epsilon_{Ga3:s}+2t_{Ga3:s,Ga3:s}\end{pmatrix}\\\
&\rightarrow\begin{pmatrix}-7.498&0&0&0\\\ 0&-7.498&0&0\\\ 0&0&-7.222&0.216\\\
0&0&0.216&-7.222\end{pmatrix}\end{split}$ (3)

which is written in the basis
$\left(\left|Ga1:s\right>,\left|Ga2:s\right>,\left|Ga3:s\right>,\left|Ga4:s\right>\right)$
of Ga-s Wannier functions. Taking values from Table SIII, we have included the
numerical value of the block in eV. Since Ga1:s and Ga2:s are equivalent
tetrahedral sites while Ga3:s and Ga4:s are equivalent octahedral sites, each
have the same onsite energy and coupling. Notice that the onsite energy
$\epsilon_{Ga:1s}$ and $\epsilon_{Ga:3s}$ are modified by coupling to the same
Ga sites in the neighboring unit cells ($t_{Ga1:s,Ga1:s}$ and
$t_{Ga3:s,Ga3:s}$). Furthermore, the octahedral sites are weakly coupled to
each other.

In simplifying the description of the valence band, we neglect the coupling
between O-p Wannier functions. As a result each $H_{O:px}$, $H_{O:py}$, and
$H_{O:pz}$ block is diagonal and can be written as:

$H_{O:px}=H_{O:py}=H_{O:pz}=\epsilon_{O:p}\times\mathbb{I}_{6\times
6}\rightarrow-12\times\mathbb{I}_{6\times 6}$ (4)

which is written in the basis
$\left(\right.\left|O1:p_{l}\right>$,$\left|O2:p_{l}\right>$,$\left|O3:p_{l}\right>$,$\left|O4:p_{l}\right>$,$\left|O5:p_{l}\right>$,$\left|O6:p_{l}\right>\left.\right)$,
of O-p Wannier functions where $l=x,y,z$ . Here, $\mathbb{I}_{6\times 6}$ is
the 6 $\times$ 6 identity matrix.

The coupling between Ga-s and O-px and O-pz are described by $H_{Ga:s,O:px}$
and $H_{Ga:s,O:pz}$ which evaluate to (hopping strengths connected by
symmetries are shown with the same symbol):

$H_{Ga:s,O:px}=\begin{pmatrix}t_{Ga1:s,O1:px}&0&2t_{Ga3:s,O1:px}&0\\\
0&-t_{Ga1:s,O1:px}&0&2t_{Ga4:s,O2:px}\\\
t_{Ga1:s,O3:px}&0&0&-t_{Ga3:s,O4:px}\\\ 0&t_{Ga2:s,O4:px}&t_{Ga3:s,O4:px}&0\\\
0&-2t_{Ga1:s,O6:px}&t_{Ga3:s,O5:px}&0\\\
2t_{Ga1:s,O6:px}&0&0&-t_{Ga3:s,O5:px}\\\
\end{pmatrix}\rightarrow\begin{pmatrix}0.651&0&0.710&0\\\ 0&-0.651&0&-0.703\\\
2.592&0&0&-2.877\\\ 0&-2.595&2.877&0\\\ 0&3.499&-2.965&0\\\
-3.496&0&0&2.965\\\ \end{pmatrix}$ (5)
$H_{Ga:s,O:pz}=\begin{pmatrix}t_{Ga1:s,O1:pz}&0&2t_{Ga3:s,O1:pz}&0\\\
0&-t_{Ga1:s,O1:pz}&0&2t_{Ga4:s,O2:pz}\\\
t_{Ga1:s,O3:px}&0&2t_{Ga3:s,O3:pz}&-t_{Ga3:s,O4:pz}\\\
0&-t_{Ga2:s,O4:px}&t_{Ga3:s,O4:pz}&-2t_{Ga4:s,O4:pz}\\\
0&-2t_{Ga1:s,O6:pz}&t_{Ga3:s,O5:pz}&0\\\
2t_{Ga1:s,O6:pz}&0&0&-t_{Ga3:s,O5:pz}\\\
\end{pmatrix}\rightarrow\begin{pmatrix}3.466&0&-3.307&0\\\ 0&-3.466&0&3.314\\\
2.592&0&3.345&-1.101\\\ 0&2.595&1.101&-3.326\\\ 0&0.885&-0.621&0\\\
-0.880&0&0&0.621\\\ \end{pmatrix}$ (6)

Here, terms appear as pairs with opposite signs, manifesting the two-fold
symmetry of the monoclinic structure. $H_{Ga:s,O:px}$ and $H_{Ga:s,O:pz}$
looks formally similar but are not exactly the same. Again, the matrix
coupling Ga-s to O-py vanishes at the $\Gamma$ point but will show up away
from it. The complete set of matrix elements can be found in Table SIII. The
eigenvalues of the 22 $\times$ 22 matrix correspond to the energy eigenstates
at the given k-point, and the eigenvectors correspond to the wavefunctions. To
generate the band structure, one generates the matrix at sampled k-points
along the path and solve for the eigenvalues. The complexity of numerical
eigensolvers depends on the matrix dimension as $O(n^{3})$. As a comparison,
our DFT calculation relies on around 5000 planar waves as a basis, so the
tight-binding model provides a significant speedup. In our experience, a
tight-binding band structure can be generated on a personal computer in
seconds.

Figure 2:  Electronic structure of $\beta$-Ga2O3, $\alpha$-Ga2O3, and
$\alpha$-Al2O3. (a,b) – The first Brillouin zone of the monoclinic $\beta$ and
rhombohedral $\alpha$ phases. (c,f,i) – The DFT band structure and orbital-
projected DOS of $\beta$-Ga2O3, $\alpha$-Ga2O3, and $\alpha$-Al2O3. The DFT
bandgaps have been tuned to the experimental bandgaps via a scissor cut. Color
coordination indicates the orbital character of the bands and projected DOS.
The blue dots indicate the valence band maximum. (d,g,j) – The tight-binding
band structure and DOS (red) plotted over the DFT data (black). The sharp peak
in the tight-binding DOS at the top of the valence band is due to the lack of
O-p to O-p coupling in the tight-binding model. (e,h,k) – DFT (black) and
tight-binding (red) band structure near the $\Gamma$ point. The reciprocal
space vectors and high-symmetry points are given in Table SII.

##  IV Tight-Binding Bandstructure

Fig. 2(d,g,j) shows the tight-binding band structure and DOS superimposed on
the DFT results. In all three phases, the experimental bandgaps and conduction
band effective masses are reproduced by the adjustment of parameters mentioned
above. Moreover, the tight-binding model gives a satisfactory description of
the parabolic to linear dispersion, the second conduction band at the $\Gamma$
point, and the flat valence band edge states.

The slope of the linear dispersion of the lowest conduction band away from
$\Gamma$, and the energy of the second conduction band at $\Gamma$ are
slightly overestimated due to the absence of interaction from higher
Ga-p(Al-p) bands. This phenomenon is well-known in the tight-binding
description of Si, Ge, and III-V semiconductorsVogl, Hjalmarson, and Dow
(1983) where it is due to the lack of interactions from higher energy bands.
The lack of O-p to O-p coupling leaves some of the O-p states non-dispersive.
As a result, the tight-binding model cannot describe the small difference
between the valence band at $\Gamma$ and the DFT predicted valence band
maximum between $\Gamma$ and T in the $\beta$ phase and $\Gamma$ and $S_{0}$
in the $\alpha$ phases. This leaves a large peak in the DOS at the valence
band edge (see Fig. 2(d-f)) which could be resolved in future models that
include O-p to O-p coupling.

##  V Conclusion

We have derived a minimal tight-binding model for $\beta$-Ga2O3,
$\alpha$-Ga2O3, and $\alpha$-Al2O3 using selectively-localized Wannier
functions as a basis. The Wannier functions reflect the local symmetry of the
atomic sites and the tight-binding model satisfactorily reproduces the
electronic structure throughout the Brillouin zone. By constraining the
hopping parameters, we have fit the isotropic conduction band effective mass,
suggesting low field transport experiments can be described by the tight-
binding model. By including the higher energy Ga-s(Al-s) conduction bands, the
tight-binding model can also describe the parabolic to linear dispersion which
suggests application to high-field electronic transport as well. Finally, by
reproducing the experimental bandgap, we expect that the tight-binding model
can simulate the major features in optical absorption including the primary
optical transition from the O-p valence bands to the Ga-s(Al-s) conduction
bands, and the optical transition from the valence band to the higher energy
conduction bands. We expect the tight-binding models to aid in the description
and development of electronic and optical devices utilizing bulk,
nanostructured, heterostructured, and strained variants of $\beta$-Ga2O3,
$\alpha$-Ga2O3, and $\alpha$-Al2O3.

###### Acknowledgements.

This research project was conducted using computational resources at the
Maryland Advanced Research Computing Center (MARCC) and was supported by the
National Science Foundation (Platform for the Accelerated Realization,
Analysis and Discovery of Interface Materials (PARADIM)) under Cooperative
Agreement No. DMR-1539918.

## Data Availability Statement

The data that supports the findings of this study are available within the
article and its supplementary material.

## References

  * Higashiwaki  _et al._ (2012) M. Higashiwaki, K. Sasaki, A. Kuramata, T. Masui, and S. Yamakoshi, “Gallium oxide (Ga2O3) metal-semiconductor field-effect transistors on single-crystal $\beta$-Ga2O3 (010) substrates,” Applied Physics Letters 100, 013504 (2012).
  * Kuramata  _et al._ (2016) A. Kuramata, K. Koshi, S. Watanabe, Y. Yamaoka, T. Masui, and S. Yamakoshi, “High-quality $\beta$-Ga2O3 single crystals grown by edge-defined film-fed growth,” Japanese Journal of Applied Physics 55, 1202A2 (2016).
  * Gogova  _et al._ (2014) D. Gogova, G. Wagner, M. Baldini, M. Schmidbauer, K. Irmscher, R. Schewski, Z. Galazka, M. Albrecht, and R. Fornari, “Structural properties of Si-doped $\beta$-Ga2O3 layers grown by MOVPE,” Journal of Crystal Growth 401, 665–669 (2014).
  * Okumura  _et al._ (2014) H. Okumura, M. Kita, K. Sasaki, A. Kuramata, M. Higashiwaki, and J. S. Speck, “Systematic investigation of the growth rate of $\beta$-Ga2O3 (010) by plasma-assisted molecular beam epitaxy,” Applied Physics Express 7, 095501 (2014).
  * Jinno  _et al._ (2021) R. Jinno, C. S. Chang, T. Onuma, Y. Cho, S. T. Ho, D. Rowe, M. C. Cao, K. Lee, V. Protasenko, D. G. Schlom, D. A. Muller, H. G. Xing, and D. Jena, “Crystal orientation dictated epitaxy of ultrawide-bandgap 5.4-to 8.6-eV c-(AlGa)2O3 on m-plane sapphire,” Science Advances 7, 1–8 (2021), arXiv:2007.03415 .
  * Roy  _et al._ (2021) S. Roy, A. Bhattacharyya, P. Ranga, H. Splawn, J. Leach, and S. Krishnamoorthy, “High-k Oxide Field-Plated Vertical (001) $\beta$-Ga2O3 Schottky Barrier Diode with Baliga’s Figure of Merit Over 1 GW/cm2,” IEEE Electron Device Letters 3106, 1–1 (2021).
  * Kalarickal  _et al._ (2021) N. K. Kalarickal, Z. Xia, H.-l. Huang, W. Moore, Y. Liu, M. Brenner, J. Hwang, and S. Rajan, “$\beta$-(Al${}_{0}.18$Ga${}_{0}.82$)$2$O3 /Ga2O3 Double Heterojunction Transistor With Average Field of 5.5 MV/cm,” IEEE Electron Device Letters 42, 899–902 (2021).
  * Pearton  _et al._ (2018) S. J. Pearton, J. Yang, P. H. Cary, F. Ren, J. Kim, M. J. Tadjer, and M. A. Mastro, “A review of Ga2O3 materials, processing, and devices,” Applied Physics Reviews 5, 011301 (2018).
  * Meneghini, Meneghesso, and Zanoni (2017) M. Meneghini, G. Meneghesso, and E. Zanoni,  _Power GaN Devices: Materials, Applications and Reliability_ (Springer International Publishing, 2017).
  * Harrison (1989) W. A. Harrison,  _Electronic Structure and the Properties of Solids \- The Physics of the Chemical Bond_ (Dover Publications, 1989).
  * Datta (2005) S. Datta,  _Quantum Transport: Atom to Transistor_ (Cambridge University Press, 2005).
  * Lake  _et al._ (1997) R. Lake, G. Klimeck, R. C. Bowen, and D. Jovanovic, “Single and multiband modeling of quantum electron transport through layered semiconductor devices,” Journal of Applied Physics 81, 7845–7869 (1997).
  * Nehari  _et al._ (2005) K. Nehari, N. Cavassilas, J. L. Autran, M. Bescond, D. Munteanu, and M. Lannoo, “Influence of band-structure on electron ballistic transport in silicon nanowire MOSFET’s: An atomistic study,” Proceedings of ESSDERC 2005: 35th European Solid-State Device Research Conference 2005, 229–232 (2005).
  * Popov (2004) V. N. Popov, “Curvature effects on the structural, electronic and optical properties of isolated single-walled carbon nanotubes within a symmetry-adapted non-orthogonal tight-binding model,” New Journal of Physics 6 (2004), 10.1088/1367-2630/6/1/017.
  * Sols  _et al._ (1989) F. Sols, M. MacUcci, U. Ravaioli, and K. Hess, “Theory for a quantum modulated transistor,” Journal of Applied Physics 66, 3892–3906 (1989).
  * Graf and Vogl (1995) M. Graf and P. Vogl, “Electromagnetic fields and dielectric response in empirical tight-binding theory,” Physical Review B 51, 4940–4949 (1995).
  * Pi and Delerue (2013) X. Pi and C. Delerue, “Tight-binding calculations of the optical response of optimally p-doped Si nanocrystals: A model for localized surface plasmon resonance,” Physical Review Letters 111, 1–5 (2013).
  * Jirauschek and Kubis (2014) C. Jirauschek and T. Kubis, “Modeling techniques for quantum cascade lasers,” Applied Physics Reviews 1 (2014), 10.1063/1.4863665, arXiv:1412.3563 .
  * Lee  _et al._ (2019) J. Lee, S. Ganguli, A. K. Roy, and S. C. Badescu, “Density functional tight binding study of $\beta$-Ga2O3: Electronic structure, surface energy, and native point defects,” The Journal of Chemical Physics 150, 174706 (2019).
  * Souza, Marzari, and Vanderbilt (2002) I. Souza, N. Marzari, and D. Vanderbilt, “Maximally localized Wannier functions for entangled energy bands,” Physical Review B - Condensed Matter and Materials Physics 65, 1–13 (2002), arXiv:0108084 [cond-mat] .
  * Gresch  _et al._ (2018) D. Gresch, Q. Wu, G. W. Winkler, R. Häuselmann, M. Troyer, and A. A. Soluyanov, “Automated construction of symmetrized Wannier-like tight-binding models from ab initio calculations,” Physical Review Materials 2, 1–15 (2018), arXiv:1805.12148 .
  * Souza  _et al._ (2012) I. Souza, N. Marzari, D. Vanderbilt, M. Higashiwaki, K. Sasaki, A. Kuramata, T. Masui, S. Yamakoshi, J. Jung, A. H. Macdonald, H. Peelaers, C. G. Van de Walle, D. Gresch, Q. Wu, G. W. Winkler, R. Häuselmann, M. Troyer, A. A. Soluyanov, C. Franchini, R. Kováčik, M. Marsman, S. Sathyanarayana Murthy, J. He, C. Ederer, G. Kresse, S. Fang, R. Kuate Defo, S. N. Shirodkar, S. Lieu, G. A. Tritsaris, E. Kaxiras, S. Carr, S. Fang, H. C. Po, A. Vishwanath, and E. Kaxiras, “Ab initio tight-binding Hamiltonian for transition metal dichalcogenides,” Physical Review B - Condensed Matter and Materials Physics 1, 1–15 (2012), arXiv:1907.06282 .
  * Giannozzi  _et al._ (2009) P. Giannozzi, S. Baroni, N. Bonini, M. Calandra, R. Car, C. Cavazzoni, D. Ceresoli, G. L. Chiarotti, M. Cococcioni, I. Dabo, A. Dal Corso, S. de Gironcoli, S. Fabris, G. Fratesi, R. Gebauer, U. Gerstmann, C. Gougoussis, A. Kokalj, M. Lazzeri, L. Martin-Samos, N. Marzari, F. Mauri, R. Mazzarello, S. Paolini, A. Pasquarello, L. Paulatto, C. Sbraccia, S. Scandolo, G. Sclauzero, A. P. Seitsonen, A. Smogunov, P. Umari, and R. M. Wentzcovitch, “QUANTUM ESPRESSO: a modular and open-source software project for quantum simulations of materials,” Journal of Physics: Condensed Matter 21, 395502 (2009).
  * Perdew, Burke, and Ernzerhof (1996) J. P. Perdew, K. Burke, and M. Ernzerhof, “Generalized gradient approximation made simple,” Physical Review Letters 77, 3865–3868 (1996).
  * 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,” Physical Review Letters 100, 1–4 (2008), arXiv:0711.0156 .
  * Zade  _et al._ (2019) V. Zade, B. Mallesham, S. Shantha-Kumar, A. Bronson, and C. V. Ramana, “Interplay between Solubility Limit, Structure, and Optical Properties of Tungsten-Doped Ga 2 O 3 Compounds Synthesized by a Two-Step Calcination Process,” Inorganic Chemistry 58, 3707–3716 (2019).
  * Marezio and Remeika (1967) M. Marezio and J. P. Remeika, “Bond Lengths in the $\alpha$-Ga2O3 Structure and the High-Pressure Phase of Ga2-xFexO3,” The Journal of Chemical Physics 46, 1862–1865 (1967).
  * Ishizawa  _et al._ (1980) N. Ishizawa, T. Miyata, I. Minato, F. Marumo, and S. Iwai, “A structural investigation of $\alpha$-Al 2 O 3 at 2170 K,” Acta Crystallographica Section B Structural Crystallography and Crystal Chemistry 36, 228–230 (1980).
  * Peelaers and Van de Walle (2015) H. Peelaers and C. G. Van de Walle, “Brillouin zone and band structure of $\beta$-Ga 2 O 3,” physica status solidi (b) 252, 828–832 (2015).
  * Zhang  _et al._ (2018) Y. Zhang, A. Neal, Z. Xia, C. Joishi, J. M. Johnson, Y. Zheng, S. Bajaj, M. Brenner, D. Dorsey, K. Chabak, G. Jessen, J. Hwang, S. Mou, J. P. Heremans, and S. Rajan, “Demonstration of high mobility and quantum transport in modulation-doped $\beta$-(AlxGa1-x)2O3/Ga2O3 heterostructures,” Applied Physics Letters 112, 1–6 (2018).
  * Mohamed  _et al._ (2010) M. Mohamed, C. Janowitz, I. Unger, R. Manzke, Z. Galazka, R. Uecker, R. Fornari, J. R. Weber, J. B. Varley, and C. G. Van De Walle, “The electronic structure of $\beta$-Ga2O3,” Applied Physics Letters 97 (2010), 10.1063/1.3521255.
  * Janowitz  _et al._ (2011) C. Janowitz, V. Scherer, M. Mohamed, A. Krapf, H. Dwelk, R. Manzke, Z. Galazka, R. Uecker, K. Irmscher, R. Fornari, M. Michling, D. Schmeißer, J. R. Weber, J. B. Varley, and C. G. Van De Walle, “Experimental electronic structure of In2O3 and Ga2O3,” New Journal of Physics 13 (2011), 10.1088/1367-2630/13/8/085014.
  * Ghosh and Singisetti (2017) K. Ghosh and U. Singisetti, “Ab initio velocity-field curves in monoclinic $\beta$-Ga2O3,” Journal of Applied Physics 122 (2017), 10.1063/1.4986174.
  * Ghosh and Singisetti (2018) K. Ghosh and U. Singisetti, “Impact ionization in $\beta$-Ga2O3,” Journal of Applied Physics 124 (2018), 10.1063/1.5034120.
  * Singh  _et al._ (2020) A. Singh, O. Koksal, N. Tanen, J. McCandless, D. Jena, H. G. Xing, H. Peelaers, and F. Rana, “Intra- and inter-conduction band optical absorption processes in $\beta$-Ga2O3,” Applied Physics Letters 117, 072103 (2020), arXiv:2006.14850 .
  * Shinohara and Fujita (2008) D. Shinohara and S. Fujita, “Heteroepitaxy of corundum-structured $\alpha$-Ga2O3 thin films on $\alpha$-Ga2O3 substrates by ultrasonic mist chemical vapor deposition,” Japanese Journal of Applied Physics 47, 7311–7313 (2008).
  * French (1990) R. H. French, “Electronic Band Structure of Al2O3, with Comparison to Alon and AIN,” Journal of the American Ceramic Society 73, 477–489 (1990).
  * Pizzi  _et al._ (2020) G. Pizzi, V. Vitale, R. Arita, S. Blügel, F. Freimuth, G. Géranton, M. Gibertini, D. Gresch, C. Johnson, T. Koretsune, J. Ibañez-Azpiroz, H. Lee, J.-M. Lihm, D. Marchand, A. Marrazzo, Y. Mokrousov, J. I. Mustafa, Y. Nohara, Y. Nomura, L. Paulatto, S. Poncé, T. Ponweiser, J. Qiao, F. Thöle, S. S. Tsirkin, M. Wierzbowska, N. Marzari, D. Vanderbilt, I. Souza, A. A. Mostofi, and J. R. Yates, “Wannier90 as a community code: new features and applications,” Journal of Physics: Condensed Matter 32, 165902 (2020).
  * Wang  _et al._ (2014) R. Wang, E. A. Lazar, H. Park, A. J. Millis, and C. A. Marianetti, “Selectively localized Wannier functions,” Physical Review B - Condensed Matter and Materials Physics 90 (2014), 10.1103/PhysRevB.90.165125, arXiv:1407.5124 .
  * Coh and Vanderbilt (2013) S. Coh and D. Vanderbilt, “Python tight binding (PythTB),” (2013).
  * Vogl, Hjalmarson, and Dow (1983) P. Vogl, H. P. Hjalmarson, and J. D. Dow, “A Semi-empirical tight-binding theory of the electronic structure of semiconductors†,” Journal of Physics and Chemistry of Solids 44, 365–378 (1983).
  * Tsao  _et al._ (2018) J. Y. Tsao, S. Chowdhury, M. A. Hollis, D. Jena, N. M. Johnson, K. A. Jones, R. J. Kaplar, S. Rajan, C. G. Van de Walle, E. Bellotti, C. L. Chua, R. Collazo, M. E. Coltrin, J. A. Cooper, K. R. Evans, S. Graham, T. A. Grotjohn, E. R. Heller, M. Higashiwaki, M. S. Islam, P. W. Juodawlkis, M. A. Khan, A. D. Koehler, J. H. Leach, U. K. Mishra, R. J. Nemanich, R. C. Pilawa-Podgurski, J. B. Shealy, Z. Sitar, M. J. Tadjer, A. F. Witulski, M. Wraback, and J. A. Simmons, “Ultrawide-Bandgap Semiconductors: Research Opportunities and Challenges,” Advanced Electronic Materials 4 (2018), 10.1002/aelm.201600501.
  * Montgomery, Todorov, and Todorov (2002) M. J. Montgomery, T. N. Todorov, and T. N. Todorov, “Tight-binding simulation of current-carrying nanostructures,” (2002).
  * Paskaleva, Spassov, and Terziyska (2017) A. Paskaleva, D. Spassov, and P. Terziyska, “Electric, dielectric and optical properties of Ga2O3 grown by metal organic chemical vapour deposition,” Journal of Physics: Conference Series 794, 012017 (2017).
  * Ruch and Kino (1968) J. G. Ruch and G. S. Kino, “Transport Properties of GaAs,” Physical Review 174, 921–931 (1968).

Supplemental Information for Tight-Binding Bandstructure of $\beta-$ and
$\alpha-$ phase Ga2O3 and Al2O3

Y. Zhang1 M. Liu2 D. Jena1,2 and G. Khalsa2

1\. Department of Electrical and Computer Engineering and 2. Department of
Materials Science and Engineering, Cornell University, Ithaca NY 14853, USA

In the Supplementary Tables we include additional information for the tight-
binding models of the monoclinic $\beta$-Ga2O3, and the rhombohedral
$\alpha$-Ga2O3 and $\alpha$-Al2O3 phases. Table SI gives the structural
information from density functional theory in scaled coordinates. Table SII
shows the high-symmetry points in reciprocal space for the monoclinic and
rhombohedral structures. Table SIII-SV give the explicit tight-binding models
for each phase.

TABLE SI: Structural parameters of $\beta$-Ga2O3, $\alpha$-Ga2O3, and
$\alpha$-Al2O3. The Cartesian coordinate of a site $i$ is given by
$\vec{d}_{i}=A_{1}\vec{a_{1}}+A_{2}\vec{a_{2}}+A_{3}\vec{a_{3}}$

lattice |  $\beta$-Ga2O3 |
---|---|---
vector (Å) | $x$ | $y$ | $z$ |
$\vec{a_{1}}$ | 6.115 | 1.620 | 0.000 |
$\vec{a_{2}}$ | -6.115 | 1.620 | 0.000 |
$\vec{a_{3}}$ | -1.374 | 0.000 | 5.635 |
site |  |  |  |
coordinate | $A_{1}$ | $A_{2}$ | $A_{3}$ |
Ga1 | 0.090 | -0.090 | 0.795 |
Ga2 | 0.910 | 0.090 | 0.205 |
Ga3 | 0.659 | 0.341 | 0.314 |
Ga4 | 0.341 | -0.341 | 0.686 |
O1 | 0.165 | -0.165 | 0.109 |
O2 | 0.835 | 0.165 | 0.891 |
O3 | 0.173 | -0.173 | 0.563 |
O4 | 0.827 | 0.173 | 0.437 |
O5 | 0.496 | -0.496 | 0.256 |
O6 | 0.504 | 0.496 | 0.743 |

lattice |  $\alpha$-Ga2O3 |
---|---|---
vector (Å) | $x$ | $y$ | $z$ |
$\vec{a_{1}}$ | 2.491 | 1.438 | 4.478 |
$\vec{a_{2}}$ | -2.491 | 1.438 | 4.478 |
$\vec{a_{3}}$ | 0.000 | -2.877 | 4.478 |
site |  |  |  |
coordinate | $A_{1}$ | $A_{2}$ | $A_{3}$ |
Ga1 | 0.179 | 0.179 | 0.179 |
Ga2 | 0.678 | 0.678 | -0.322 |
Ga3 | 0.821 | -0.179 | -0.179 |
Ga4 | 0.322 | 0.322 | -0.678 |
O1 | 0.553 | -0.053 | 0.250 |
O2 | 0.941 | 0.250 | -0.441 |
O3 | 0.250 | 0.559 | -0.059 |
O4 | 1.059 | -0.250 | -0.559 |
O5 | 0.750 | -0.559 | 0.059 |
O6 | 0.441 | 0.059 | -0.250 |

lattice |  $\alpha$-Al2O3 |
---|---|---
vector (Å) | $x$ | $y$ | $z$ |
$\vec{a_{1}}$ | 2.403 | 1.387 | 4.372 |
$\vec{a_{2}}$ | -2.403 | 1.387 | 4.372 |
$\vec{a_{3}}$ | 0.000 | -2.774 | 44.372 |
site |  |  |  |
coordinate | $A_{1}$ | $A_{2}$ | $A_{3}$ |
Al1 | 0.818 | -0.162 | -0.177 |
Al2 | 0.671 | 0.659 | -0.311 |
Al3 | 0.330 | 0.318 | -0.664 |
Al4 | 0.174 | 0.182 | 0.164 |
O1 | 0.750 | -0.558 | 0.057 |
O2 | 0.942 | 0.250 | -0.442 |
O3 | 1.058 | -0.250 | -0.558 |
O4 | 0.442 | 0.058 | -0.250 |
O5 | 0.558 | -0.058 | 0.250 |
O6 | 0.250 | 0.558 | -0.058 |

TABLE SII: High-symmetry points in the first Brilloin zone of the $\beta$ and
$\alpha$ phase. Coordinates are calculated by
$\vec{k}=B_{1}\vec{b_{1}}+B_{2}\vec{b_{2}}+B_{3}\vec{b_{3}}$, in which
$\vec{b_{1}}$, $\vec{b_{2}}$, and $\vec{b_{3}}$ are the reciprocal lattice
vector. $\vec{b_{1}}$, $\vec{b_{2}}$, and $\vec{b_{3}}$ can be calculated by
$\begin{pmatrix}\vec{b_{1}}&\vec{b_{2}}&\vec{b_{3}}\end{pmatrix}=(\begin{pmatrix}\vec{a_{1}}&\vec{a_{2}}&\vec{a_{3}}\end{pmatrix}^{-1})^{T}$

$\beta$ phase
---
| $B_{1}$ | $B_{2}$ | $B_{3}$
$\Gamma$ | 0.0000 | 0.0000 | 0.0000
$C$ | 0.2662 | 0.2662 | 0.0000
$C_{2}$ | -0.2662 | 0.7338 | 0.0000
$Y_{2}$ | -0.5000 | 0.5000 | 0.0000
$M_{2}$ | -0.5000 | 0.5000 | 0.5000
$D$ | -0.2580 | 0.7419 | 0.5000
$D_{2}$ | 0.2580 | 0.2580 | 0.0000
$A$ | 0.0000 | 0.0000 | 0.5000
$L_{2}$ | 0.0000 | 0.5000 | 0.5000
$V_{2}$ | 0.0000 | 0.5000 | 0.0000

$\alpha$ phase
---
| $B_{1}$ | $B_{2}$ | $B_{3}$
$\Gamma$ | 0.0000 | 0.0000 | 0.0000
$T$ | 0.5000 | 0.5000 | 0.5000
$H_{2}$ | 0.7641 | 0.2358 | 0.5000
$H_{0}$ | 0.5000 | -0.2358 | 0.2358
$L$ | 0.5000 | 0.0000 | 0.0000
$S_{0}$ | 0.3679 | -0.3679 | 0.0000
$S_{2}$ | 0.6320 | 0.0000 | 0.3679
$F$ | 0.5000 | 0.0000 | 0.5000

TABLE SIII: Tight-binding parameters of $\beta$-Ga2O3. $\vec{R}$ is written in
integer multiples of lattice vectors $[ijk]$ and can be translated to
Cartesian coordinates by $\vec{R}=i\vec{a_{1}}+j\vec{a_{2}}+k\vec{a_{3}}$,
where $a_{1}$, $a_{2}$, and $a_{3}$ are the lattice vectors.

onsite energy |  |  |
---|---|---|---
site $i$ |  $\epsilon_{i}$ (eV) |  |  |
Ga1, Ga2 | 4.95 |  |  |
Ga3, Ga4 | 4.52 |  |  |
all O | 0 |  |  |
hopping parameters |  |  |
site $i$ | site $j$ | $\vec{R}$ |  $\tilde{t}_{ij}$ (eV) | $\phi_{ij}$
Ga1 s | Ga1 s | [1,1,0] | 0.224 | $\pi$
Ga2 s | Ga2 s | [1,1,0] | 0.224 | $\pi$
Ga3 s | Ga3 s | [1,1,0] | 0.129 | 0
Ga3 s | Ga4 s | [0,0,0] | 0.108 | 0
Ga3 s | Ga4 s | [1,1,0] | 0.108 | 0
Ga4 s | Ga4 s | [1,1,0] | 0.129 | 0
Ga1 s | O1 px | [0,0,1] | 0.652 | 0
Ga1 s | O1 pz | [0,0,1] | 3.467 | 0
Ga1 s | O3 px | [0,0,0] | 2.592 | 0
Ga1 s | O3 pz | [0,0,0] | 2.592 | $\pi$
Ga1 s | O6 px | [0,0,0] | 1.748 | $\pi$
Ga1 s | O6 px | [-1,-1,0] | 1.748 | $\pi$
Ga1 s | O6 py | [0,0,0] | 3.029 | 0
Ga1 s | O6 py | [-1,-1,0] | 3.029 | $\pi$
Ga1 s | O6 pz | [0,0,0] | 0.440 | $\pi$
Ga1 s | O6 pz | [-1,-1,0] | 0.440 | $\pi$
Ga2 s | O2 px | [0,0,-1] | 0.652 | $\pi$
Ga2 s | O2 pz | [0,0,-1] | 3.467 | $\pi$
Ga2 s | O4 px | [0,0,0] | 2.596 | $\pi$
Ga2 s | O4 pz | [0,0,0] | 2.596 | 0
Ga2 s | O5 px | [0,0,0] | 1.750 | 0
Ga2 s | O5 px | [1,1,0] | 1.750 | 0
Ga2 s | O5 py | [0,0,0] | 3.030 | $\pi$
Ga2 s | O5 py | [1,1,0] | 3.030 | 0
Ga2 s | O5 pz | [0,0,0] | 0.443 | 0
Ga2 s | O5 pz | [1,1,0] | 0.443 | 0

hopping parameters continued |  |
---|---|---
site $i$ | site $j$ | $\vec{R}$ |  $\tilde{t}_{ij}$ (eV) | $\phi_{ij}$
Ga3 s | O1 px | [0,0,0] | 0.355 | 0
Ga3 s | O1 px | [1,1,0] | 0.355 | 0
Ga3 s | O1 py | [0,0,0] | 2.484 | $\pi$
Ga3 s | O1 py | [1,1,0] | 2.484 | 0
Ga3 s | O1 pz | [0,0,0] | 1.654 | $\pi$
Ga3 s | O1 pz | [1,1,0] | 1.654 | $\pi$
Ga3 s | O3 py | [0,0,0] | 2.162 | $\pi$
Ga3 s | O3 py | [1,1,0] | 2.162 | 0
Ga3 s | O3 pz | [0,0,0] | 1.673 | 0
Ga3 s | O3 pz | [1,1,0] | 1.673 | 0
Ga3 s | O4 px | [0,0,0] | 2.878 | 0
Ga3 s | O4 pz | [0,0,0] | 1.102 | 0
Ga3 s | O5 px | [0,1,0] | 2.965 | $\pi$
Ga3 s | O5 pz | [0,1,0] | 0.622 | $\pi$
Ga4 s | O2 px | [0,0,0] | 0.352 | $\pi$
Ga4 s | O2 px | [-1,-1,0] | 0.352 | $\pi$
Ga4 s | O2 py | [0,0,0] | 2.484 | 0
Ga4 s | O2 py | [-1,-1,0] | 2.484 | $\pi$
Ga4 s | O2 pz | [0,0,0] | 1.657 | 0
Ga4 s | O2 pz | [-1,-1,0] | 1.657 | 0
Ga4 s | O3 px | [0,0,0] | 2.878 | $\pi$
Ga4 s | O3 pz | [0,0,0] | 1.102 | $\pi$
Ga4 s | O4 py | [0,0,0] | 2.164 | 0
Ga4 s | O4 py | [-1,-1,0] | 2.164 | $\pi$
Ga4 s | O4 pz | [0,0,0] | 1.663 | $\pi$
Ga4 s | O4 pz | [-1,-1,0] | 1.663 | $\pi$
Ga4 s | O6 px | [0,-1,0] | 2.965 | 0
Ga4 s | O6 pz | [0,-1,0] | 0.622 | 0

TABLE SIV: Tight-binding parameters of $\alpha$-Ga2O3. $\vec{R}$ is written in
integer multiples of lattice vectors $[ijk]$ and can be translated to
Cartesian coordinates by $\vec{R}=i\vec{a_{1}}+j\vec{a_{2}}+k\vec{a_{3}}$,
where $a_{1}$, $a_{2}$, and $a_{3}$ are the lattice vectors.

onsite energy |  |  |
---|---|---|---
site $i$ |  $\epsilon_{i}$ (eV) |  |  |
all Ga | -5.48 |  |  |
all O | -10.5 |  |  |
hopping parameters |  |  |
site $i$ | site $j$ | $\vec{R}$ | $\tilde{t}_{ij}$ | $\phi_{ij}$
Ga1 s | Ga3 s | [0,0,0] | 0.012 | 0
Ga1 s | Ga3 s | [-1,1,0] | 0.013 | 0
Ga1 s | Ga3 s | [-1,0,1] | 0.013 | 0
Ga1 s | Ga4 s | [0,0,1] | 0.014 | $\pi$
Ga1 s | Ga4 s | [0,-1,1] | 0.010 | $\pi$
Ga1 s | Ga4 s | [-1,0,1] | 0.010 | $\pi$
Ga2 s | Ga3 s | [0,0,0] | 0.010 | $\pi$
Ga2 s | Ga3 s | [0,1,0] | 0.014 | $\pi$
Ga2 s | Ga3 s | [-1,1,0] | 0.010 | $\pi$
Ga2 s | Ga4 s | [0,0,1] | 0.013 | 0
Ga2 s | Ga4 s | [0,1,0] | 0.013 | 0
Ga2 s | Ga4 s | [1,0,0] | 0.013 | 0
Ga1 s | O1 px | [0,0,0] | 1.983 | 0
Ga1 s | O1 py | [0,0,0] | 0.105 | 0
Ga1 s | O1 pz | [0,0,0] | 1.860 | 0
Ga1 s | O2 px | [-1,0,1] | 0.910 | $\pi$
Ga1 s | O2 py | [-1,0,1] | 1.766 | $\pi$
Ga1 s | O2 pz | [-1,0,1] | 1.864 | 0
Ga1 s | O3 px | [0,0,0] | 1.076 | $\pi$
Ga1 s | O3 py | [0,0,0] | 1.673 | 0
Ga1 s | O3 pz | [0,0,0] | 1.867 | 0
Ga1 s | O4 px | [-1,0,1] | 1.151 | 0
Ga1 s | O4 py | [-1,0,1] | 2.402 | $\pi$
Ga1 s | O4 pz | [-1,0,1] | 1.248 | $\pi$
Ga1 s | O5 px | [-1,1,0] | 2.654 | $\pi$
Ga1 s | O5 py | [-1,1,0] | 0.216 | 0
Ga1 s | O5 pz | [-1,1,0] | 1.246 | $\pi$
Ga1 s | O6 px | [0,0,0] | 1.509 | 0
Ga1 s | O6 py | [0,0,0] | 2.197 | 0
Ga1 s | O6 pz | [0,0,0] | 1.242 | $\pi$
Ga2 s | O1 px | [0,1,-1] | 1.511 | $\pi$
Ga2 s | O1 py | [0,1,-1] | 2.198 | 0
Ga2 s | O1 pz | [0,1,-1] | 1.249 | $\pi$
Ga2 s | O2 px | [0,0,0] | 2.656 | 0
Ga2 s | O2 py | [0,0,0] | 0.216 | 0
Ga2 s | O2 pz | [0,0,0] | 1.242 | $\pi$
Ga2 s | O3 px | [0,0,0] | 1.142 | $\pi$
Ga2 s | O3 py | [0,0,0] | 2.404 | $\pi$
Ga2 s | O3 pz | [0,0,0] | 1.245 | $\pi$
Ga2 s | O4 px | [0,1,0] | 1.074 | 0
Ga2 s | O4 py | [0,1,0] | 1.664 | 0
Ga2 s | O4 pz | [0,1,0] | 1.866 | 0
Ga2 s | O5 px | [0,1,0] | 0.904 | 0
Ga2 s | O5 py | [0,1,0] | 1.759 | $\pi$
Ga2 s | O5 pz | [0,1,0] | 1.872 | 0
Ga2 s | O6 px | [0,1,0] | 1.980 | $\pi$
Ga2 s | O6 py | [0,1,0] | 0.096 | 0
Ga2 s | O6 pz | [0,1,0] | 1.870 | 0

hopping parameters continued |  |
---|---|---
site $i$ | site $j$ | $\vec{R}$ |  $\tilde{t}_{ij}$ (eV) | $\phi_{ij}$
Ga3 s | O1 px | [0,0,0] | 1.508 | $\pi$
Ga3 s | O1 py | [0,0,0] | 2.191 | $\pi$
Ga3 s | O1 pz | [0,0,0] | 1.252 | 0
Ga3 s | O2 px | [0,0,0] | 1.151 | $\pi$
Ga3 s | O2 py | [0,0,0] | 2.402 | 0
Ga3 s | O2 pz | [0,0,0] | 1.251 | 0
Ga3 s | O3 px | [1,-1,0] | 2.656 | 0
Ga3 s | O3 py | [1,-1,0] | 0.199 | $\pi$
Ga3 s | O3 pz | [1,-1,0] | 1.253 | 0
Ga3 s | O4 px | [0,0,0] | 0.910 | 0
Ga3 s | O4 py | [0,0,0] | 1.779 | 0
Ga3 s | O4 pz | [0,0,0] | 1.851 | $\pi$
Ga3 s | O5 px | [0,0,0] | 1.084 | 0
Ga3 s | O5 py | [0,0,0] | 1.681 | $\pi$
Ga3 s | O5 pz | [0,0,0] | 1.853 | $\pi$
Ga3 s | O6 px | [0,0,0] | 1.995 | $\pi$
Ga3 s | O6 py | [0,0,0] | 0.098 | $\pi$
Ga3 s | O6 pz | [0,0,0] | 1.856 | $\pi$
Ga4 s | O1 px | [0,0,-1] | 1.989 | 0
Ga4 s | O1 py | [0,0,-1] | 0.098 | $\pi$
Ga4 s | O1 pz | [0,0,-1] | 1.859 | $\pi$
Ga4 s | O2 px | [-1,0,0] | 1.077 | $\pi$
Ga4 s | O2 py | [-1,0,0] | 1.677 | $\pi$
Ga4 s | O2 pz | [-1,0,0] | 1.862 | $\pi$
Ga4 s | O3 px | [0,0,-1] | 0.908 | $\pi$
Ga4 s | O3 py | [0,0,-1] | 1.775 | 0
Ga4 s | O3 pz | [0,0,-1] | 1.856 | $\pi$
Ga4 s | O4 px | [-1,1,0] | 2.653 | $\pi$
Ga4 s | O4 py | [-1,1,0] | 0.213 | $\pi$
Ga4 s | O4 pz | [-1,1,0] | 1.253 | 0
Ga4 s | O5 px | [0,1,-1] | 1.151 | 0
Ga4 s | O5 py | [0,1,-1] | 2.405 | 0
Ga4 s | O5 pz | [0,1,-1] | 1.249 | 0
Ga4 s | O6 px | [0,0,0] | 1.510 | 0
Ga4 s | O6 py | [0,0,0] | 2.196 | $\pi$
Ga4 s | O6 pz | [0,0,0] | 1.254 | 0

TABLE SV: Tight-binding parameters of $\alpha$-Al2O3. $\vec{R}$ is written in
integer multiples of lattice vectors $[ijk]$ and can be translated to
Cartesian coordinates by $\vec{R}=i\vec{a_{1}}+j\vec{a_{2}}+k\vec{a_{3}}$,
where $a_{1}$, $a_{2}$, and $a_{3}$ are the lattice vectors.

onsite energy (eV) |  |  |
---|---|---|---
all Al | 7.00 |  |  |
all O | 0.00 |  |  |
hopping parameters |  |  |
site $i$ | site $j$ | $\vec{R}$ | $\tilde{t}_{ij}$ | $\phi_{ij}$
Al1 s | Al2 s | [0,-1,0] | 1.019 | $\pi$
Al1 s | Al3 s | [0,0,0] | 0.429 | 0
Al1 s | Al3 s | [0,0,1] | 0.389 | 0
Al1 s | Al3 s | [1,0,0] | 0.405 | 0
Al2 s | Al4 s | [0,0,0] | 0.406 | 0
Al2 s | Al4 s | [1,0,0] | 0.415 | 0
Al3 s | Al4 s | [0,0,-1] | 1.024 | $\pi$
Al1 s | O1 px | [0,0,0] | 1.592 | 0
Al1 s | O1 py | [0,0,0] | 2.614 | $\pi$
Al1 s | O1 pz | [0,0,0] | 2.646 | $\pi$
Al1 s | O2 px | [0,0,0] | 1.008 | $\pi$
Al1 s | O2 py | [0,0,0] | 2.321 | 0
Al1 s | O2 pz | [0,0,0] | 1.738 | 0
Al1 s | O3 px | [0,0,0] | 1.553 | 0
Al1 s | O3 py | [0,0,0] | 2.708 | 0
Al1 s | O3 pz | [0,0,0] | 2.685 | $\pi$
Al1 s | O4 px | [0,0,0] | 3.249 | $\pi$
Al1 s | O4 pz | [0,0,0] | 2.871 | $\pi$
Al1 s | O5 px | [0,0,0] | 1.392 | $\pi$
Al1 s | O5 py | [0,0,0] | 2.058 | $\pi$
Al1 s | O5 pz | [0,0,0] | 1.717 | 0
Al1 s | O6 px | [1,-1,0] | 2.352 | 0
Al1 s | O6 py | [1,-1,0] | 0.250 | $\pi$
Al1 s | O6 pz | [1,-1,0] | 1.597 | 0
Al2 s | O1 px | [0,1,0] | 1.578 | 0
Al2 s | O1 py | [0,1,0] | 2.783 | $\pi$
Al2 s | O1 pz | [0,1,0] | 2.994 | 0
Al2 s | O2 px | [0,0,0] | 2.511 | 0
Al2 s | O2 py | [0,0,0] | 0.365 | 0
Al2 s | O2 pz | [0,0,0] | 1.704 | $\pi$
Al2 s | O3 px | [0,1,0] | 1.473 | 0
Al2 s | O3 py | [0,1,0] | 2.567 | 0
Al2 s | O3 pz | [0,1,0] | 2.657 | 0
Al2 s | O4 px | [0,1,0] | 2.998 | $\pi$
Al2 s | O4 pz | [0,1,0] | 2.722 | 0
Al2 s | O5 px | [0,1,-1] | 1.381 | $\pi$
Al2 s | O5 py | [0,1,-1] | 1.901 | 0
Al2 s | O5 pz | [0,1,-1] | 1.554 | $\pi$
Al2 s | O6 px | [0,0,0] | 1.144 | $\pi$
Al2 s | O6 py | [0,0,0] | 2.297 | $\pi$
Al2 s | O6 pz | [0,0,0] | 1.769 | $\pi$

hopping parameters continued |  |
---|---|---
site $i$ | site $j$ | $\vec{R}$ |  $\tilde{t}_{ij}$ (eV) | $\phi_{ij}$
Al3 s | O1 px | [0,1,-1] | 0.924 | 0
Al3 s | O1 py | [0,1,-1] | 2.264 | 0
Al3 s | O1 pz | [0,1,-1] | 1.620 | 0
Al3 s | O2 px | [-1,0,0] | 1.644 | $\pi$
Al3 s | O2 py | [-1,0,0] | 2.699 | $\pi$
Al3 s | O2 pz | [-1,0,0] | 2.806 | $\pi$
Al3 s | O3 px | [-1,1,0] | 2.385 | $\pi$
Al3 s | O3 py | [-1,1,0] | 0.189 | $\pi$
Al3 s | O3 pz | [-1,1,0] | 1.636 | 0
Al3 s | O4 px | [0,0,0] | 1.476 | 0
Al3 s | O4 py | [0,0,0] | 2.108 | $\pi$
Al3 s | O4 pz | [0,0,0] | 1.798 | 0
Al3 s | O5 px | [0,0,-1] | 3.176 | 0
Al3 s | O5 pz | [0,0,-1] | 2.857 | $\pi$
Al3 s | O6 px | [0,0,-1] | 1.451 | $\pi$
Al3 s | O6 py | [0,0,-1] | 2.631 | 0
Al3 s | O6 pz | [0,0,-1] | 2.596 | $\pi$
Al4 s | O1 px | [-1,1,0] | 2.510 | $\pi$
Al4 s | O1 py | [-1,1,0] | 0.197 | 0
Al4 s | O1 pz | [-1,1,0] | 1.727 | $\pi$
Al4 s | O2 px | [-1,0,1] | 1.468 | $\pi$
Al4 s | O2 py | [-1,0,1] | 2.660 | $\pi$
Al4 s | O2 pz | [-1,0,1] | 2.715 | 0
Al4 s | O3 px | [-1,0,1] | 0.974 | 0
Al4 s | O3 py | [-1,0,1] | 2.175 | $\pi$
Al4 s | O3 pz | [-1,0,1] | 1.606 | $\pi$
Al4 s | O4 px | [0,0,0] | 1.535 | 0
Al4 s | O4 py | [0,0,0] | 2.022 | 0
Al4 s | O4 pz | [0,0,0] | 1.687 | $\pi$
Al4 s | O5 px | [0,0,0] | 3.038 | 0
Al4 s | O5 pz | [0,0,0] | 2.702 | 0
Al4 s | O6 px | [0,0,0] | 1.619 | $\pi$
Al4 s | O6 py | [0,0,0] | 2.746 | 0
Al4 s | O6 pz | [0,0,0] | 2.903 | 0


# Rejuvenation engineering in metallic glasses by complementary stress and
structure modulation

D. Şopu1,2<EMAIL_ADDRESS>F .Spieckermann3
<EMAIL_ADDRESS>X.L. Bian4<EMAIL_ADDRESS>S.
Fellner1 J. Wright5 M. Cordill1 C. Gammer1 G. Wang4 M. Stoica5 J. Eckert1,3
1Erich Schmid Institute of Materials Science, Austrian Academy of Sciences,
Jahnstraße 12, A-8700 Leoben, Austria
2Institut für Materialwissenschaft, Fachgebiet Materialmodellierung,
Technische Universität Darmstadt, Otto-Berndt-Straße 3, D-64287 Darmstadt,
Germany
3Department of Materials Science, Chair of Materials Physics,
Mountanuniversität Leoben, Jahnstraße 12, A-8700 Leoben, Austria
4Laboratory for Microstructures, Institute of Materials, Shanghai University,
200444, Shanghai, China
5European Synchrotron Radiation Facility (ESRF), 38042 Grenoble, France
6Laboratory of Metal Physics and Technology, Department of Materials, ETH
Zurich, 8093 Zurich, Switzerland

###### Abstract

Residual stress engineering is very widely used in the design of new advanced
lightweight materials. For metallic glasses the attention has been on
structural changes and rejuvenation processes. High energy scanning X-ray
diffraction strain mapping reveals large elastic fluctuations in metallic
glasses after deformation under triaxial compression. Microindentation
hardness mapping hints to a competing hardening-softening mechanism after
compression and further reveals the complementary effects of stress and
structure modulation. Transmission electron microscopy proves that structure
modulation under room temperature deformation relates to the shear band
formation that closely correlates to the distribution of elastic
heterogeneities. Molecular dynamics simulations provide an atomistic
understanding of the complex shear band activity in notched metallic glasses
and the related fluctuations in the strain/stress heterogeneity. Thus, future
focus should be given to stress engineering and elastic heterogeneity that
together with structure modulation may allow to design metallic glasses with
enhanced ductility and strain hardening ability.

##  I Introduction

Given the metastable nature of metallic glasses (MGs) and the related rugged
energy landscape with a large diversity of energy barriers makes the
properties of this class of material to be easily tunable Isner2006 ; fan2017
. Nowadays there is a large number of methods available to alter the structure
and corresponding properties of MGs saida2013recovery ; kuchemann2018energy ;
wakeda2015controlled ; pan2018extreme ; park_2008 ; concustell2009structural ;
dmowski2010 ; ding2019ultrafast ; feng2018rejuvenation . Rejuvenation and
relaxation are competing processes in controlling the energy of MGs towards
higher and lower levels resulting, in the end, to a wide range of possible
glassy states concustell2009structural ; mahmoud2021identification .
Rejuvenation of MGs is usually referred to increased disorder resulting in
excess free volume and a raise in energy dmowski2010 ; Liu2010 . Nevertheless,
rejuvenation is not always associated with structural softening and hardness
reduction. Rejuvenation processes done by combining deformation protocols and
thermal treatments revealed a density and hardness increase of some MGs
(ge2017unusual, ; wang2017high, ; miyazaki2016prediction, ). The observed
hardening behavior was associated to confined microplasticity and local
accumulation of irreversible compressive strains Deng2012 ; Zhao2017 .
Furthermore, deformation by cold-rolling and shot-peening results in a
transition from softening to hardening behavior that was associated with the
accumulation of compressive residual stresses overcompensating shear softening
and the corresponding free volume generation Stolpe2014 ; zhang2006 . Highly
rejuvenated MGs obtained under triaxial compression revealed both, strain
hardening and enhanced ductility pan2018extreme ; Pan2020 . The exceptionally
efficient strain-hardening and suppression of shear-banding were associated
with structural relaxation of the extremely rejuvenated MG Pan2020 .

While structural rejuvenation received much attention, there has been less
interest in recent years in exploring stress engineering in MGs. Due to their
very low thermal conductivities and the high cooling rates involved in the
quenching process, MGs undergo differential cooling and large thermal
gradients are generated across the sample leading to residual stress buildup
Ustundag1998 ; Wang2011 . Residual stresses and stress field gradients affect
shear band dynamics and the overall plastic deformation Kosiba2020 ; zhang2006
, hardness Wang2011 ; Haag2014  and even magnetic properties of MGs Wang2012a
. Recently, molecular dynamics (MD) simulations of stress modulated MG
heterostructures revealed that strain hardening together with tensile
ductility can be attained by only modulating the internal residual stress
without altering the structure Yuan2021 .

Here, we focus our attention on clarifying the key factors underlying the
rejuvenation process under triaxial compression and distinguish between
structural and elastic fluctuations. Using a similar setup proposed as by Pan
et al. Pan2020  at the ESRF Extremely Brilliant Source, strain mapping of
extremely deformed MG samples was performed. This novel approach allows to
distinguish hydrostatic from deviatoric strains and to quantify structural and
stress heterogeneities ShakurShahabi2016 ; scudino2019 . Transmission electron
microscopy (TEM) is also applied to highlight intense shear band activity
during the room temperature deformation. Hardness mapping with
microindentation is used to emphasize the competing effects of stress and
structure modulation. Finally, an atomistic understanding of the confined
deformation mechanisms in a notched MG is provided by MD simulations.

Figure 1: Plane strain distribution in the Zr61Ti2Cu25Al12 MG deformed under
triaxial compression. The deviatoric plane strains $\epsilon_{11}$ (a),
$\epsilon_{12}$ (b), and $\epsilon_{22}$ (c) determined from the ellipticity
of the first amorphous diffraction maximum. (d) The volumetric strain $\Delta
q/\Delta q_{0}$ is represented to show the change in free volume. Figure 2:
HAADF-STEM images of two random selected-areas. Two HAADF-STEM images were
taken from lamellae with thickness of $\approx$ 100 nm between the notch area
and along the shear front (black rectangles marked with A and B in the
volumetric strain map). For each cross-section HAADF-STEM images show
alternating contrast changes, which are related to presence of shear bands in
the material. Arrows pointing at the shear bands are a guidance for the eye.
Figure 3: Eigenvalue directions and values determined from the deviatoric
strains. The flow lines for the strain directions $\epsilon_{1}$ and
$\epsilon_{2}$ reveal turbulent behavior around the notches and at along the
shear fronts due to the change from positive to negative strains.

##  II Results

The plane strain distribution can be experimentally determined with high
spatial resolution in MGs binkowski2015 ; ShakurShahabi2016 ; scudino2019 . In
Fig. 1 we show for the first time the antisymmetric nature of the strain
evolution in the MG sample at the notch with unprecedented resolution. As
shown in Figs. 1a and c, the elastic strains parallel and perpendicular to the
loading direction ($\epsilon_{11}$ and $\epsilon_{22}$) display a highly
heterogeneous profile that seems to be stabilized even after releasing the
external load. This behavior indicates a highly inhomogeneous deformation
mechanism, i.e., shear banding ShakurShahabi2016 ; scudino2019 ; Kosiba2020
and, hence, it further confirms that, at low temperature, the shear bands are
the only plastic carriers in MGs. HAADF-STEM images of random selected-areas
along the shear fronts (of highest volumetric strain) and between the notch
area highlight shear band activity (see Fig. 2).

The $\epsilon_{11}$ component shows localization of dilatation at the notches
and along the two shear fronts that extend even away from the notches.
Contraction can be observed between the notch area and ahead of the
intersection point of the two shear fronts. The $\epsilon_{22}$ component
shows exactly an inverse trend to $\epsilon_{11}$. The strong variations in
the shear strain ($\epsilon_{12}$) further highlight the presence of a complex
interpenetrating network of shear bands (Fig. 1b) ShakurShahabi2016 .
Variations in the elastic strain, from compressive to tensile states, are not
only present across the two shear fronts but can be also seen between the
notches. The latter Eigenvector evaluation also allows to determine volumetric
strains (free volume) that are in good agreement with the reciprocal space
dilatation evaluations (Fig. 1d). The highest volumetric strain can be seen at
the intersection point of the two shear fronts, where a high density of shear
bands concentrates, and shows a descending gradient towards the two notches.

The variations in the Eigenvectors of the strain tensor indicates the
direction of the principal strain axes and allows for the localization of
shear bands in MGs scudino2018 . The Eigenvectors exhibit an antisymmetric
sigmoidal profile consisting of large angular variation at the position of
SBs. The evaluation of the Eigenvectors in Fig. 3 highlights a turbulent
strain distribution with a localization of these turbulent areas close to the
notches and to the area of highest dilatation. Such a turbulent strain
distribution correlates to the change in the sign of the strain from
compressive to tensile strain. The sigmoidal profile is stronger at the
position of the two main shear fronts but weaker variations in the sign of the
strain can be also observed between the notches highlighting the presence of
shear activity.

To distinguish between structural and stress fluctuations in the extremely
deformed MG sample, hardness mapping over half of the sample was performed
with microindentation. The hardness map indicates a pronounced local
heterogeneity across the sample (see Fig. 4). As expected, the hardness
decreases monotonically with increasing free volume Yoo2009 . This correlates
to the increase in the volumetric strain, especially at the intersection point
of the two shear fronts of highest dilatation (highlighted in Fig. 1d).
However, below the two shear fronts where the volumetric strain is still
higher, as compared to the undeformed sample, the hardness shows the highest
values. The shear induced softening and the related large volumetric strain do
not decrease the hardness as we would have expected. The only reasonable
explanation for such a questionable hardness mapping could be the presence of
residual stresses and stress field gradients that overcompensate the shear
softening. Such a behavior could already predict that besides the structural
fluctuations (free volume) also the deviatoric strains (stress
heterogeneities) may play a major role in the hardness behavior and in general
in the mechanical properties of MGs zhang2006 .

Figure 4: Hardness mapping with microindentation. The hardness map indicates a
high local heterogeneity over the sample and hints at an unsymmetrical
hardening/softening after compression.

MD simulations were used to provide an atomistic understanding of the complex
deformation mechanisms between the two notches. A simplified rectangular
sample with two notches was deformed. When the yield commences four shear
fronts initiate from the two notches (as highlighted by the atomic strain in
Fig. 5a, upper panel). During loading only two shear bands become dominant and
define the shear direction. Even with continued deformation of the sample,
these two shear bands do not become critical and do not transect the entire
sample (see the lower panel in Fig. 5a). Although the other two shear bands do
not mature, their strain fields perturb the percolation and further
propagation of the dominant shear bands. A better view of the shear band
interaction mechanism is given by the sign of the rotation angle, generally
used to identify vortex-like (rotating) units in glassy structures (Fig. 5b)
Sopu2020 ; Sopu2021 . The sign of rotation defines the shear direction and has
high implication on the dynamics of shear band interaction. When two shear
bands intersect, vortexes of opposite rotation directions come into contact
and the further movement of the shear bands is hindered Sopu2021 . This
process blocks the initial shear front and changes the shear band dynamics and
morphology, ultimately leading to shear band blocking, branching and
multiplication Sopu2020 . Shear band multiplication and blocking processes are
particularly important as they play a decisive role in enhancing the
plasticity of MGs, what could explain the extended compressive plastic flow in
the central notched region pan2018extreme .

However, structural variation is not only confined to the narrow shear bands.
When analyzing the volumetric strain one can see that the local density does
not increase only at the position of the shear bands (Fig. 5c). The entire
area at the notch, between the four shear fronts, shows larger variations with
respect to the structure away from the notches. The complex shear band network
confines the material and changes the local strain/stress state. In fact, this
could be responsible for the observed homogeneous rejuvenation of MG when
deformed under triaxial compression Pan2020 . Indeed, when analyzing the
strain tensor ($\epsilon_{12}$) distribution obtained from experiments (Fig.
1b) and MD simulations (Fig. 5d) a striking resemblance can be observed. The
variation in the elastic strain from compressive to tensile states closely
follows the shear band front. Moreover, at the position where the shear band
branches, the sign of the strain shows a sinusoidal strain profile (tensile
vs. compressive), highlighting the formation of multiple shear bands.

##  III Discussion

From the above results it can be stated that triaxial deformation of MGs
induces both plane (directional) and deviatoric (volumetric) strains. Hence,
the reason for the observed extreme rejuvenation corresponds to both,
structural softening and residual stresses imparted during the deformation
process. Additionally, the stress heterogeneity is strongly modulated and
correlated with the intense shear band activity, as indicated by the shear
fluctuations (change in the sign from positive to negative) in the area
between the two notches (Fig. 1). It is well accepted that shear banding is
not only confined to a narrow thin layer of about 20 nm but extends over a
longer range (tens of micrometers) across the shear band Maass2011 ; liu2018 ;
scudino2018 . Experimental observations of high-energy XRD across an
individual shear band revealed that the strain at opposite sides of a shear
band changes sign from compressive to tensile scudino2019 ; scudino2018 ;
Kosiba2020 . Likewise, the strains/stresses are created in an inhomogeneous
manner during triaxial compression (see Fig. 1). Moreover, these fluctuations
are even there after unloading, cutting and polishing the highly deformed
sample. Previously it has been stated that residual stresses are released
after the small sample is machined within the notches of a pre-deformed rod
Pan2020 . There, the degree of rejuvenation was quantified using the
exothermic heat of relaxation ($\Delta H_{rel}$). The measured $\Delta
H_{rel}$ was mainly related to structural changes, and the volume
accumulation, and a only small fraction of this energy corresponds to residual
stresses Stolpe2014 . Consequently, the contributions of stress
heterogeneity/fluctuations and residual stresses to rejuvenation and strain-
hardening were underestimated Pan2020 . Nevertheless, the presence of residual
stresses is reflected in the observed large change in density that exceeds
0.2$\%$ (Fig. 1) johnson2005 ; wang_2006 . Indeed, the relative change in
density as a function of cooling rate is usually lower than 0.07 $\%$ Hu2001 ;
Fu2013 . Similar variations in density were found in cold-rolled metallic
glasses Haruyama2013 . Large variations in density up to 0.16 $\%$ have been
observed in heavily deformed MGs Vempati2012 . In this case, the density
fluctuations correspond to both, enhanced disorder (structural rejuvenation)
and the stress constraint introduced by the residual stresses, respectively.

Figure 5: MD simulations of the deformation mechanisms in the notched MG
alloy. Structural variations in the simulated Cu64Zr36 MG at two loading
stages, at an early stage 1 of deformation of shear band formation (upper
panels) and a later stage 2 when the shear bands interact (lower panels).
Evolution of the atomic von Mises strain (a), the rotation angle (b), the
volumetric strain (c) and the strain tensor ($\epsilon_{12}$) (d).

The effects of compressive and tensile residual stresses after mechanical
loading are especially important when considering that they can overcompensate
or undercompensate the softening caused by shear banding and free volume
generation Launey2008 ; Yuan2021 . These competing effects may help to explain
why the hardness map does not perfectly correlate to the volumetric strain
(free volume), especially at the position between the two notches (see Fig.
1c). Here, the volumetric strain shows higher values (Fig. 1c) and hence,
lower hardness values are expected. Contrarily, the microhardnesses map
indicates hardening in the region between the notches (Fig. 4). In this area
the elastic strains $\epsilon_{11}$ and $\epsilon_{22}$ display large
fluctuations (Figs. 1a, b and d). Since the sample is symmetric perpendicular
to the loading direction also $\epsilon_{33}$ is tensile. The residual tensile
stress can be seen as a back stress Zhu2019 ; Yuan2021  that must be overcome
by the indent to penetrate the system and to switch locally to a compressive
stress state. Namely, the contact pressure underneath the indent that is
required to initiate the yield significantly increases in these regions
Wang2011 ; Chen2008 .

##  IV Conclusion

In conclusion, this research has significantly advanced our understanding of
MG rejuvenation under constrained deformation protocols. Based on high energy
scanning X-ray diffraction strain mapping, STEM and microindentation one can
distinguish between structural and elastic fluctuations, the two key factors
for the observed extreme rejuvenation in triaxial compression. STEM
characterization shows that structural rejuvenation under room temperature
deformation relates to shear-induced softening and dilatation (large
volumetric strain). An interpenetrating network of multiple shear bands not
only increases the free volume but also results in large variations in the
elastic strain. High energy scanning X-ray diffraction strain mapping reveals
large elastic fluctuations corresponding to the intense shear band activity.
The complex superimposed spatial stress and structure heterogeneities reflect
in the weak correlation between hardness and structural softening and support
the competing hardening-softening mechanism. Hence, the classical description
of rejuvenation should be reconsidered and future focus should be also given
to stress-driven modulation. Stress and structural engineering can be used
simultaneously to enhance rejuvenation beyond the limits known so far and,
consequently, can help to design MGs with enhanced ductility and strain
hardening ability.

##  V Acknowledgments

The authors acknowledge financial support by the Deutsche
Forschungsgemeinschaft (DFG) through Grant No. SO 1518/1-1, the European
Research Council under the ERC Advanced Grant INTELHYB (grant
ERC-2013-ADG-340025), Austrian Science Fund (FWF): Y1236-N37 and the China
Scholarship Council (CSC, 201806220096). We acknowledge the European
Synchrotron Radiation Facility (ESRF) for provision of synchrotron radiation
facilities and we would like to thank Dr. Carlotta Giacobbe for assistance and
support in using beamline ID 11. Calculations for this research were conducted
on the Lichtenberg high performance computer of the Technische Universität
Darmstadt.

##  VI Methods

###  VI.1 Preparation of metallic glass samples

Metallic glass with a nominal composition of Zr61Ti2Cu25Al12 (at.$\%$) was
prepared by arc-melting a mixture of high-purity ($>$ 99.9$\%$) metals at
least four times in a titanium-gettered high-purity argon atmosphere.
Cylindrical samples with a diameter of 4.5 mm and a length of 70 mm were
fabricated by suction casting into a copper mold. The structure of the glassy
phase of the as-cast MG samples was ascertained by X-ray diffraction (XRD)
using a Rigaku Dmax-2550 diffractometer with Cu-K$\alpha$ radiation generated
at 40 kV. Specimens with a length/diameter ratio of 2 were cut from the rods
by a diamond saw, and the surface roughness of the two parallel ends were
carefully polished. Notched specimens with a circumferential notch dimension
(notch diameter, d $\times$ notch height, h) of 2.5 $\times$ 0.5 mm2 were
produced using a 0.5 mm diameter round-head grooving knife at a feed rate of
0.06 mm/r in a lathe machine with a running speed of 500 r/min followed by
fine polishing and final cleaning in an ultrasonic bath. Then, the notched
specimens were compressed along the cylindrical axis at room temperature using
MTS CMT5205 device at a strain rate of 10-4s-1 until the width of the notch
was reduced by $\approx$ 40$\%$ plastic flow. Due to the notch constrain the
central zone is deformed under triaxial compression, suppressing relaxation
Pan2020 . The deformed specimens were then machined along the central lateral
and longitudinal planes. The longitudinal cross-sections were mechanically
polished to mirror finish for the following high energy XRD and
nanoindentation experiments.

###  VI.2 Structural characterization

High energy scanning XRD was used to perform strain mapping on the cross-
section of the deformed specimen, which was carried out at the nanostation of
the beamline ID11 at ESRF. A wavelength of 0.01919 $nm$ was used and the beam-
size was focused to 500 $\times$ 500 $nm$. The region of interest was scanned
with 3500 points with a spacing of 1 $\mu m$ each in the fast scanning
direction and 65 with a spacing of 50 $\mu m$ in the slow scanning direction
leading to a total of 227500 two-dimensional diffraction patterns recorded
with an exposure time of 10 $ms$ with an Eiger CdTe 4M fast pixel array
detector. In the fast direction 50 frames were summed to increase the signal-
to-noise ratio. The newly upgraded Extremely Brilliant Source (EBS) of ESRF
allowed measurement of half of the sample in about 2 hours.

The data was calibrated using a CeO2 reference and integrated using the pyFAI
software. For the evaluation of plane strains, the integration was performed
in 72 azimuthal sections. The volumetric strains were evaluated using the
average center of mass of the first strong diffraction maximum $q_{COM}/q_{0}$
of the glassy material. The fitting procedure is reproduced for a
representative diffraction pattern in Fig. S1.

For the plane strains $\epsilon_{11}$, $\epsilon_{12}$ and $\epsilon_{12}$,
the normalized center of mass of the first diffraction maximum
$q(\phi)/q_{COM}$ was fitted using the following quation for the azimuthal
angle $\phi$ similar to the established procedure ShakurShahabi2016

$\phi=\epsilon_{11}cos(\phi)^{2}+\epsilon_{12}\mathrm{cos}(\phi)sin(\phi)+\epsilon_{22}sin(\phi)^{2}.$
(1)

The high flux of EBS and new fast pixel array detectors with single photon
counting ability allow correlation of the mesoscale information occurring in
the plastic deformation of the notched MG sample with the formation of
individual shear bands and their interactions. Strain mapping allows
identification and visualization of shear bands in the bulk material in 2D.

Thin lamella for scanning transmission electron microscopy (STEM) were
prepared with a final thickness of 100 $nm$ by focused ion-beam milling using
a Zeiss Auriga dual-beam workstation. One lamella was taken from the middle
part between the two notches of the specimen. A second lamella was taken from
a position inclined by 45∘ with respect to the notch direction. High Angle
Annular Dark Field (HAADF) STEM images were recorded for each lamella with a
spot size of 1 $nm$.

Macroscopic Vickers hardness maps were recorded using a universal fully-
automatic DuraScan70G5 laboratory ZwickRoell hardness tester. The areal
distribution of indentation points was recorded in a grid of 200 $\mu m$ x 200
$\mu m$ with a total number of 522 indents in a matrix of 26x21. To avoid the
influence of the elastic-plastic stress field, the indents were spaced three
diagonals apart. The indentation depth of each indent was approximately 15
$\mu m$, calculated using the 136∘ opening angle of the pyramid-shaped diamond
indenter. A load of 9,807 $N$ was applied for each impression, which
corresponds to the small load range of HV1.

###  VI.3 Molecular dynamics simulations

A qualitative and elemental atomistic picture of the confined deformation in a
notched MG was provided by means of classical molecular dynamics simulations
carried out with the LAMMPS software plimpton1995 . Because no appropriate
interatomic potentials exist for Zr61Ti2Cu25Al12, the Cu64Zr36 MG was used as
a prototype material since a reliable EAM potential was made available by
Mendelev et al. Mendelev2007 . The amorphous Cu64Zr36 alloy structures
predicted by the potential are in good agreement with X-ray diffraction data
and the potential can be reliably used to simulate the structure and
properties of amorphous Cu-Zr alloys Sopu2020 ; Sopu2021 .

The starting liquid structure was created by randomly distributing 968000
atoms in a box of $111.0\times 55.5\times 2.5\;nm$ with periodic boundary
conditions. The amorphous structure was obtained by quenching from the melt
starting from 2000 K. More details about quenching procedure can be found in
ref. Moitzi2020 . The cooled sample was replicated twice along the
$z$-direction and two notches were created along the $y$-direction. Open
boundaries in $y$-direction and periodic boundaries in the $x$\- and
$z$-directions were applied. The time step was set to 2 fs.

Uniaxial tensile simulations along the $z$-direction were performed on the
notched MG sample under a constant engineering strain rate of
4$\times$10${}^{7}~{}s^{-1}$. The stress along $x$-direction was relaxed to
allow for lateral contraction using Nose/Hoover barostat plimpton1995 . A low
temperature of 50 K was chosen to compensate for the very high strain rate and
guarantee a regime of localized deformation and shear band formation. The
deformation mechanisms and the related structural fluctuations at two
different strain levels were analyzed using the OVITO software Stukowski2010
by calculating and visualizing the von Mises strains, simultaneously occurring
directional rotation fields, volumetric strains and the strain tensors. All
these micromechanical quantities were evaluated for each atom in the system
from the relative motion of the neighboring atoms within a cutoff range of 1
$nm$.

## References

  * (1) Isner, B. A. & Lacks, D. J.  Generic rugged landscapes under strain and the possibility of rejuvenation in glasses.  _Phys. Rev. Lett._ 96, 025506 (2006).  URL https://link.aps.org/doi/10.1103/PhysRevLett.96.025506.
  * (2) Fan, Y., Iwashita, T. & Egami, T.  Energy landscape-driven non-equilibrium evolution of inherent structure in disordered material.  _Nature Communications_ 8, 15417 (2017).  URL https://doi.org/10.1038/ncomms15417.
  * (3) Saida, J., Yamada, R. & Wakeda, M.  Recovery of less relaxed state in zr-al-ni-cu bulk metallic glass annealed above glass transition temperature.  _Applied Physics Letters_ 103, 221910 (2013).
  * (4) Küchemann, S. _et al._ Energy storage in metallic glasses via flash annealing.  _Advanced Functional Materials_ 28, 1805385 (2018).
  * (5) Wakeda, M., Saida, J., Li, J. & Ogata, S.  Controlled rejuvenation of amorphous metals with thermal processing.  _Scientific reports_ 5, 1–8 (2015).
  * (6) Pan, J. _et al._ Extreme rejuvenation and softening in a bulk metallic glass.  _Nature communications_ 9, 1–9 (2018).
  * (7) Park, K. W.  Effect of the atomic packing density on the structural change rate of amorphous alloys under elastostatic stress.  _Metals and Materials International_ 14, 159–163 (2008).
  * (8) Concustell, A., Méar, F., Surinach, S., Baró, M. & Greer, A.  Structural relaxation and rejuvenation in a metallic glass induced by shot-peening.  _Philosophical Magazine Letters_ 89, 831–840 (2009).
  * (9) Dmowski, W. _et al._ Structural rejuvenation in a bulk metallic glass induced by severe plastic deformation.  _Acta Materialia_ 58, 429 – 438 (2010).  URL http://www.sciencedirect.com/science/article/pii/S135964540900617X.
  * (10) Ding, G. _et al._ Ultrafast extreme rejuvenation of metallic glasses by shock compression.  _Science advances_ 5, 6249 (2019).
  * (11) Feng, S. _et al._ Rejuvenation by weakening the medium range order in zr46cu46al8 metallic glass with pressure preloading: A molecular dynamics simulation study.  _Materials & Design_ 158, 248–255 (2018).
  * (12) Mahmoud, Z. H. _et al._ identification of rejuvenation and relaxation regions in a zr-based metallic glass induced by laser shock peening.  _Journal of Materials Research and Technology_ 11, 2015–2020 (2021).
  * (13) Liu, J. W., Cao, Q. P., Chen, L. Y., Wang, X. D. & Jiang, J. Z.  Shear band evolution and hardness change in cold-rolled bulk metallic glasses.  _Acta Materialia_ 58, 4827–4840 (2010).  URL https://www.sciencedirect.com/science/article/pii/S1359645410002934.
  * (14) Ge, T. _et al._ Unusual energy state evolution in ce-based metallic glass under high pressure.  _J. Appl. Phys._ 121, 205109 (2017).
  * (15) Wang, C. _et al._ High stored energy of metallic glasses induced by high pressure.  _Appl. Phys. Lett._ 110, 111901 (2017).
  * (16) Miyazaki, N., Wakeda, M., Wang, Y. J. & Ogata, S.  Prediction of pressure-promoted thermal rejuvenation in metallic glasses.  _Npj Comput. Mater._ 2, 1–9 (2016).
  * (17) Deng, C. & Schuh, C. A.  Atomistic mechanisms of cyclic hardening in metallic glass.  _Appl. Phys. Lett._ 100, 251909 (2012).  URL https://doi.org/10.1063/1.4729941.
  * (18) Zhao, D., Zhao, H., Zhu, B. & Wang, S.  Investigation on hardening behavior of metallic glass under cyclic indentation loading via molecular dynamics simulation.  _Applied Surface Science_ 416, 14–23 (2017).  URL https://www.sciencedirect.com/science/article/pii/S0169433217311522.
  * (19) Stolpe, M., Kruzic, J. J. & Busch, R.  Evolution of shear bands, free volume and hardness during cold rolling of a zr-based bulk metallic glass.  _Acta Materialia_ 64, 231–240 (2014).  URL https://www.sciencedirect.com/science/article/pii/S1359645413007982.
  * (20) Zhang, Y., Wang, W. H. & Greer, A. L.  Making metallic glasses plastic by control of residual stress.  _Nat Mater_ 5, 857–860 (2006).  URL http://dx.doi.org/10.1038/nmat1758.
  * (21) Pan, J., Ivanov, Y. P., Zhou, W. H., Li, Y. & Greer, A. L.  Strain-hardening and suppression of shear-banding in rejuvenated bulk metallic glass.  _Nature_ 578, 559–562 (2020).  URL https://doi.org/10.1038/s41586-020-2016-3.
  * (22) Ustundag, E. _et al._ Residual stresses in bulk metallic glasses due to differential cooling or thermal tempering.  _MRS Online Proceedings Library_ 554, 431–436 (1998).  URL https://doi.org/10.1557/PROC-554-431.
  * (23) Wang, L., Bei, H., Gao, Y. F., Lu, Z. P. & Nieh, T. G.  Effect of residual stresses on the hardness of bulk metallic glasses.  _Acta Materialia_ 59, 2858–2864 (2011).  URL https://www.sciencedirect.com/science/article/pii/S1359645411000310.
  * (24) Kosiba, K. _et al._ Guiding shear bands in bulk metallic glasses using stress fields: A perspective from the activation of flow units.  _Phys. Rev. B_ 102, 134113 (2020).  URL https://link.aps.org/doi/10.1103/PhysRevB.102.134113.
  * (25) Haag, F., Beitelschmidt, D., Eckert, J. & Durst, K.  Influences of residual stresses on the serrated flow in bulk metallic glass under elastostatic four-point bending - a nanoindentation and atomic force microscopy study.  _Acta Materialia_ 70, 188–197 (2014).  URL https://www.sciencedirect.com/science/article/pii/S1359645414000743.
  * (26) Wang, J.-Q. _et al._ Rapid degradation of azo dye by fe-based metallic glass powder.  _Adv. Funct. Mater._ 22, 2567–2570 (2012).  URL https://doi.org/10.1002/adfm.201103015.
  * (27) Yuan, X., Şopu, D. & Eckert, J.  Origin of strain hardening in monolithic metallic glasses.  _Phys. Rev. B_ 103, L140107 (2021).  URL https://link.aps.org/doi/10.1103/PhysRevB.103.L140107.
  * (28) Shakur Shahabi, H. _et al._ Mapping of residual strains around a shear band in bulk metallic glass by nanobeam x-ray diffraction.  _Acta Materialia_ 111, 187–193 (2016).  URL https://www.sciencedirect.com/science/article/pii/S1359645416301884.
  * (29) Scudino, S.  Mechanism of shear banding during cold rolling of a bulk metallic glass.  _Journal of Alloys and Compounds_ 773, 883–889 (2019).  URL https://www.sciencedirect.com/science/article/pii/S0925838818335606.
  * (30) Binkowski, I. _et al._ Sub-micron strain analysis of local stick-slip motion of individual shear bands in a bulk metallic glass.  _Applied Physics Letters_ 107, 221902 (2015).  URL http://dx.doi.org/10.1063/1.4936388.  eprint http://dx.doi.org/10.1063/1.4936388.
  * (31) Scudino, S. & Şopu, D.  Strain distribution across an individual shear band in real and simulated metallic glasses.  _Nano Lett._ 18, 1221–1227 (2018).  URL https://doi.org/10.1021/acs.nanolett.7b04816.
  * (32) Yoo, B.-G., Kim, Y.-J., Oh, J.-H., Ramamurty, U. & Jang, J.-i.  On the hardness of shear bands in amorphous alloys.  _Scripta Materialia_ 61, 951–954 (2009).  URL https://www.sciencedirect.com/science/article/pii/S1359646209004953.
  * (33) Şopu, D., Scudino, S., Bian, X. L., Gammer, C. & Eckert, J.  Atomic-scale origin of shear band multiplication in heterogeneous metallic glasses.  _Scripta Materialia_ 178, 57–61 (2020).  URL http://www.sciencedirect.com/science/article/pii/S1359646219306499.
  * (34) Şopu, D. _et al._ From elastic excitations to macroscopic plasticity in metallic glasses.  _Applied Materials Today_ 22, 100958 (2021).  URL https://www.sciencedirect.com/science/article/pii/S2352940721000238.
  * (35) Maaß, R., Klaumünzer, D. & Löffler, J. F.  Propagation dynamics of individual shear bands during inhomogeneous flow in a zr-based bulk metallic glass.  _Acta Materialia_ 59, 3205–3213 (2011).  URL http://www.sciencedirect.com/science/article/pii/S1359645411000760.
  * (36) Liu, C. & Maaß, R.  Elastic fluctuations and structural heterogeneities in metallic glasses.  _Adv. Funct. Mater._ 28, 1800388 (2018).  URL https://doi.org/10.1002/adfm.201800388.
  * (37) Johnson, W. L. & Samwer, K.  A universal criterion for plastic yielding of metallic glasses with a $(t/{T}_{g}{)}^{2/3}$ temperature dependence.  _Phys. Rev. Lett._ 95, 195501 (2005).  URL http://link.aps.org/doi/10.1103/PhysRevLett.95.195501.
  * (38) Wang, W. H.  Correlations between elastic moduli and properties in bulk metallic glasses.  _Journal of Applied Physics_ 99, 93506 (2006).
  * (39) Hu, X., Ng, S. C., Feng, Y. P. & Li, Y.  Cooling-rate dependence of the density of ${\mathrm{pd}}_{40}{\mathrm{ni}}_{10}{\mathrm{cu}}_{30}{\mathrm{p}}_{20}$ bulk metallic glass.  _Phys. Rev. B_ 64, 172201 (2001).  URL https://link.aps.org/doi/10.1103/PhysRevB.64.172201.
  * (40) Fu, W., Sun, Y. & Zhang, W.  The effect of cooling rate on microstructure and mechanical properties of zr-based bulk metallic glasses.  _Advances in Materials Science and Engineering_ 2013, 826758 (2013).  URL https://doi.org/10.1155/2013/826758.
  * (41) Haruyama, O. _et al._ Characterization of free volume in cold-rolled zr55cu30ni5al10 bulk metallic glasses.  _Acta Materialia_ 61, 3224–3232 (2013).  URL https://www.sciencedirect.com/science/article/pii/S1359645413001286.
  * (42) Vempati, U. K., Valavala, P. K., Falk, M. L., Almer, J. & Hufnagel, T. C.  Length-scale dependence of elastic strain from scattering measurements in metallic glasses.  _Phys. Rev. B_ 85, 214201 (2012).  URL https://link.aps.org/doi/10.1103/PhysRevB.85.214201.
  * (43) Launey, M. E., Busch, R. & Kruzic, J. J.  Effects of free volume changes and residual stresses on the fatigue and fracture behavior of a zr-ti-ni-cu-be bulk metallic glass.  _Acta Materialia_ 56, 500–510 (2008).  URL https://www.sciencedirect.com/science/article/pii/S135964540700688X.
  * (44) Zhu, Y. & Wu, X.  Perspective on hetero-deformation induced (hdi) hardening and back stress.  _Materials Research Letters_ 7, 393–398 (2019).  URL https://doi.org/10.1080/21663831.2019.1616331.
  * (45) Chen, M. W.  Mechanical behavior of metallic glasses: Microscopic understanding of strength and ductility.  _Annual Review of Materials Research_ 38, 445–469 (2008).
  * (46) Plimpton, S.  Fast parallel algorithms for short-range molecular-dynamics.  _J. Comput. Phys._ 117, 1–19 (1995).
  * (47) Mendelev, M. I., Sordelet, D. J. & Kramer, M. J.  Using atomistic computer simulations to analyze x-ray diffraction data from metallic glasses.  _J. Appl. Phys._ 102, 043501 (2007).  URL http://link.aip.org/link/?JAP/102/043501/1.
  * (48) Moitzi, F. _et al._ Chemical bonding effects on the brittle-to-ductile transition in metallic glasses.  _Acta Materialia_ 188, 273–281 (2020).  URL http://www.sciencedirect.com/science/article/pii/S135964542030094X.
  * (49) Stukowski, A.  Visualization and analysis of atomistic simulation data with ovito-the open visualization tool.  _Modell. Simul. Mater. Sci. Eng._ 18, 015012 (2010).