# Structural disorder by octahedral tilting in inorganic halide perovskites:
New insight with Bayesian optimization

Jingrui Li<EMAIL_ADDRESS>Electronic Materials Research Laboratory,
Key Laboratory of the Ministry of Education and International Center for
Dielectric Research, School of Electronic Science and Engineering &
International Joint Laboratory for Micro/Nano Manufacturing and Measurement
Technology, Xi’an Jiaotong University, Xi’an 710049, China  Fang Pan
Electronic Materials Research Laboratory, Key Laboratory of the Ministry of
Education and International Center for Dielectric Research, School of
Electronic Science and Engineering & International Joint Laboratory for
Micro/Nano Manufacturing and Measurement Technology, Xi’an Jiaotong
University, Xi’an 710049, China  Guo-Xu Zhang  MIIT Key Laboratory of Critical
Materials Technology for New Energy Conversion and Storage, School of
Chemistry and Chemical Engineering, Harbin Institute of Technology, Harbin
150001, China  Zenghui Liu  Electronic Materials Research Laboratory, Key
Laboratory of the Ministry of Education and International Center for
Dielectric Research, School of Electronic Science and Engineering &
International Joint Laboratory for Micro/Nano Manufacturing and Measurement
Technology, Xi’an Jiaotong University, Xi’an 710049, China  Hua Dong  Key
Laboratory for Physical Electronics and Devices of the Ministry of Education
and Shaanxi Key Lab of Information Photonic Technique, School of Electronic
Science and Engineering, Xi’an Jiaotong University, Xi’an 710049, China  Dawei
Wang  School of Microelectronics and Key Lab of Micro-Nano Electronics and
System Integration of Xi’an City, Xi’an Jiaotong University, Xi’an 710049,
China  Zhuangde Jiang  State Key Laboratory for Manufacturing Systems
Engineering & International Joint Laboratory for Micro/Nano Manufacturing and
Measurement Technology, Xi’an Jiaotong University, Xi’an 710049, China  Wei
Ren  Electronic Materials Research Laboratory, Key Laboratory of the Ministry
of Education and International Center for Dielectric Research, School of
Electronic Science and Engineering & International Joint Laboratory for
Micro/Nano Manufacturing and Measurement Technology, Xi’an Jiaotong
University, Xi’an 710049, China  Zuo-Guang Ye  Department of Chemistry and 4D
LABS, Simon Fraser University, Burnaby, British Columbia V5A 1S6, Canada
Milica Todorović<EMAIL_ADDRESS>Department of Mechanical and
Materials Engineering, University of Turku, FI-20014 Turku, Finland  Patrick
Rinke  Department of Applied Physics, Aalto University, P.O.Box 11100,
FI-00076 AALTO, Finland

###### Abstract

Structural disorder is common in metal-halide perovskites and important for
understanding the functional properties of these materials. First-principles
methods can address structure variation on the atomistic scale, but they are
often limited by the lack of structure-sampling schemes required to
characterize the disorder. In this work, structural disorder in the benchmark
inorganic halide perovskites $\text{CsPbI}_{3}$ and $\text{CsPbBr}_{3}$ is
computationally studied in terms of the three octahedral-tilting angles. The
consequent variation in energetics and properties are described by three-
dimensional potential-energy surfaces (PESs) and property landscapes,
delivered by Bayesian Optimization Structure Search method with integrated
density-functional-theory (DFT) calculations. The rapid convergence of the PES
with about 200 DFT data points in three-dimensional searches demonstrates the
power of active learning and strategic sampling with Bayesian optimization.
Further analysis indicates that disorder grows with increasing temperature,
and reveals that the materials band gap at finite temperatures is a
statistical mean over disordered structures.

##  I Introduction

Structural disorder in materials has become an important topic in both
experimental and computational materials science [1, 2, 3, 4, 5, 6, 7, 8].
Disorder phenomena are very common in emergent functional materials, as ion-
mixing or doping strategies are widely applied to obtain high performance and
stability. Notable examples include ferroelectric oxide perovskite solid
solutions [9, 10, 11], multiple-cation metal halide perovskites for solar
cells and light-emitting diodes [12, 13], and kesterite photovoltaic materials
[14, 15]. Structural disorder can also occur in pure materials with perfect
stoichiometry, if the material exhibits a series of stable structures with
similar thermodynamic free energies, and the energy barriers separating free-
energy minima are easy to overcome. Consequently, properties of complex
functional materials often arise from the thermal population of a number of
low-energy structures.

Cesium lead iodide ($\text{CsPbI}_{3}$) is an example of the latter category
that has received increasing attention in recent years. It is a promising
photovoltaic material as $\text{CsPbI}_{3}\,$-based perovskite solar cells
have a reported power conversion efficiency above $21\%$ and good stability
[16]. It is generally believed that $\text{CsPbI}_{3}$ adopts the cubic
$\upalpha$ phase (space group $Pm\bar{3}m$) at high temperatures, which
successively converts into the tetragonal $\upbeta$ phase ($P4/mbm$) at
$533\leavevmode\nobreak\ \text{K}$ then the orthorhombic $\upgamma$ phase
($Pnma$) at $448\leavevmode\nobreak\ \text{K}$ upon cooling [17]. The
structural details of the $\upbeta$ phase are not well established. In a joint
experimental-theoretical study, Marronnier et al. analyzed the phonon
instabilities of this phase [18]. Jinnouchi et al. related the phase
transition to the change in effective radius of $\text{Cs}^{+}$ by thermal
fluctuations [19]. Klarbring suggested that the macroscopic tetragonal phase
consists of dynamically fluctuating orthorhombic structures [20]. Yang et al.
further showed that some other low-symmetry structures are also involved [21].
Dynamical disorder has previously been found in other isostructural inorganic
halide perovskites [22, 20, 23, 24]. For example, a recent density-functional-
theory (DFT) and molecular dynamics study indicated that the disorder in cubic
$\text{CsPbBr}_{3}$ is closely related to the octahedral-tilting dynamics
[24]. The nature of disorder in the high-temperature phases of inorganic
halide perovskites remains unclear to-date, thus it is important to tackle
this issue given the link between the atomic structure and functional
properties of materials. Better knowledge about the involved structures,
energetics, and mechanisms of structural fluctuations is crucial to understand
the atomistic origin of functional properties and to design the next-
generation high-performance perovskites.

To gain theoretical insight into structural disorder and its energetics, it is
useful to study the system’s multidimensional potential energy surface (PES)
[25, 26], i.e., the total energy as a function of the relevant degrees of
freedom (DOFs) for disorder. Energy differences governing structural disorder
may be very small, so we must employ accurate DFT methods which have
celebrated great successes in modeling molecules and condensed matter systems
[27]. DFT can supply one- and two-dimensional PESs with an equispaced-grid
approach (such as in Ref. [20]). Yet PES computations quickly become
intractable as the number of DOFs and thus the number of grid points grows.
Attempts to construct accurate high-dimensional PESs based on first-principles
calculations can be traced back several decades, when a series of
interpolation schemes were developed to study polyatomic chemical reactions
[28, 29]. However, it is not easy to directly apply this approach to complex-
materials problems such as solid states and interfaces. Among the many atomic
DOFs of a complex material, there are usually several “principal” DOFs that
play a decisive role in structural energetics, while other DOFs are dependent
or of secondary importance. As a result, such a PES could be much more
complicated than that of a polyatomic molecular system. It might have multiple
maxima and minima, making simple interpolation approaches difficult.

In the past decade, machine learning (ML) has made an impact in the field of
computational materials science [30, 31, 32, 33, 34]. ML and DFT can be
combined to approximate PESs [35, 36] and to accelerate structure search with
improved sampling schemes [37, 38, 39, 40, 41] and accelerated force
evaluation [42, 43, 44, 45]. Force fields can of course provide the PESs in
terms of all atomic DOFs yet typically lack the accuracy. Modern machine-
learned force fields or interatomic potentials [46, 47, 48, 49, 50, 51, 52,
53] also provide DFT accuracy, but they take a long time to train in materials
with many chemical species. We are here, however, interested in PESs that we
can visualize, which is not possible for high-dimensional PESs of force
fields.

Recently, some of us have developed the Bayesian Optimization Structure Search
(BOSS) approach, an ML-based structure search scheme for accelerated and
unbiased PES computation [54]. BOSS couples state-of-the-art DFT calculations
with the active learning Bayesian optimization (BO) technique. It employs
Gaussian-process to fit a PES surrogate model to DFT data points, and refines
this model until convergence is achieved by acquiring further data with a
smart sampling strategy. In such a way, BOSS can construct a complete high-
dimensional PES using a relatively modest number of energy data points. BOSS
has already been applied to solve problems such as conformer search for
organic molecules [54, 55] and adsorption of organic molecules at
semiconductor surfaces [54], to resolve different organic adsorbate types and
film growth at metallic surfaces [56, 57, 58], and to identify the interface
geometries of inorganic or organic materials at perovskite surfaces [59, 60,
61]. The application of BO in materials science is not limited to structure
search, e.g., it was used to train the force fields to study phase transition
in hybrid perovskites [19], to select optimal ML hyperparameters [62, 63], or
to efficiently solicit experimental data [64, 65].

In this work, we apply BOSS to study disorder in the benchmark systems
$\text{CsPbI}_{3}$ and the isostructural bromide perovskite
($\text{CsPbBr}_{3}$). The latter is a standard green-emitting material for
perovskite light-emitting diodes [66], another of today’s key materials in
perovskite optoelectronics. We select the three octahedral-tilting angles [67]
as the principal DOFs. Previous experimental and theoretical studies [68, 69,
22, 20, 70, 24] have established that they are key to determining the
materials structure. One of the main objectives of this work is to compute a
realistic and accurate three-dimensional (3D) PES that best reflects the
structures inside the perovskite materials. The smart data sampling by BOSS
enables us to include all three octahedral-tilting angles as variables for the
3D PES. We thus go beyond the previously established two-dimensional (2D) view
of Klarbring [20], which nevertheless is considered as reference for our
study. In addition we use BOSS to generate a 3D profile of the important band
gap property. This enables us to gain insight into the effects of structural
disorder and temperature on functional properties. Important information about
the structural disorder in both $\text{CsPbI}_{3}$ and $\text{CsPbBr}_{3}$ is
obtained by gathering all these findings.

##  II Bayesian optimization structure search scheme

###  II.1 Octahedral tilting and the disorder search spaces

The disorder search spaces for the benchmark systems $\text{CsPbI}_{3}$ and
$\text{CsPbBr}_{3}$ are spanned by the three octahedral-tilting angles. We
adopted the Glazer notation $a^{r_{a}}b^{r_{b}}c^{r_{c}}$ of octahedral
tilting in perovskites [67] throughout this manuscript. Here, $x^{r_{x}}$
denotes the tilt around the lattice vector $x$ ($=a,b,c$) of the quasi-cubic
cell, and $r_{x}=+,-,0$ indicate the in-phase, out-of-phase, and zero tilting
modes (sometimes shortened as “tilts”), respectively (Figure 1a). The
corresponding tilting angle is denoted by $\theta_{x}$ hereafter. Different
combinations of tilting modes and angles around the three lattice vectors
result in different tilting patterns. The perovskite lattice accordingly
exhibits different types of distortions. Originally, Glazer had derived 23
tilting patterns [67, 71] and later Woodward revised this framework into 15
based on experimental results and space-group analysis [72, 73, 74]. In this
work we considered all possible tilting patterns (listed in Table 1) that are
compatible with both Glazer’s and Woodward’s conventions. A recent study of
$\text{CsSnI}_{3}\,$, which is isostructural and closely related to our
systems, indicates that this is necessary [75]. Results of DFT full relaxation
for both systems with different tilting patterns are provided in Section S1 in
Supporting Information (SI).

FIG 1: Definition of octahedral-tilting configuration spaces for $\text{CsPbX}_{3}$ ($\text{X}=\text{I},\text{Br}$) used in this work. (a) Tilting modes along a lattice vector: in-phase tilt (left, Glazer notation “$+$”) and out-of-phase tilt (right, “$-$”), with $\text{Pb}^{2+}$ cation colored in blue and $\text{X}^{-}$ anions in red, octahedral rotation directions indicated by yellow arrows; (b) the four configuration spaces representing different categories of octahedral-tilting patterns; and (c) the total configuration space. Spaces in (b) and octants in (c) are colored to indicate the correspondence between them, e.g., the “mmm” space in (b) and the “mmm” octant in (c) are both colored in green. TABLE 1: Attribution of Glazer tilting patterns to different categories (disorder search spaces). Tilting | ppp | ppm | pmm | mmm
---|---|---|---|---
$a^{0}a^{0}a^{0}$ | ✓ | ✓ | ✓ | ✓
$a^{+}b^{0}b^{0}$ | ✓ | ✓ | ✓ |
$a^{-}b^{0}b^{0}$ |  | ✓ | ✓ | ✓
$a^{+}a^{+}c^{0}$, $a^{+}b^{+}c^{0}$ | ✓ | ✓ |  |
$a^{+}b^{-}c^{0}$ |  | ✓ | ✓ |
$a^{-}a^{-}c^{0}$, $a^{-}b^{-}c^{0}$ |  |  | ✓ | ✓
$a^{+}a^{+}a^{+}$, $a^{+}a^{+}c^{+}$, $a^{+}b^{+}c^{+}$ | ✓ |  |  |
$a^{+}a^{+}c^{-}$, $a^{+}b^{+}c^{-}$ |  | ✓ |  |
$a^{+}b^{-}b^{-}$, $a^{+}b^{-}c^{-}$ |  |  | ✓ |
$a^{-}a^{-}a^{-}$, $a^{-}a^{-}c^{-}$, $a^{-}b^{-}c^{-}$ |  |  |  | ✓

The Glazer tilting patterns listed in Table 1 can be divided into four
categories: those with three in-phase tilts, with two in-phase and one out-of-
phase tilts, with one in-phase and two out-of-phase tilts, and with three out-
of-phase tilts. For each of them, the three octahedral-tilting angles as
coordinates span a three-dimensional configuration space (Figure 1b). We
denote them by $a^{+}b^{+}c^{+}$, $a^{+}b^{+}c^{-}$, $a^{+}b^{-}c^{-}$, and
$a^{-}b^{-}c^{-}$, or simply by ppp, ppm, pmm, and mmm, respectively, with “p”
standing for plus (“$+$”) and “m” for minus (“$-$”). For each of them, the
“negative tilt” $(-x)^{r_{x}}$ describes the same tilt as the positive one,
$(+x)^{r_{x}}$ (i.e., $r_{x}$ mode by $\theta_{x}$ angle), but in the opposite
direction [20]. We should note two important features of these configuration
spaces: (a) tilting patterns with zero tilt(s) belong to more than one spaces,
e.g., $a^{+}b^{+}c^{0}$ belongs to both ppp and ppm, and (b) there are
symmetry-equivalent points for each data (except $a^{0}a^{0}a^{0}$). As a
result, any of the eight octants of a configuration space, e.g., octant I with
$\theta_{a}\geqslant 0,\theta_{b}\geqslant 0,\theta_{c}\geqslant 0$, contains
all information for this category. A brief discussion about symmetry images of
data entries that represent different tilting patterns are supplied in Section
S2 in SI.

Figure 1c illustrates another way to construct the 3D disorder configuration
space. All tilting patterns listed in Table 1 are included in one space by
assigning the positive and negative values of each tilting-angle coordinate to
in-phase and out-of-phase tilts, respectively. We thus refer to it as the
“total configuration space” hereafter. The color correspondence between Figure
1b and c indicates that each of the eight octants of this total space can
represent a tilting-pattern category. For example, points in the mmm octant,
i.e., $\theta_{a}\leqslant 0,\theta_{b}\leqslant 0,\theta_{c}\leqslant 0$, of
this total space and points in the octant I of the mmm space represent
identical structures. In addition, either ppm or pmm has two symmetry-
equivalent octants in the total space.

###  II.2 BOSS protocol

BOSS is an ML method that accelerates structure search via strategic sampling
of the PES. In general, it can be used to rapidly explore $N$-dimensional
domains and build surrogate models for any simulated properties. A more in-
depth description of the BOSS scheme can be found in Refs. [76, 77, 54]. Here
we only outline the search principle of BOSS with highlighting the aspects
that are different to earlier applications.

First, the structures were sampled separately in the ppp, ppm, pmm, and mmm
search spaces to ensure that there are enough data in all four tilting-pattern
categories. The 3D surrogate models for visualization and analysis, in
contrast, were constructed within the total space to avoid ambiguous data at
the subspaces (planes and lines) shared by different tilting-pattern
categories. Second, for each arrangement of tilting angles sampled, we
performed a full structural relaxation of all other DOFs. We opted for this
approach instead of the building-block models in previous BOSS applications to
release any strains or structural artefacts introduced by the rearrangement of
perovskite octohedra. Third, we guaranteed that the surrogate model was
constructed and data acquisition was performed based on a set of DFT data
points that reflect the structural symmetry of materials as alluded to. To
this end, the sampling query additionally returned the data of all symmetry-
equivalent points in the search space after each DFT calculation.

The BOSS workflow in this study is as follows. First, DFT calculations of
tilted structures were separately sampled in each of the ppp, ppm, pmm, and
mmm search spaces:

  1. 1.

The database was initialized with three non-equivalent structures together
with their symmetry images;

  2. 2.

A 3D surrogate model was fitted to all existing data using Gaussian process
regression;

  3. 3.

Based on the surrogate model, an acquisition function was calculated and used
to determine the next sample point in the configuration space;

  4. 4.

The database was updated with the newly sampled point and its symmetry images;

  5. 5.

Steps 2 – 4 were repeated until the parameters of the most-stable structure
did not change within 10 consecutive iterations (variation of total energy
within $0.5\leavevmode\nobreak\ \text{meV}$ per perovskite unit and all
tilting angles within $0.5\text{\textdegree}$).

Then the sampled data were merged into one data set according to the symmetry
of the total configuration space. 3D surrogate models were fitted to both DFT-
calculated total energy and band gap data to build the PES and band gap
landscape, respectively. The surrogate model for the band gap was obtained in
the same way from the same set of computations.

Once the PES was evaluated, we analyzed it to extract the locations of all
local minima and their predicted energies. For all minima, we removed the
fixed tilting angle restriction and optimized the structures fully to evaluate
the quality of this small approximation. For any pair of minima, we utilized
the surrogate model to compute the minimum energy path (MEP) between them and
determine the energies and structures of transition states.

For data sampling, we used the exploratory Lower Confidence Bound [78]
acquisition function which balances exploitation (search for better points
near the predicted global minimum) and exploration (acquiring data in the less
visited regions of the search space). Such heuristic approaches are suitable
for an unbiased exploration of the entire configuration space, not only
limited to the vicinity of extrema.

For each sampled structure, we carried out DFT calculations for
$\text{CsPbX}_{3}$ using a $2\times 2\times 2$ supercell model that can host
all Glazer tilting patterns. The Pb–X framework of the system was initialized
according to the desired tilting pattern. The Cs locations and lattice
parameters were relaxed during structural optimization, while the Pb and X
fractional coordinates were frozen so that the tilting pattern was (at least
approximately) maintained. We chose the Perdew-Burke-Ernzerhof exchange-
correlation functional for solids (PBEsol) [79] and analytical stress tensor
[80] implemented in the all-electron numeric-atom-centered orbital code FHI-
aims [81, 82, 83, 84], as PBEsol describes the lattice constants of halide
perovskites well with moderate computational cost [85, 86, 87]. Single-point
band gap calculations were performed for each optimized structure using a
hybrid functional PBEsol0 that contains $25\%$ exact exchange. Scalar
relativistic effects were included by means of the zero-order regular
approximation [88], while spin-orbit coupling (SOC) was further included in
PBEsol0 calculations. Standard FHI-aims tier-2 basis sets were used in
combination with a $\Gamma$-centered $4\times 4\times 4$ $k$-point mesh. The
results of all relevant calculations of this work are available from the NOMAD
(Novel Materials Discovery) repository [89].

##  III Results and discussion

###  III.1 Octahedral-tilting PESs

The surrogate models for the energy landscapes produced by 3D octahedral
tilting allow us insight into the energetics of disorder in $\text{CsPbI}_{3}$
and $\text{CsPbBr}_{3}\,$. We monitored the convergence of the models with
structural sampling to ensure that the landscapes were well converged and
reliable. Figure 2a and b illustrate important 2D cross-sections of the 3D PES
of $\text{CsPbI}_{3}\,$, which was constructed by merging the first $23$,
$75$, $76$, and $34$ DFT calculations ($809$, $795$, $797$, and $793$
symmetry-unfolded points) in the ppp, ppm, pmm, and mmm search spaces,
respectively. The total-space database consists of $203$ non-equivalent DFT
calculations which are symmetrized into $982$ data points. The 3D PES of
$\text{CsPbBr}_{3}$ (Figure 2c and d) was fitted to $195$ DFT calculations
($937$ data points) collected from $38$, $65$, $70$, and $27$ DFT-calculation
entries ($815$, $799$, $797$, and $793$ data points) in ppp, ppm, pmm, and
mmm, respectively. Some further 2D cross-sections of the 3D PESs of both
benchmark systems can be found in Section S3 in SI.

FIG 2: 2D cross-sections of 3D PESs of (a,b) $\text{CsPbI}_{3}$ and (c,d)
$\text{CsPbBr}_{3}\,$. Total energies, $E_{\text{total}}\,$, are given in
their difference to the $a^{0}a^{0}a^{0}$ structure in meV per
$\text{CsPbX}_{3}$ unit. For each material: (a,c) cross-sections
$a^{+}b^{-}b^{-}$, $a^{+}b_{\text{m}}^{-}c^{-}$, and
$a_{\text{m}}^{+}b^{-}c^{-}$ within pmm (subscript “m” in the Glazer tilting
notation means that the corresponding variable is kept constant at its value
at the PES minimum), with white dots indicating the global $E_{\text{total}}$
minima in the total space; and (b,d) border planes $a^{+}b^{+}c^{0}$ (shared
by ppp and ppm), $a^{+}b^{0}c^{-}$ (shared by ppm and pmm), and
$a^{0}b^{-}c^{-}$ (shared by pmm and mmm), with white dots corresponding to
the lowest-energies in ppp, ppm, and mmm, respectively. The axis names of each
plot are given at the top-left corner, e.g., $\theta_{b}\sim\theta_{a}$
meaning $\theta_{b}$ and $\theta_{a}$ (both in °) for the vertical and
horizontal axes, respectively. All figures have the same $x$ and $y$ scales.
Note the different colormap scales for different materials.

At the outset, it is important to verify the quality of the PES surrogate
models. The $a^{+}b^{-}b^{-}$ and $a^{+}b_{\text{m}}^{-}c^{-}$ cross-sections
of $\text{CsPbI}_{3}$ (Figure 2a-top and middle, respectively) agree with the
previous reported 2D PESs of the same material very well [20]. We note that
$76$ DFT calculations sufficed to converge the 3D PES of pmm in our approach,
while in Ref. [20], $289$ DFT data points sampled over an equispaced grid in a
quadrant of the 2D configuration space supplied a 2D PES. This comparison
demonstrates the validity, accuracy, and the remarkable efficiency of our
full-dimensional (3D) approach aided by Bayesian optimization.

The global minima (including symmetry-equivalent structures) in each search
space are marked by white dots in Figure 2. Table 2 features the tilting angle
locations and PES values of all minima, as well the results of subsequent full
structure relaxations. The total energies derived from the surrogate model
($E_{\text{total}}$) are only $<3.5\leavevmode\nobreak\ \text{meV}$ per
perovskite unit higher than the DFT full-relaxation results
($E_{\text{f.r.}}$), further confirming the accuracy of the 3D PES surrogate
models. This difference is attributed to the constraint of the tilting
patterns in DFT calculations. Since the difference is very small, we conclude
that the three tilting angles are suitable principal DOFs.

TABLE 2: Structural data of global minima of all search spaces (the global minimum of the total space is marked by ∗). Positive and negative tilting angles (in °) indicate in-phase and out-of-phase tilts, respectively. $E_{\text{total}}$ is the total energy evaluated based on the 3D PES surrogate models. Also listed are the total energy of corresponding structures fully relaxed with DFT ($E_{\text{f.r.}}$). Both energies are given in their difference to the $a^{0}a^{0}a^{0}$ structure in meV per $\text{CsPbX}_{3}$ unit. Search space |  ∗pmm∗ | ppp | ppm | mmm
---|---|---|---|---
$\text{CsPbI}_{3}$
$\theta_{a}$ | $12.2$ | $15.1$ | $14.5$ | $-14.8$
$\theta_{b}$ | $-9.1$ | $0{\color[rgb]{1,1,1}\definecolor[named]{pgfstrokecolor}{rgb}{1,1,1}\pgfsys@color@gray@stroke{1}\pgfsys@color<EMAIL_ADDRESS>| $0{\color[rgb]{1,1,1}\definecolor[named]{pgfstrokecolor}{rgb}{1,1,1}\pgfsys@color@gray@stroke{1}\pgfsys@color<EMAIL_ADDRESS>| $-3.5$
$\theta_{c}$ | $-9.1$ | $0{\color[rgb]{1,1,1}\definecolor[named]{pgfstrokecolor}{rgb}{1,1,1}\pgfsys@color@gray@stroke{1}\pgfsys@color<EMAIL_ADDRESS>| $-6.2$ | $0{\color[rgb]{1,1,1}\definecolor[named]{pgfstrokecolor}{rgb}{1,1,1}\pgfsys@color@gray@stroke{1}\pgfsys@color<EMAIL_ADDRESS>
$E_{\text{total}}$ | $-129.2$ | $-91.3$ | $-103.5$ | $-93.2$
Space group | $Pnma$ | $P4/mbm$ | $Cmcm$ | $C2/m$
$E_{\text{f.r.}}$ | $-132.6$ | $-91.5$ | $-104.1$ | $-93.7$
In Figure 2 | (a) | (b)-$a^{+}b^{+}c^{0}$ | (b)-$a^{+}b^{0}c^{-}$ | (b)-$a^{0}b^{-}c^{-}$
$\text{CsPbBr}_{3}$
$\theta_{a}$ | $11.1$ | $14.4$ | $13.4$ | $-14.3$
$\theta_{b}$ | $-8.3$ | $0{\color[rgb]{1,1,1}\definecolor[named]{pgfstrokecolor}{rgb}{1,1,1}\pgfsys@color@gray@stroke{1}\pgfsys@color<EMAIL_ADDRESS>| $0{\color[rgb]{1,1,1}\definecolor[named]{pgfstrokecolor}{rgb}{1,1,1}\pgfsys@color@gray@stroke{1}\pgfsys@color<EMAIL_ADDRESS>| $-2.7$
$\theta_{c}$ | $-8.3$ | $0{\color[rgb]{1,1,1}\definecolor[named]{pgfstrokecolor}{rgb}{1,1,1}\pgfsys@color@gray@stroke{1}\pgfsys@color<EMAIL_ADDRESS>| $-6.2$ | $0{\color[rgb]{1,1,1}\definecolor[named]{pgfstrokecolor}{rgb}{1,1,1}\pgfsys@color@gray@stroke{1}\pgfsys@color<EMAIL_ADDRESS>
$E_{\text{total}}$ | $-93.9$ | $-71.5$ | $-77.5$ | $-72.5$
Space group | $Pnma$ | $P4/mbm$ | $Cmcm$ | $C2/m$
$E_{\text{f.r.}}$ | $-97.3$ | $-71.5$ | $-78.3$ | $-73.2$
In Figure 2 | (c) | (d)-$a^{+}b^{+}c^{0}$ | (d)-$a^{+}b^{0}c^{-}$ | (d)-$a^{0}b^{-}c^{-}$

The lowest-energy structures of both materials exhibit an $a^{+}b^{-}b^{-}$
tilting ($\upgamma$ phase, space group $Pnma$) which is found in pmm (Figure
2a and c). The 2D cross-sections of the 3D PESs in this space are
qualitatively similar for these two materials. The main difference is that the
global minimum of $\text{CsPbI}_{3}$ is lower than $\text{CsPbBr}_{3}$ by
$35.3\leavevmode\nobreak\ \text{meV}$ per unit. In addition, we can observe
less steep low-energy basins with $\text{CsPbBr}_{3}$ than
$\text{CsPbI}_{3}\,$. This can be rationalized by structural consideration:
the Goldschmidt tolerance factor of $\text{CsPbBr}_{3}$ ($0.862$) is larger
than that of $\text{CsPbI}_{3}$ ($0.851$) due to the smaller radius of the
$\text{Br}^{-}$ anion. $\text{CsPbI}_{3}$ therefore has a stronger tendency of
octahedral tilting to stabilize the structure, resulting in more energy gain
compared to the untilted structure and a steeper PES in the vicinity of the
minimum.

Based on the 3D PES surrogate models, we can use the nudged elastic band (NEB)
method to easily evaluate the barrier of the MEPs that connect a global
minimum and its symmetry image. Important results are presented in Table 3
(see Section S4 of SI for more detail). The smallest barriers are $25.7$ and
$16.4\leavevmode\nobreak\ \text{meV}$ for $\text{CsPbI}_{3}$ and
$\text{CsPbBr}_{3}\,$, respectively, associated with the
$a^{+}b^{-}b^{-}\rightarrow a^{+}b^{-}c^{0}\rightarrow a^{+}b^{-}(-b)^{-}$
transition path. The lower barrier of $\text{CsPbBr}_{3}$ is a direct result
of the shallower basin in the vicinity of the global minima. The MEPs linking
$a^{+}b^{-}b^{-}$ and $(-a)^{+}b^{-}b^{-}$ have higher barriers, indicating
that larger energy is required for reversing the direction of an in-phase tilt
than for an out-of-phase tilt. The transition $a^{+}b^{-}b^{-}\rightarrow
a^{+}(-b)^{-}(-b)^{-}$ follows two consecutive MEPs, each associated with the
direction reversal of an out-of-phase tilt. A simultaneous direction reversal
of the two out-of-phase tilts is not favored. Our findings for
$\text{CsPbI}_{3}$ that (a) the lowest barrier is associated with a
transition-state structure with an $a^{+}b^{-}c^{0}$ tilting, and (b) its
energy is lower than what is found in the 2D approach by a few meV (see
Section S4 of SI), are in good agreement with Klarbring’s results [20]. The
advantage of our BOSS approach is that we can directly evaluate the MEP within
the full-dimensional (3D) PES, so that we neither need to calculate several 2D
PESs nor will miss the path with the lowest barrier.

TABLE 3: MEPs connecting the global-minimum structure $a^{+}b^{-}b^{-}$ and its symmetry equivalent structures and the associated barriers.  | Barrier $[\text{meV}]$
---|---
Structure-variation MEP | $\text{CsPbI}_{3}$ | $\text{CsPbBr}_{3}$
$a^{+}b^{-}b^{-}\rightarrow a^{+}b^{-}c^{0}\rightarrow a^{+}b^{-}(-b)^{-}$ | $25.7$ | $16.4$
$a^{+}b^{-}b^{-}\rightarrow a^{0}b^{-}c^{-}\rightarrow(-a)^{+}b^{-}b^{-}$ | $35.9$ |
$a^{+}b^{-}b^{-}\rightarrow a^{0}b^{-}b^{-}\rightarrow(-a)^{+}b^{-}b^{-}$ |  | $23.1$
$a^{+}b^{-}b^{-}\rightarrow a^{+}b^{-}c^{0}\rightarrow a^{+}b^{-}(-b)^{-}\rightarrow a^{+}b^{0}(-c)^{-}\rightarrow a^{+}(-b)^{-}(-b)^{-}$ | $25.7$ | $16.4$

The other three search spaces are naturally of less interest than pmm as they
do not contain the overall most stable structure. For an in-depth insight of
the tilting modes, we nevertheless briefly discuss them based on Figure 2b and
d. Table 2 already indicates that the white dots in the top, middle, and
bottom of Figure 2b or d are the global minima in the ppp, ppm, and mmm search
spaces, respectively. The lowest-energy structure of ppp has an
$a^{+}b^{0}b^{0}$ tilting (space group $P4/mbm$) which is generally considered
as the structure of the recently reported tetragonal $\upbeta$ phase [17, 18,
90]. This structure is very localized in terms of the single in-phase tilting
angle, which is obviously larger than the in-phase tilting angle of the
overall most stable ($a^{+}b^{-}b^{-}$) structure (see Table 2). This
structure also belongs to ppm but is not the most stable one within this
search space. The global minimum of ppm is instead located at
$a^{+}b^{0}c^{-}$ whose energy is lower than the $a^{+}b^{0}b^{0}$ structure
by $12.2\leavevmode\nobreak\ \text{meV}$. It is the transition-state structure
of the $a^{+}b^{-}b^{-}\rightarrow a^{+}(-b)^{-}b^{-}$ MEP as alluded to. The
local minimum of mmm exhibits an $a^{-}b^{-}c^{0}$ tilting with one large
($>10\text{\textdegree}$) and one small ($<5\text{\textdegree}$) out-of-phase
tilt. The barrier of the amplitude-exchange pathway,
$a^{0}b^{-}c^{-}\rightarrow a^{0}c^{-}b^{-}$, is rather low (see Section S4 of
SI). Finally, none of the $a^{-}b^{0}b^{0}$ (space group $I4/mcm$) structures
is a global minimum of any search space. This might explain why it has never
been experimentally observed.

Summarizing Figures 2b and d provides insight of structural disorder in these
two benchmark systems. The in-phase and out-of-phase tilts have different
characters. An in-phase tilt around one lattice vector “excludes” in-phase
tilts around any other lattice vectors; while an out-of-phase tilt can
stabilize the structure with an in-phase tilt, and can co-exist and exchange
energy with another out-of-phase tilt. We can thus conclude that the out-of-
phase tilting mode is more active than the in-phase mode in
$\text{CsPbX}_{3}\,$. The motion of out-of-phase tilts transfers vibrational
energy from one lattice vector to another, giving rise to noticeable dynamical
structural disorder.

###  III.2 Configuration space distribution of structures at finite
temperatures

Based on the total-energy profiles in the full-dimensional (3D) configuration
space, we can analyze the distribution of the structural disorder of
$\text{CsPbI}_{3}$ and $\text{CsPbBr}_{3}$ at finite temperatures. Figure 3
sketches the total energy data listed in Table 2 and facilitates the
prediction of how structures evolve and disorder grows in both systems with
increasing temperature.

FIG 3: Energies of the most stable structures (labeled by A–D) in the four
search spaces (red triangles and blue circles for $\text{CsPbI}_{3}$ and
$\text{CsPbBr}_{3}\,$, respectively).

From the BOSS surrogate model of PES (Figure 2), we calculated the probability
density of finding $\text{CsPbI}_{3}$ or $\text{CsPbBr}_{3}$ in a structure
that is labeled by point $\bm{\theta}$ (a shortened notation of tilting
angles) in the total configuration space at temperature $T$ using the
Boltzmann distribution

$\displaystyle\rho(\bm{\theta},T)$
$\displaystyle\propto\text{e}^{-[E(\bm{\theta})-E_{\min}]/k_{\text{B}}T}\,.$

Here $E(\bm{\theta})$ and $E_{\min}$ denote the energies of structure
$\bm{\theta}$ and the overall most stable structure, respectively.
$k_{\text{B}}$ is the Boltzmann constant. In this work, we have not normalized
the total probability but rather kept the maximum probability density always
at $1$. We can thus set a criterion of $\rho$ ($0.5$ in this work) to help
investigate the distribution of thermal population of disordered structures.
We are especially interested in two features: (a) the region within which
$\rho>0.5$ (the larger the region, the more extensive the disorder), and (b)
whether the $\rho>0.5$ distribution emerges at some particular regions in the
configuration space (meaning noticeable population at the corresponding
tilting patterns). Representative results are illustrated in Figure 4 (see
also Section S5 of SI).

FIG 4: 2D cross-sections of unnormalized distribution probability density of
(a,b) $\text{CsPbI}_{3}$ and (c,d) $\text{CsPbBr}_{3}$ in the search spaces of
tilting angles at different temperatures. Shown for each material: (a,c)
$a^{+}b^{-}b^{-}$ cross-section from pmm, and (b,c) $a^{\pm}b^{\pm}c^{0}$
cross-section from the total configuration space. The thick yellow lines in
each plot mark a relative probability of $0.5\,$.

Figure 4a shows that $\text{CsPbI}_{3}$ is quite strongly bound in its most
stable $a^{+}b^{-}b^{-}$ ($Pnma$, $\upgamma$ phase) structure at
$300\leavevmode\nobreak\ \text{K}$. As temperature increases, the $\rho>0.5$
regions grow gradually and reach the $a^{+}b^{-}c^{0}$ plane first at
$450\leavevmode\nobreak\ \text{K}$ (Figure 4b). This corresponds to the A
$\rightarrow$ B transition in Figure 3, which requires an energy of
$25.7\leavevmode\nobreak\ \text{meV}$. At this point, the structure has not
yet gained enough energy to further travel across a series of
$a^{+}b^{+}c^{-}$ structures (B $\rightarrow$ C in Figure 3), but will rather
be “reflected” by this $a^{+}b^{-}c^{0}$ border back to pmm. It can either go
back to the initial $a^{+}b^{-}b^{-}$ structure, or fall into another
equivalent $a^{+}b^{-}(-b)^{-}$ minimum with the direction of one out-of-phase
tilt reversed around the lattice vector. Our observation of the
$a^{+}b^{-}b^{-}\rightarrow a^{+}b^{-}c^{0}\rightarrow a^{+}b^{-}(-b)^{-}$
mechanism is in good agreement with the previous study of the structural
disorder for the observed $\upbeta$-$\text{CsPbI}_{3}$ by Klarbring [20].

As temperature further increases beyond $450\leavevmode\nobreak\ \text{K}$,
noticeable population of $\text{CsPbI}_{3}$ structure at $a^{+}b^{0}b^{0}$
(Figure 4a) and $a^{-}b^{-}c^{0}$ (Figure 4b) is observed at
$650\leavevmode\nobreak\ \text{K}$. The former corresponds to “C” in Figure 3,
with only one single in-phase tilt while both out-of-phase tilts are largely
suppressed. The latter corresponds to “D” in Figure 3. The structure can
easily fluctuate between in-phase and out-of-phase tilts at this temperature,
too. The direction reversal of the in-phase tilt, $a^{+}\rightarrow(-a)^{+}$,
becomes highly probable at $800\leavevmode\nobreak\ \text{K}$. The motion of
the less active in-phase tilt gives rise to an extremely high level of
structural disorder.

$\text{CsPbBr}_{3}$ exhibits generally similar features (Figure 4c and d), yet
noticeable differences can be observed. Because of the shallower PES basins,
the probability distribution of $\text{CsPbBr}_{3}$ at
$300\leavevmode\nobreak\ \text{K}$ is much broader than $\text{CsPbI}_{3}\,$.
First, $\rho>0.5$ can already be observed at $a^{+}b^{-}c^{0}$ at room
temperature due to a much lower activation energy ($16.4\leavevmode\nobreak\
\text{meV}$). At $400\leavevmode\nobreak\ \text{K}$, the $a^{+}b^{0}b^{0}$
($P4/mbm$) structure is populated. At this temperature, we also observe
$a^{+}b^{0}\leftrightarrow a^{0}b^{-}$ tilting-mode transitions. Finally, at
$450\leavevmode\nobreak\ \text{K}$ the $a^{+}\leftrightarrow(-a)^{+}$ motion
sets in. All of these findings indicate a much more disordered structure in
$\text{CsPbBr}_{3}$ than its iodide analog at the same temperature.

Table 4 summarizes the important phenomena derived from Figure 4 and Figure
S2. Both benchmark systems exhibit a three-stage mechanism of structure
variation as temperature increases. Starting from the most stable
$a^{+}b^{-}b^{-}$ structure, we can first observe the motion of one out-of-
phase tilting mode, $a^{+}b^{-}b^{-}\leftrightarrow a^{+}b^{-}c^{0}$. This
process is closely related to the $\upgamma\rightarrow\upbeta$ phase
transition. In our evaluation, the temperature of this stage is very close to
the phase-transition temperature for $\text{CsPbI}_{3}$ ($450$ vs.
$448\leavevmode\nobreak\ \text{K}$ [17, 18]) but clearly lower than that for
$\text{CsPbBr}_{3}$ ($300$ vs. $361\leavevmode\nobreak\ \text{K}$ [91, 69]).
Further increase of temperature activates the motion of another out-of-phase
tilt and the fluctuation between the out-of-phase and in-phase tilts. Finally,
direction reversal of any mode and switch between any two modes become
prevalent, giving rise to the highly disordered, effectively cubic structure
with pronounced lattice vibration. Our evaluation of this stage’s temperature
is a bit higher than the $\upbeta\rightarrow\upalpha$ phase-transition
temperature for $\text{CsPbBr}_{3}$ ($450$ vs. $403\leavevmode\nobreak\
\text{K}$ [91, 69]) and much higher than that for $\text{CsPbI}_{3}$ ($800$
vs. $533\leavevmode\nobreak\ \text{K}$ [17, 18]). The focus of this work is to
understand the nature of structural disorder of these emergent and important
optoelectronic materials at the atomic scale. The PES database created in this
paper has laid a solid groundwork for the direct simulation of the
$\upgamma\rightarrow\upbeta\rightarrow\upalpha$ phase transition, which is
important but beyond the scope of this work.

TABLE 4: Important tilting patterns: the (approximate) temperature at which they become noticeably populated, and the tilting-mode transition to which they closely correspond. | Relevant transition |  $T$ of emergence $[\text{K}]$
---|---|---
Tilting pattern | of tilting modes | $\text{CsPbI}_{3}$ | $\text{CsPbBr}_{3}$
$a^{+}b^{-}c^{0}$ | $a^{+}b^{-}b^{-}\rightarrow a^{+}b^{-}(-b)^{-}$ | $450$ | $<300$
$a^{+}b^{0}b^{0}$ | $\,\,a^{+}b^{-}c^{0}\rightarrow a^{+}(-b)^{-}c^{0}$ | $650$ | $\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ 400$
$a^{-}b^{0}b^{0}$ | $\,\,\,a^{+}b^{0}b^{0}\rightarrow a^{0}b^{-}a^{0}$ | $650$ | $\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ 400$
$a^{-}a^{-}c^{0}$ | $\,\,\,a^{-}b^{0}b^{0}\rightarrow a^{0}b^{-}a^{0}$ | $800$ | $\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ 450$

As a concluding note, we emphasize that the zero-tilting structure
$a^{0}a^{0}a^{0}$, which is generally regarded as the $\upalpha$ (cubic)
phase, is very rarely populated even at $>800\leavevmode\nobreak\ \text{K}$.
This is simply because it is energetically much less favorable than the most
stable $a^{+}b^{-}b^{-}$ structure. Our results indicate that the cubic phase
must be understood as a dynamical mixture of low-symmetry structures, a result
that was found by a recent DFT and molecular dynamics study for
$\text{CsPbBr}_{3}$ [24]. In this phase, plenty of domains or nano-regions of
low-symmetry structures (such as $a^{+}b^{-}b^{-}$) occur with equal
probability; the octahedra of each nano-region are dynamically fluctuating to
be adaptive to their neighbors, thus giving rise to an effectively cubic
average structure.

###  III.3 Band gap distributed over disordered structures

With the sampled DFT calculations we are able to fit surrogate models for
materials properties of interest other than PESs. The band gap is one of the
most important electronic-structure properties for optoelectronic materials.
Here we focus on how band gaps of the benchmark systems vary with octahedral-
tilting angles.

Figure 5 provides representative 2D cross-sections of 3D surrogate models that
are fitted to the band gap data calculated with DFT (PBEsol0+SOC). For
$\text{CsPbI}_{3}\,$, we observe that the band gap increases monotonically
with tilting angle. While for $\text{CsPbBr}_{3}$ the overall trend is the
same, small deviations are observed when all three tilting angles are small.

FIG 5: Band gaps of (a) $\text{CsPbI}_{3}$ and (b) $\text{CsPbBr}_{3}$ in the
search space of tilting angles (shown are 2D cross-sections $a^{+}b^{-}b^{-}$,
$a^{+}b_{\text{m}}^{-}c^{-}$, and $a^{\pm}b^{\pm}c^{0}$ of the 3D surrogate
models).

It it noteworthy that the cubic $a^{0}a^{0}a^{0}$ structure has an obviously
lower band gap than all minima structures in different search spaces (white
dots in Figure 5) for each benchmark system. This is due to the maximal
overlap of atomic orbitals as a result of the zero tilt, which broadens the
valence and conduction bands thus reducing the gap. As already alluded to,
such a cubic structure is not seen in experiments, because it is too high in
energy. The experimental band gaps of $\text{CsPbI}_{3}$ and
$\text{CsPbBr}_{3}$ are $\sim 1.7$ [92, 93, 17] and $\sim
2.3\leavevmode\nobreak\ \text{eV}$ [69, 94], respectively. DFT band gaps of
all minima structures fall into those regions (see also Section S6 in SI),
indicating that the PBEsol0+SOC approach with $\alpha=0.25$ is a suitable
choice for both materials.

Based on the band gap surrogate model, we calculated the temperature-dependent
average band gap of each system by

$\displaystyle\langle E_{\text{gap}}(T)\rangle$
$\displaystyle=\frac{\int\rho(\bm{\theta},T)E_{\text{gap}}(\bm{\theta})\operatorname{d}\bm{\theta}}{\int\rho(\bm{\theta},T)\operatorname{d}\bm{\theta}},$

where the integrations run over the whole configuration space,
$E_{\text{gap}}(\bm{\theta})$ is the BOSS-predicted band gap at $\bm{\theta}$,
the weighting factor $\rho(\bm{\theta},T)$ is the distribution density
displayed in Figure 4. Results are illustrated in the left column of Figure 6.
The right column of Figure 6 displays the search space-resolved band gap vs.
total energy relationship based on the BOSS surrogate model data.

FIG 6: Temperature-dependence of the average band gaps (left) and the band-gap
vs. total-energy relationship (right) of (a) $\text{CsPbI}_{3}$ and (b)
$\text{CsPbBr}_{3}\,$. In the right plots, only the BOSS surrogate-model data
within the low-energy range are shown. Data from pmm, ppm, mmm, and ppp are
colored in blue, green, red, and yellow, respectively. The tilting pattern of
the data entry corresponding to the local minimum of each tilting-pattern
category is highlighted.

At $0\leavevmode\nobreak\ \text{K}$, both systems are in their
$a^{+}b^{-}b^{-}$ global-minimum structures (with some small-amplitude zero-
point vibrations). The calculated band gaps are accordingly $1.90$ and
$2.33\leavevmode\nobreak\ \text{eV}$ for $\text{CsPbI}_{3}$ and
$\text{CsPbBr}_{3}\,$, respectively. As temperature increases, structural
disorder grows first mainly within pmm (blue in Figure 6-right). The average
band gap decreases as the contribution of ppm (green) increases, since the
band gaps of ppm structures within the low-energy range are generally smaller
than pmm. A further temperature increase enlarges the contribution from both
mmm and ppp to the average band gap. The band gap range of mmm (red) is above
ppp (yellow). For $\text{CsPbI}_{3}$ (Figure 6a), the band gap range of the
union of mmm and ppp approximately superposes ppm. Consequently, the average
band gap stabilizes between $300$ and $700\leavevmode\nobreak\ \text{K}$
($1.81$ and $1.82\leavevmode\nobreak\ \text{eV}$, respectively). While for
$\text{CsPbBr}_{3}$ (Figure 6b) at $T>300\leavevmode\nobreak\ \text{K}$, the
band gap range of mmm$\oplus$ppp is slightly higher than ppm, thus resulting
in a slowly increasing average band gap with temperature ($2.28$ and
$2.32\leavevmode\nobreak\ \text{eV}$ at $300$ and $700\leavevmode\nobreak\
\text{K}$, respectively).

For $\text{CsPbBr}_{3}\,$, the calculated average band gap at room temperature
or above agrees well with experiments (ranging between
$2.25\leavevmode\nobreak\ \text{eV}$ [69] and $2.39\leavevmode\nobreak\
\text{eV}$ [94]). Our data reproduce the experimentally observed trend that
the band gap increases with temperature above $300\leavevmode\nobreak\
\text{K}$ [94]. For $\text{CsPbI}_{3}\,$, our model somehow overestimates the
band gap (reported experimental band gaps are $1.67\leavevmode\nobreak\
\text{eV}$ for the $\upgamma$ phase [17] and $1.73\leavevmode\nobreak\
\text{eV}$ for the $\upalpha$ phase [92, 93]). The accuracy of band gap
estimation certainly depends on the computational settings (e.g., the amount
of exact exchange in a hybrid functional), for which we need benchmark studies
in the future. Based on our full-dimensional PESs, however, we can already
propose valuable strategies to design next-generation optoelectronic
materials. For example, smaller tilting angles are desired to tune the band
gap of photovoltaic material $\text{CsPbI}_{3}$ to the optimal $\sim
1.35\leavevmode\nobreak\ \text{eV}$ according to the Shockley-Queisser theory
[95]. Doping with organic monovalent A-site cations that are larger than
$\text{Cs}^{+}$ would be a solution. While for $\text{CsPbBr}_{3}$ in green-
light emission, we are interested in materials with minimal structural
disorder for high monochromaticity. A-site mixing with cations smaller than
$\text{Cs}^{+}$ and X-site mixing with larger anions would be promising
options.

##  IV Conclusions

This work offers new computational insight into the structural disorder in
benchmark systems $\text{CsPbI}_{3}$ and $\text{CsPbBr}_{3}\,$. With
relatively few first-principles (DFT) calculations sampled with the BOSS
scheme, we have obtained accurate landscapes for both total energy and band
gap, each as a function of octahedral-tilting angles. We could infer how the
distribution of crystal structures evolves within the three-dimensional
octahedral-tilting configuration space with increasing temperature. We then
evaluated the statistical mean of band gap, an important materials property
for optoelectronic application, over the disordered structures. Our work
demonstrates that novel cross-disciplinary machine-learning computational
materials science tools such as BOSS can aid the description of structural
disorder phenomena, which is important to make recommendation for rational
design of advanced functional materials such as halide perovskites.

## Acknowledgments

We thank Jari Järvi, Joakim Löfgren, Jarno Laakso, Azimatu Fangnon, Zhaoxin
Wu, and Nan Zhang for the fruitful discussions. We acknowledge the computing
resources by Xi’an Jiaotong University’s HPC platform, the Hefei Advanced
Computing Center, the Aalto Science-IT project, and the CSC-IT Center for
Science. This work was supported by the National Natural Science Foundation of
China under Grant Nos. 62281330043, 11974268, 21503057, 12111530061, and
51902244, the Fundamental Research Funds for the Central Universities (Grant
Nos. HIT.NSRIF.2017032 and xzy012021025), the China Postdoctoral Science
Foundation (Grant No. 2018M643632), the Natural Science Foundation of Shaanxi
Province of China (Grant No. 2023-YBGY-447, 2019JQ-389), the Natural Sciences
& Engineering Research Council of Canada (NSERC, Grant No. RGPIN-2017-06915),
the European Union’s Horizon 2020 research and innovation program under Grant
Agreement No. 676580 [The Novel Materials Discovery (NOMAD)], and the Academy
of Finland (Grant Nos. 316601, 334532 and 305632).

## References

  * Nayak et al. [2012] P. K. Nayak, G. Garcia-Belmonte, A. Kahn, J. Bisquert, and D. Cahen, Energy Environ. Sci. 5, 6022 (2012).
  * Cairns and Goodwin [2013] A. B. Cairns and A. L. Goodwin, Chem. Soc. Rev. 42, 4881 (2013).
  * Keen and Goodwin [2015] D. A. Keen and A. L. Goodwin, Nature 521, 303 (2015).
  * Rhodes et al. [2019] D. Rhodes, S. H. Chae, R. Ribeiro-Palau, and J. Hone, Nat. Mater. 18, 541 (2019).
  * Dragoe and Bérardan [2019] N. Dragoe and D. Bérardan, Science 366, 573 (2019).
  * Deringer [2020] V. L. Deringer, J. Phys. Energy 2, 041003 (2020).
  * Simonov and Goodwin [2020] A. Simonov and A. L. Goodwin, Nat. Rev. Chem. 4, 657 (2020).
  * Deringer et al. [2021] V. L. Deringer, N. Bernstein, G. Csányi, C. Ben Mahmoud, M. Ceriotti, M. Wilson, D. A. Drabold, and S. R. Elliott, Nature 589, 59 (2021).
  * Ye [2009] Z.-G. Ye, MRS Bullet. 34, 277 (2009).
  * Zhang et al. [2015] S. Zhang, F. Li, X. Jiang, J. Kim, J. Luo, and X. Geng, Prog. Mater. Sci. 68, 1 (2015).
  * Sun and Cao [2014] E. Sun and W. Cao, Prog. Mater. Sci. 65, 124 (2014).
  * Abdi-Jalebi et al. [2018] M. Abdi-Jalebi, Z. Andaji-Garmaroudi, A. J. Pearson, G. Divitini, S. Cacovich, B. Philippe, H. Rensmo, C. Ducati, R. H. Friend, and S. D. Stranks, ACS Energy Lett. 3, 2671 (2018).
  * Correa-Baena et al. [2019] J.-P. Correa-Baena, T. M. Brenner, J. Snaider, S. Sun, X. Li, M. A. Jensen, N. T. P. Hartono, L. Nienhaus, S. Wieghold, J. R. Poindexter, et al., Science 363, 627 (2019).
  * Shin et al. [2017] D. Shin, B. Saparov, and D. B. Mitzi, Adv. Energy Mater. 7, 1602366 (2017).
  * Giraldo et al. [2019] S. Giraldo, Z. Jehl, M. Placidi, V. Izquierdo-Roca, A. Pérez-Rodríguez, and E. Saucedo, Adv. Mater. 31, 1806692 (2019).
  * Tan et al. [2022] S. Tan, B. Yu, Y. Cui, F. Meng, C. Huang, Y. Li, Z. Chen, H. Wu, J. Shi, Y. Luo, et al., Angew. Chem. Int. Ed. 61, e202201300 (2022).
  * Stoumpos and Kanatzidis [2015] C. C. Stoumpos and M. G. Kanatzidis, Acc. Chem. Res. 48, 2791 (2015).
  * Marronnier et al. [2018] A. Marronnier, G. Roma, S. Boyer-Richard, L. Pedesseau, J.-M. Jancu, Y. Bonnassieux, C. Katan, C. C. Stoumpos, M. G. Kanatzidis, and J. Even, ACS Nano 12, 3477 (2018).
  * Jinnouchi et al. [2019] R. Jinnouchi, J. Lahnsteiner, F. Karsai, G. Kresse, and M. Bokdam, Phys. Rev. Lett. 122, 225701 (2019).
  * Klarbring [2019] J. Klarbring, Phys. Rev. B 949, 104105 (2019).
  * Chen et al. [2020] L. Chen, B. Xu, Y. Yang, and L. Bellaiche, Adv. Funct. Mater. 30, 1909496 (2020).
  * Patrick et al. [2015] C. E. Patrick, K. W. Jacobsen, and K. S. Thygesen, Phys. Rev. B 92, 201205 (2015).
  * Yang et al. [2020] R. X. Yang, J. M. Skelton, E. L. da Silva, J. M. Frost, and A. Walsh, J. Phys. Chem. Lett. 152, 024703 (2020).
  * Zhu et al. [2022] X. Zhu, S. Caicedo-Dávila, C. Gehrmann, and D. A. Egger, ACS Appl. Mater. Interfaces 14, 22973 (2022).
  * Kirkpatrick et al. [1983] S. Kirkpatrick, C. G. Gellatt, and M. P. Vecchi, Science 220, 671 (1983).
  * Goedecker et al. [2005] S. Goedecker, W. Hellmann, and T. Lenosky, Phys. Rev. Lett. 95, 055501 (2005).
  * Jones [2015] R. O. Jones, Rev. Mod. Phys. 87, 897 (2015).
  * Bisseling et al. [1986] R. H. Bisseling, R. Kosloff, and D. Kosloff, Comput. Phys. Commun. 39, 313 (1986).
  * Wu et al. [2004] T. Wu, H.-J. Werner, and U. Manthe, Science 306, 2227 (2004).
  * Schleder et al. [2019] G. R. Schleder, A. C. M. Padilha, C. Mera Acosta, M. Costa, and A. Fazzio, J. Phys.: Mater. 2, 032001 (2019).
  * Schmidt et al. [2019] J. Schmidt, M. R. G. Marques, S. Botti, and M. A. L. Marques, npj Comput. Mater. 5, 83 (2019).
  * Himanen et al. [2019] L. Himanen, A. Geurts, A. S. Foster, and P. Rinke, Adv. Sci. 6, 1900808 (2019).
  * Deringer et al. [2019] V. L. Deringer, M. A. Caro, and G. Csányi, Adv. Mater. 31, 1902765 (2019).
  * Hart et al. [2021] G. L. W. Hart, T. Mueller, C. Toher, and S. Curtarolo, Nat. Rev. Mater. 6, 730 (2021).
  * Bartók et al. [2010] A. P. Bartók, M. C. Payne, R. Kondor, and G. Csányi, Phys. Rev. Lett. 104, 136403 (2010).
  * Behler [2014] J. Behler, J. Phys.: Condens. Matter 26, 183001 (2014).
  * d’Avezac and Zunger [2008] M. d’Avezac and A. Zunger, Phys. Rev. B 78, 064102 (2008).
  * Wang et al. [2012] Y. Wang, J. Lv, L. Zhu, and Y. Ma, Comput. Phys. Commun. 183, 2063 (2012).
  * Bhattacharya et al. [2013] S. Bhattacharya, S. V. Levchenko, L. M. Ghiringhelli, and M. Scheffler, Phys. Rev. Lett. 111, 135501 (2013).
  * Yamashita et al. [2018] T. Yamashita, N. Sato, H. Kino, T. Miyake, K. Tsuda, and T. Oguchi, Phys. Rev. Materials 2, 013803 (2018).
  * Jørgensen et al. [2018] M. S. Jørgensen, U. F. Larsen, K. W. Jacobsen, and B. Hammer, J. Phys. Chem. A 122, 1504 (2018).
  * Bisbo and Hammer [2020] M. K. Bisbo and B. Hammer, Phys. Rev. Lett. 124, 086102 (2020).
  * Mortensen et al. [2020] H. L. Mortensen, S. A. Meldgaard, M. K. Bisbo, M.-P. V. Christiansen, and B. Hammer, Phys. Rev. B 102, 075427 (2020).
  * Kaappa et al. [2021] S. Kaappa, E. G. del Río, and K. W. Jacobsen, Phys. Rev. B 103, 174114 (2021).
  * Wanzenböck et al. [2022] R. Wanzenböck, A. Marco, B. Sebastian, B. Florian, J. Carrete, and G. K. H. Madsen, Digital Discovery 1, 703 (2022).
  * Rupp et al. [2012] M. Rupp, A. Tkatchenko, K.-R. Müller, and O. A. von Lilienfeld, Phys. Rev. Lett. 108, 058301 (2012).
  * Botu et al. [2017] V. Botu, R. Batra, J. Chapman, and R. Ramprasad, J. Phys. Chem. C 121, 511 (2017).
  * von Lilienfeld and Burke [2020] O. A. von Lilienfeld and K. Burke, Nature Commun. 11, 4895 (2020).
  * Sauceda et al. [2020] H. E. Sauceda, S. Chmiela, I. Poltavsky, K.-R. Müller, and A. Tkatchenko, in _Machine Learning Meets Quantum Physics_ , edited by K. T. Schütt, S. Chmiela, O. A. von Lilienfeld, A. Tkatchenko, K. Tsuda, and K.-R. Müller (Springer International Publishing, Cham, 2020), p. 277.
  * Vandermause et al. [2020] J. Vandermause, S. B. Torrisi, S. Batzner, Y. Xie, L. Sun, A. M. Kolpak, and B. Kozinsky, npj Comput. Mater. p. 20 (2020).
  * Timmermann et al. [2021] J. Timmermann, Y. Lee, C. G. Staacke, J. T. Margraf, C. Scheurer, and K. Reuter, J. Chem. Phys. 155, 244107 (2021).
  * Unke et al. [2021] O. T. Unke, S. Chmiela, H. E. Sauceda, M. Gastegger, I. Poltavsky, K. T. Schütt, A. Tkatchenko, and K.-R. Müller, Chem. Rev. 121, 10142 (2021).
  * Westermayr et al. [2022] J. Westermayr, S. Chaudhuri, A. Jeindl, O. T. Hofmann, and R. J. Maurer, Digital Discovery 1, 463 (2022).
  * Todorović et al. [2019] M. Todorović, M. U. Gutman, J. Corander, and P. Rinke, npj Comput. Mater. 5, 35 (2019).
  * Fang et al. [2021] L. Fang, E. Makkonen, M. Todorovć, P. Rinke, and X. Chen, J. Chem. Theory Comput. 17, 1955 (2021).
  * Egger et al. [2020] A. T. Egger, L. Hörmann, A. Jeindl, M. Scherbela, V. Obersteiner, M. Todorović, P. Rinke, and O. T. Hofmann, Adv. Sci. 7, 2000992 (2020).
  * Järvi et al. [2020a] J. Järvi, P. Rinke, and M. Todorović, Beilstein J. Nanotechnol. 11, 1577 (2020a).
  * Järvi et al. [2020b] J. Järvi, B. Alldritt, O. Krejčí, M. Todorović, P. Liljeroth, and P. Rinke, Adv. Funct. Mater. 31, 2010853 (2020b).
  * Fangnon et al. [2022] A. Fangnon, M. Dvorak, V. Havu, M. Todorović, J. Li, and P. Rinke, ACS Appl. Mater. Interfaces 14, 12758 (2022).
  * Chen et al. [2022] J. Chen, H. Dong, J. Li, X. Zhu, J. Xu, F. Pan, R. Xu, J. Xi, B. Jiao, X. Hou, et al., ACS Energy Lett. 7, 3685 (2022).
  * [61] C. Zhang, Q. Liao, J. Chen, B. Li, C. Xu, K. Wei, G. Du, Y. Wang, D. Liu, J. Deng, et al., Adv. Mater. pp. (in press, doi:10.1002/adma.202209422) (????).
  * Stuke et al. [2021] A. Stuke, P. Rinke, and M. Todorović, Mach. Learn.: Sci. Technol. 2, 035022 (2021).
  * Laakso et al. [2022] J. Laakso, M. Todorović, J. Li, G.-X. Zhang, and P. Rinke, Phys. Rev. Materials 6, 113801 (2022).
  * Sun et al. [2021] S. Sun, A. Tiihonen, F. Oviedo, Z. Liu, J. Thapa, Y. Zhao, N. T. P. Hartono, A. Goyal, T. Heumueller, C. Batali, et al., Matter 4, 1305 (2021).
  * Löfgren et al. [2022] J. Löfgren, D. Tarasov, T. Koitto, P. Rinke, M. Balakshin, and M. Todorović, ACS Sustain. Chem. Eng. 10, 9469 (2022).
  * Sim et al. [2019] K. Sim, T. Jun, J. Bang, H. Kamioka, J. Kim, H. Hiramatsu, and H. Hosono, Appl. Phys. Rev. 6, 031402 (2019).
  * Glazer [1972] A. M. Glazer, Acta Crystallogr. B28, 3384 (1972).
  * Chung et al. [2012] I. Chung, J.-H. Song, J. Im, J. Androulakis, C. D. Malliakas, H. Li, A. J. Freeman, J. T. Kenney, and M. G. Kanatzidis, J. Am. Chem. Soc. 134, 8579 (2012).
  * Stoumpos et al. [2013] C. C. Stoumpos, C. D. Malliakas, J. A. Peters, Z. Liu, M. Sebastian, J. Im, T. C. Chasapis, A. C. Wibowo, D. Y. Chung, A. J. Freeman, et al., Cryst. Growth Des. 13, 2722 (2013).
  * Gehrmann and Egger [2019] C. Gehrmann and D. A. Egger, Nature Commun. 10, 3141 (2019).
  * Glazer [1975] A. M. Glazer, Acta Crystallogr. A31, 756 (1975).
  * Woodward [1997a] P. M. Woodward, Acta Crystallogr. B53, 32 (1997a).
  * Woodward [1997b] P. M. Woodward, Acta Crystallogr. B53, 44 (1997b).
  * Islam et al. [2013] M. A. Islam, J. M. Rondinelli, and J. E. Spanier, J. Phys.: Condens. Matter 25, 175902 (2013).
  * Xie et al. [2020] N. Xie, J. Zhang, S. Raza, N. Zhang, X. Chen, and D. Wang, J. Phys.: Condens. Matter 32, 315901 (2020).
  * Rasmussen and Williams [2006] C. E. Rasmussen and C. K. I. Williams, _Gaussian Processes for Machine Learning_ (The MIT Press, Cambridge, MA, 2006).
  * Gutman and Corander [2016] M. U. Gutman and J. Corander, J. Mach. Learn. Res 17, 1 (2016).
  * [78] E. Brochu, V. M. Cora, and N. de Freitas, arXiv:1012.2599 (2010).
  * 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, Phys. Rev. Lett. 100, 136406 (2008).
  * Knuth et al. [2015] F. Knuth, C. Carbogno, V. Atalla, V. Blum, and M. Scheffler, Comput. Phys. Commun. 190, 33 (2015).
  * Blum et al. [2009] V. Blum, R. Gehrke, F. Hanke, P. Havu, V. Havu, X. Ren, K. Reuter, and M. Scheffler, Comput. Phys. Commun. 180, 2175 (2009).
  * Havu et al. [2009] V. Havu, V. Blum, P. Havu, and M. Scheffler, J. Comput. Phys. 228, 8367 (2009).
  * Ren et al. [2012] X. Ren, P. Rinke, V. Blum, J. Wieferink, A. Tkatchenko, A. Sanfilippo, K. Reuter, and M. Scheffler, New J. Phys. 14, 053020 (2012).
  * Levchenko et al. [2015] S. V. Levchenko, X. Ren, J. Wieferink, R. Johanni, P. Rinke, V. Blum, and M. Scheffler, Comput. Phys. Commun. 192, 60 (2015).
  * Yang et al. [2017] R. X. Yang, J. M. Skelton, E. L. da Silva, J. M. Frost, and A. Walsh, J. Phys. Chem. Lett. 8, 4720 (2017).
  * Bokdam et al. [2017] M. Bokdam, J. Lahnsteiner, B. Ramberger, T. Schäfer, and G. Kresse, Phys. Rev. Lett. 119, 145501 (2017).
  * Seidu et al. [2021] A. Seidu, M. Dvorak, P. Rinke, and J. Li, J. Chem. Phys. 154, 074712 (2021).
  * van Lenthe et al. [1993] E. van Lenthe, E. J. Baerends, and J. G. Sneijders, J. Chem. Phys. 99, 4597 (1993).
  * [89] See https://dx.doi.org/10.17172/NOMAD/2021.09.15-1.
  * Wang et al. [2019] Y. Wang, M. Ibrahim Dar, L. K. Ono, T. Zhang, M. Kan, Y. Li, L. Zhang, X. Wang, Y. Yang, X. Gao, et al., Science 365, 591 (2019).
  * Hirotsu et al. [2013] S. Hirotsu, J. Harada, M. Iizumi, and K. Gesi, J. Phys. Soc. Jpn. 37, 1393 (2013).
  * Eperon et al. [2014] G. E. Eperon, S. D. Stranks, C. Menelaou, M. B. Johnston, L. M. Herz, and H. J. Snaith, Energy Environ. Sci. 7, 982 (2014).
  * Eperon et al. [2015] G. E. Eperon, G. M. Paternò, R. J. Sutton, A. Zampetti, A. A. Haghighirad, F. Cacialli, and H. Snaith, J. Mater. Chem. A 3, 19688 (2015).
  * Mannino et al. [2020] G. Mannino, I. Deretzis, E. Smecca, A. La Magna, A. Alberti, D. Ceratti, and D. Cahen, J. Phys. Chem. Lett. 11, 2490 (2020).
  * Shockley and Queisser [1961] W. Shockley and H.-J. Queisser, J. Appl. Phys. 32, 510 (1961).

Supporting Information


# Enhancing Channel Estimation in Quantized Systems with a Generative Prior
††thanks: The authors acknowledge the financial support by the Federal
Ministry of Education and Research of Germany in the program of “Souverän.
Digital. Vernetzt.”. Joint project 6G-life, project identification number:
16KISK002.

Benedikt Fesl, Aziz Banna, and Wolfgang Utschick

School of Computation, Information and Technology, Technical University of
Munich, Germany
Email: {benedikt.fesl, aziz.banna<EMAIL_ADDRESS>
###### Abstract

Channel estimation in quantized systems is challenging, particularly in low-
resolution systems. In this work, we propose to leverage a Gaussian mixture
model (GMM) as generative prior, capturing the channel distribution of the
propagation environment, to enhance a classical estimation technique based on
the expectation-maximization (EM) algorithm for one-bit quantization. Thereby,
a maximum a posteriori (MAP) estimate of the most responsible mixture
component is inferred for a quantized received signal, which is subsequently
utilized in the EM algorithm as side information. Numerical results
demonstrate the significant performance improvement of our proposed approach
over both a simplistic Gaussian prior and current state-of-the-art channel
estimators. Furthermore, the proposed estimation framework exhibits
adaptability to higher resolution systems and alternative generative priors.

###### Index Terms:

Generative prior, Gaussian mixture, one-bit quantization, channel estimation,
expectation-maximization.

©This work has been submitted to the IEEE for possible publication. Copyright
may be transferred without notice, after which this version may no longer be
accessible.

##  I Introduction

Utilizing a generative prior for inverse problems in wireless communications
has attracted considerable attention recently due to the high potential of
improved estimation performance, e.g., using GMMs [1], mixture of factor
analyzers (MFAs) [2], variational autoencoders (VAEs) [3], generative
adversarial networks (GANs) [4], or diffusion models (DMs) [5, 6]. These works
promise a significant performance improvement over classical techniques and
learning-based regression neural networks (NNs) by learning the channel
distribution for a whole radio propagation environment, capturing valuable
prior information. A general limitation is the common assumption of infinite
precision analog-to-digital converters (ADCs) at the receiver. The primary
motivation for deploying low-resolution ADCs is the significant improvement in
energy efficiency, particularly crucial in massive multiple-input multiple-
output (MIMO) systems. However, utilizing a generative prior for low-
resolution quantization at the receiver, introducing a challenging inverse
problem due to the pronounced nonlinearity introduced by the low-resolution
ADCs, is not very well studied.

An exception is the work in [7], where conditionally Gaussian latent
generative models are utilized to parameterize a conditional Bussgang
estimator. Although this work offers promising performance gains for channel
estimation, it is limited by the linearity of the parameterized Bussgang
estimator, based on the observation that the Bussgang estimator deteriorates
from the mean square error (MSE)-optimal conditional mean estimator (CME),
especially for higher dimensions [8].

A different approach was discussed in [9, 10], focusing on the MAP
optimization problem instead of the minimum mean square error (MMSE)
estimator, which is solved with an EM algorithm. It was shown that explicit
knowledge of the prior distribution is essential, and the case of a simple
Gaussian prior is discussed, leading to closed-form solutions. However, a
simplistic Gaussian prior considerably underrepresents the channel
distribution of a whole propagation environment.

Contributions: We propose a novel estimation framework for quantized systems,
enhancing the EM algorithm from [9, 10] by utilizing a conditionally Gaussian
latent model as generative prior for learning the channel distribution of a
whole radio propagation environment, providing valuable information for the
subsequent estimation task. Due to the conditional Gaussianity of the
leveraged model, we get an analytic closed-form solution for the corresponding
M-step of the EM algorithm, resulting in reduced computational complexity.
While our study concentrates on a GMM as generative prior and the extreme case
of one-bit quantization, the proposed estimation framework exhibits
adaptability for higher resolution systems and alternative generative priors,
e.g., MFAs or VAEs. We discuss different structural constraints on the
covariances of the GMM, i.e., Toeplitz or circulant structure, resulting in
reduced memory overhead and computational complexity. Numerical results
demonstrate significant performance gains compared to using a simplistic
Gaussian model as prior and state-of-the-art channel estimators in quantized
systems.

Notation: The $i$-th element of a vector $\bm{x}$ is denoted as $[\bm{x}]_{i}$
and its real and imaginary parts are given as $\bm{x}_{\Re}$ and
$\bm{x}_{\Im}$, respectively. The standard Gaussian cumulative distribution
function (CDF) is denoted by
${\operatorname{\Phi}}(x)=\int_{-\infty}^{x}\tfrac{1}{\sqrt{2\pi}}\exp(-\tfrac{x^{2}}{2}){\operatorname{d}}x$.

##  II System and Channel Model

Consider an uplink transmission of $P$ pilot signals from single-antenna
mobile terminals (MTs) to an $N$-antenna base station (BS), operating one-bit
ADCs. The quantized receive signal is
$\bm{R}=Q(\bm{Y})=Q(\bm{h}\bm{a}^{{\operatorname{T}}}+\bm{N})$, where
$\bm{R}=[\bm{r}_{1},\dots,\bm{r}_{P}]\in\mathbb{C}^{N\times P}$, $\bm{Y}$ is
the unquantized receive signal, $\bm{h}\in\mathbb{C}^{N}$ denotes the wireless
channel, and $\bm{a}\in\mathbb{C}^{P}$ is the pilot vector with the entries,
cf. [7],
$[\tilde{\bm{a}}]_{i}=\beta_{i}\exp\left({\operatorname{j}}\frac{\pi}{2P}(i-1)\right),~{}i\in\\{1,\dots,P\\}$,
where $\beta_{i}=\frac{1}{2}+\frac{i-1}{2(P-1)}$ is the amplitude spacing. The
pilot vector is subsequently normalized to fulfill the power constraint
$\|\bm{a}\|_{2}^{2}=P$. Furthermore,
$\bm{N}=[\bm{n}_{1},\dots,\bm{n}_{P}]\in\mathbb{C}^{N\times P}$ is white
Gaussian noise with
$\bm{n}_{i}\sim\mathcal{N}_{\mathbb{C}}(\bm{0},\sigma^{2}\bm{{\operatorname{I}}})$
and $Q$ denotes the complex-valued one-bit quantization function

$\displaystyle
Q(\bm{y})=\frac{1}{\sqrt{2}}\left(\operatorname{sign}(\bm{y}_{\Re})+{\operatorname{j}}\operatorname{sign}(\bm{y}_{\Im})\right).$
(1)

After column-wise vectorization, the system model is

$\displaystyle\bm{r}=Q(\bm{y})=Q(\bm{A}\bm{h}+\bm{n})\in\mathbb{C}^{NP}$ (2)

with $\bm{r}=\operatorname{vec}(\bm{R})$, $\bm{y}=\operatorname{vec}(\bm{Y})$,
$\bm{n}=\operatorname{vec}(\bm{N})$, and
$\bm{A}=\bm{a}\otimes\bm{{\operatorname{I}}}$. By normalizing the channels as
${\operatorname{E}}[\|\bm{h}\|_{2}^{2}]=N$, the signal-to-noise ratio (SNR) is
defined as $\text{SNR}=1/\sigma^{2}$.

We utilize a stochastic-geometric channel model based on the 3rd Generation
Partnership Project (3GPP) spatial channel model [11, 12], where channels are
modeled conditionally Gaussian:
$\bm{h}|\bm{\delta}\sim\mathcal{N}_{\mathbb{C}}(\bm{0},\bm{C}_{\bm{h}|\bm{\delta}})$.
The random vector $\bm{\delta}$ collects the uniform distributed angles of
arrival/departure and path gains of the main propagation clusters between a MT
and the BS. The BS employs a uniform linear array (ULA) such that the spatial
channel covariance matrix is given by

$\bm{C}_{\bm{h}|\bm{\delta}}=\int_{-\pi}^{\pi}\omega(\gamma;\bm{\delta})\bm{t}(\gamma)\bm{t}(\gamma)^{{\operatorname{H}}}{\operatorname{d}}\gamma$
(3)

with
$\bm{t}(\gamma)=[1,{\operatorname{e}}^{{\operatorname{j}}\pi\sin(\gamma)},\dots,{\operatorname{e}}^{{\operatorname{j}}\pi(N-1)\sin(\gamma)}]^{{\operatorname{T}}}$
being the array steering vector for an angle of arrival $\gamma$, and $\omega$
is a power density consisting of a sum of weighted Laplace densities whose
standard deviations describe the angle spread of the propagation clusters
[11]. To create a training dataset $\mathcal{H}$, we generate random angles
and path gains for every data sample $\bm{h}_{t}$, combined in
$\bm{\delta}_{t}$, and then draw the sample as
$\bm{h}_{t}\sim\mathcal{N}_{\mathbb{C}}(\bm{0},\bm{C}_{\bm{h}|\bm{\delta}_{t}})$,
which results in an overall non-Gaussian channel distribution [9842343] of the
training dataset $\mathcal{H}=\\{\bm{h}_{t}\\}_{t=1}^{T_{\text{train}}}.$

##  III EM Algorithm with a Gaussian Prior

In this section, we briefly summarize the EM algorithm for channel estimation
from [9, 10], where the case of a Gaussian prior is explicitly discussed. The
goal is to solve the MAP optimization problem, given as

$\displaystyle\hat{\bm{h}}_{\text{MAP}}=\operatorname*{arg\,max}_{\bm{h}\in\mathbb{C}^{N}}p(\bm{h}|\bm{r})=\operatorname*{arg\,max}_{\bm{h}\in\mathbb{C}^{N}}p(\bm{h},\bm{r}).$
(4)

However, as outlined in [9, 10], the direct optimization of 4 is intractable
since the involved densities are inaccessible even under a known (Gaussian)
prior. Thus, an EM algorithm is utilized, where the unquantized receive signal
$\bm{y}$ is considered as a latent variable. The corresponding joint
likelihood $p(\bm{r},\bm{h},\bm{y})$ is decomposed as
$p(\bm{r},\bm{h},\bm{y})=p(\bm{r}|\bm{y},\bm{h})p(\bm{y}|\bm{h})p(\bm{h})$
where $p(\bm{r}|\bm{y})$ is an indicator function since the quantization
function 1 is deterministic and
$p(\bm{y}|\bm{h})=\mathcal{N}_{\mathbb{C}}(\bm{y};\bm{A}\bm{h},\sigma^{2}\bm{{\operatorname{I}}})$.
In this section, it is further assumed that the prior distribution is a zero-
mean Gaussian, i.e.,
$p(\bm{h})=\mathcal{N}_{\mathbb{C}}(\bm{h};\bm{0},\bm{C}_{\bm{h}})$.

After some reformulations, the E-step of the EM algorithm in the $\ell$-th
iteration is shown to consist of the computation of the conditional
expectation $\hat{\bm{y}}^{(\ell)}=\mathbb{E}[\bm{y}|\bm{r},\bm{h}^{(\ell)}]$,
whose real part can be computed elementwise as, cf. [9, 10],

$\displaystyle\begin{aligned}
[\hat{\bm{y}}^{(\ell)}_{\Re}]_{i}=[\bm{b}^{(\ell)}_{\Re}]_{i}+\frac{\sigma}{\sqrt{2\pi}}\frac{[\bm{r}_{\Re}]_{i}\exp\left(-\frac{1}{\sigma^{2}}[\bm{b}_{\Re}^{(\ell)}]_{i}^{2}\right)}{{\operatorname{\Phi}}\left(\frac{2}{\sigma}[\bm{r}_{\Re}]_{i}[\bm{b}^{(\ell)}_{\Re}]_{i}\right)}\end{aligned}$
(5)

where $\bm{b}^{(\ell)}=\bm{A}\bm{h}^{(\ell)}$. The imaginary part is computed
analogously. Although the Gaussian CDF has no analytic expression, there exist
practicably feasible approximations, e.g., [13].

Furthermore, in the M-step, the updated channel estimate
$\hat{\bm{h}}^{(\ell+1)}$ is computed as [10]

$\displaystyle\hat{\bm{h}}^{(\ell+1)}=\operatorname*{arg\,min}_{\bm{h}\in\mathbb{C}^{N}}\|\bm{A}\bm{h}-\hat{\bm{y}}^{(\ell)}\|_{2}^{2}-\sigma^{2}\log
p(\bm{h}).$ (6)

In the case of a zero-mean Gaussian prior, as discussed above, the
optimization problem 6 has the closed-form analytic solution [9, 10]

$\displaystyle\hat{\bm{h}}^{(\ell+1)}=\left(\bm{A}^{{\operatorname{H}}}\bm{A}+\sigma^{2}\bm{C}_{\bm{h}}^{-1}\right)^{-1}\bm{A}^{{\operatorname{H}}}\hat{\bm{y}}^{(\ell)}.$
(7)

The E- and M-steps are repeated after properly initializing
$\hat{\bm{h}}^{(0)}$ until convergence.

##  IV Enhanced Channel Estimation with a Generative Prior

The prerequisite of the EM algorithm in Section III that the channel is
Gaussian imposes a severe limitation. This becomes particularly evident when
considering the channel distribution of a whole radio propagation environment,
generally being considerably underrepresented by a simplistic Gaussian prior.
Thus, we aim to enhance the presented EM algorithm by learning the channel
distribution with a generative model, providing a much stronger prior. A
variety of generative models exist for implicitly learning the channel
distribution through a NN that produces valid channel samples, e.g., GANs [4]
or DMs [5, 6]; however, the M-step in 6 may no longer be tractable to solve,
rendering such models difficult to use.

To circumvent this problem, we resort to conditionally Gaussian latent
variable models, i.e., GMMs, MFAs, or VAEs. These models facilitate the
inference of a latent variable that, when conditioned on the channel
distribution, produces an analytically tractable conditionally Gaussian model
for the channel. In this study, our emphasis is primarily on GMMs, although
the other models are, in principle, equally suitable.

A GMM is a probability density function (PDF) of the form

$\displaystyle
p(\bm{h})=\sum_{k=1}^{K}p(k)\mathcal{N}_{\mathbb{C}}(\bm{h};\bm{\mu}_{\bm{h}|k},\bm{C}_{\bm{h}|k})$
(8)

with the set of learnable parameters
$\\{p(k),\bm{\mu}_{\bm{h}|k},\bm{C}_{\bm{h}|k}\\}_{k=1}^{K}$ being the mixing
coefficients, means, and covariances of the corresponding mixture components.
The parameters of the GMM are fitted in an offline phase via the EM
algorithm111Note that this EM algorithm has no intended connection to the
presented algorithm in Section III. for a given training dataset $\mathcal{H}$
of channel samples, cf. [14, Ch. 9]. An essential property of GMMs is that for
a given data sample, the responsibility of each component can be computed as,
cf. [14, Ch. 9], $p(k|\bm{h})\propto
p(k)\mathcal{N}_{\mathbb{C}}(\bm{h};\bm{\mu}_{\bm{h}|k},\bm{C}_{\bm{h}|k})$.
The GMM can be described via a discrete latent variable with a categorical
distribution, which conditions on one of the $K$ components [14, Ch. 9] and,
thus, yields a conditionally Gaussian latent variable model. In [7], it was
shown that a zero-mean GMM with $\bm{\mu}_{\bm{h}|k}=\bm{0}$ for all
$k=1,\dots,K$ is sufficient to well approximate a feasible wireless channel
distribution. Since this allows for simplified expressions of the resulting
algorithm, we also use a zero-mean GMM.

Directly plugging the learned GMM distribution into 5 yields no closed-form
solution, requiring, e.g., a gradient descent technique. Solving this non-
convex optimization problem in every iteration leads to high complexity and
the possibility of reaching local minima. Thus, we instead propose to infer
the most responsible GMM component for a given pilot observation and run the
EM algorithm from Section III for the conditional Gaussian prior, leading to a
closed-form analytic solution. An estimate of the GMM components’
responsibility for a given pilot observation $\bm{r}$ from 2 can be computed
as [7]

$\displaystyle\hat{p}(k|\bm{r})=\frac{p(k)\mathcal{N}_{\mathbb{C}}(\bm{r};\bm{0},\bm{C}_{\bm{r}|k})}{\sum_{i=1}^{K}p(i)\mathcal{N}_{\mathbb{C}}(\bm{r};\bm{0},\bm{C}_{\bm{r}|i})}$
(9)

where $\bm{C}_{\bm{r}|k}$ is computed through
$\bm{C}_{\bm{y}|k}=\bm{A}\bm{C}_{\bm{h}|k}\bm{A}+\sigma^{2}\bm{{\operatorname{I}}}$
with the arcsine law [15] as

$\displaystyle\bm{C}_{\bm{r}|k}=\frac{2}{\pi}\left(\arcsin\left(\Re(\bm{R}_{\bm{y}|k})\right)+{\operatorname{j}}\arcsin\left(\Im(\bm{R}_{\bm{y}|k})\right)\right)$
(10)

where
$\bm{R}_{\bm{y}|k}=\operatorname{diag}(\bm{C}_{\bm{y}|k})^{-\frac{1}{2}}\bm{C}_{\bm{y}|k}\operatorname{diag}(\bm{C}_{\bm{y}|k})^{-\frac{1}{2}}$.

Consequently, the most responsible GMM component for a given quantized pilot
observation is the MAP estimate [16]

$\displaystyle k^{\star}=\operatorname*{arg\,max}_{k}\,\hat{p}(k|\bm{r}).$
(11)

When assuming the corresponding most responsible GMM component is the true
prior distribution, the M-step in the $\ell$-th iteration of the EM algorithm
from Section III changes to

$\displaystyle\hat{\bm{h}}^{(\ell+1)}=\operatorname*{arg\,min}_{\bm{h}\in\mathbb{C}^{N}}\|\bm{A}\bm{h}-\hat{\bm{y}}^{(\ell)}\|_{2}^{2}-\sigma^{2}\log
p(\bm{h}|k^{\star}).$ (12)

Due to the Gaussianity of
$p(\bm{h}|k^{\star})=\mathcal{N}_{\mathbb{C}}(\bm{0},\bm{C}_{\bm{h}|k^{\star}})$
by design, the closed form solution is given as

$\displaystyle\hat{\bm{h}}^{(\ell+1)}=\left(\bm{A}^{{\operatorname{H}}}\bm{A}+\sigma^{2}\bm{C}_{\bm{h}|k^{\star}}^{-1}\right)^{-1}\bm{A}^{{\operatorname{H}}}\hat{\bm{y}}^{(\ell)}.$
(13)

The E-step in 5 is unchanged since it does not depend on the prior
distribution. The offline phase, consisting of fitting the GMM once to the
underlying channel distribution of the radio propagation environment and the
online evaluation of the EM algorithm, enhanced by the generative prior, is
summarized in Algorithm 1, labeled GMM-EM.

Algorithm 1 Training and inference of the GMM-EM estimator.

0:

0: $\mathcal{H}=\\{\bm{h}_{t}\\}_{t=1}^{T_{\text{train}}}$, $K$

1: Train zero-mean GMM: $\\{p(k),\bm{C}_{\bm{h}|k}\\}_{k=1}^{K}$

1:

1: $\\{p(k),\bm{C}_{\bm{h}|k}\\}_{k=1}^{K}$, $\bm{r}$, $\sigma^{2}$, $\bm{A}$,
$L_{\text{max}}$, $\gamma$

2: Compute responsibilities $\\{\hat{p}(k|\bm{r})\\}_{k=1}^{K}$ via 9

3: Compute MAP estimate:
$k^{\star}=\operatorname*{arg\,max}_{k}\hat{p}(k|\bm{r})$

4: Initialize with least squares (LS) estimate:
$\hat{\bm{h}}^{(0)}=\bm{A}^{\dagger}\bm{r}$

5: for $\ell=0$ to $L_{\text{max}}$ do

6: E-step: Compute
$\hat{\bm{y}}^{(\ell)}=\mathbb{E}[\bm{y}|\bm{r},\hat{\bm{h}}^{(\ell)}]$ via 5

7: M-step:
$\hat{\bm{h}}^{(\ell+1)}=(\bm{A}^{{\operatorname{H}}}\bm{A}+\sigma^{2}\bm{C}_{\bm{h}|k^{\star}}^{-1})^{-1}\bm{A}^{{\operatorname{H}}}\hat{\bm{y}}^{(\ell)}$

8: if
$\|\hat{\bm{h}}^{(\ell+1)}-\hat{\bm{h}}^{(\ell)}\|_{2}\leq\gamma\|\hat{\bm{h}}^{(\ell)}\|_{2}$
then

9: break

10: end if

11: end for

###  IV-A Structural Constraints

To further incorporate structural properties to the generative prior, imposed
by the antenna structure at the BS, i.e., Toeplitz (GMM-EM toep) or circulant
(GMM-EM circ) covariances, we utilize a structured GMM as outlined in [17].
Thereby, each component’s covariance is described by the decomposition

$\displaystyle\bm{C}_{\bm{h}|k}=\bm{Q}^{{\operatorname{H}}}\operatorname{diag}(\bm{c}_{\bm{h}|k})\bm{Q},$
(14)

where $\bm{Q}$ is a (oversampled) discrete Fourier transform (DFT) matrix and
$\bm{c}_{\bm{h}|k}\in\mathbb{R}^{N}_{+}$. These structural constraints allow
for having fewer model parameters and a lower online complexity overhead of
the resulting estimator due to the usage of fast Fourier transforms (FFTs).

###  IV-B Memory and Complexity Analysis

The number of model parameters is scaling with $\mathcal{O}(KN^{2})$ or
$\mathcal{O}(KN)$ for the GMM with full or structurally constrained
covariances, respectively. For reducing the computational complexity, many
iterative computations in Algorithm 1 can be precomputed since the GMM is
fixed after the offline phase, e.g., the filter in the M-step and the inverses
$\bm{C}_{\bm{r}|k}^{-1}$ for evaluating the responsibilities 9. Thus, the
order of complexity is dominated by matrix-vector products, yielding
$\mathcal{O}(\max(L,K)PN^{2})$, where $L$ is the number of iterations of the
EM algorithm until convergence. We note that the computations of the $K$
responsibilities 9 are trivially parallelizable. Although the order of
complexity is the same for the structurally constrained versions, the number
of floating point operations (FLOPs) is drastically reduced by using FFTs.

##  V Baseline Channel Estimators

In the following, we introduce several theoretical and practical state-of-the-
art baselines that are most related and relevant for evaluating the estimation
performance of the proposed algorithm:

  * •

genie-EM: Assuming genie-knowledge of the covariance matrix
$\bm{C}_{\bm{h}|\bm{\delta}}$ from 3 for each pilot observation allows for
running the EM algorithm with the ground-truth prior by replacing
$\bm{C}_{\bm{h}}$ with $\bm{C}_{\bm{h}|\bm{\delta}}$ in the M-step 7. This
yields a theoretical bound on the best-possible estimation performance when
utilizing the EM approach.

  * •

global-EM: Since the genie-knowledge of the covariance matrix is inaccessible
in a practical scenario, we additionally evaluate the EM algorithm when
assuming the prior is Gaussian; thus, in this case, the sample covariance
matrix of the training dataset, i.e.,
$\hat{\bm{C}}=\frac{1}{T_{\text{train}}}\sum_{t=1}^{T_{\text{train}}}\bm{h}_{t}\bm{h}_{t}^{{\operatorname{H}}}$,
is used in the M-step 7.

  * •

genie-BLMMSE: Similar to the genie-EM approach, we assume genie-knowledge of
the covariance matrix $\bm{C}_{\bm{h}|\bm{\delta}}$ from 3 for each pilot
observation to evaluate a bound on the performance of the linear MMSE
estimator based on the Bussgang decomposition, cf. [7].

  * •

global-BLMMSE: In this case, instead of utilizing the ground-truth covariance
matrix, we utilize the sample covariance matrix for evaluating an
approximation of the Bussgang linear MMSE estimator [7].

  * •

EM-GM-GAMP: In [18], a compressive sensing (CS)-based channel estimator is
proposed, which is a combination of an EM algorithm for approximating the
channel PDF in the (approximately) sparse angular domain and the generalized
approximate message passing (GAMP) algorithm to solve the sparse recovery
problem.

  * •

GMM-BLMMSE: We additionally compare to the parameterized conditional Bussgang
estimator based on the GMM from [7]. Note that learning the channel
distribution via a GMM is similar to the proposed approach; however, we do not
parameterize a linear MMSE estimator.

##  VI Numerical Results

For all data-aided approaches, we utilize the same training dataset consisting
of $T_{\text{train}}=100{,}000$ data samples. As performance measure we
evaluate an estimate of the normalized MSE
$\mathbb{E}[\|\bm{h}-\hat{\bm{h}}\|_{2}^{2}]/\mathbb{E}[\|\bm{h}\|_{2}^{2}]$
using $10{,}000$ unseen channel samples. If not otherwise stated, we utilize a
GMM with $K=64$ full covariances. We set $L_{\text{max}}=1{,}000$ and
$\gamma=10^{-3}$.

Figure 1: MSE performance with one (top) and three (bottom) propagation
cluster, $N=64$ antennas, and $P=16$ pilot observations.

Fig. 1 evaluates the MSE performance over the SNR for one (top) and three
(bottom) propagation clusters with $N=64$ antennas and $P=16$ pilot
observations. It can be observed that genie-EM outperforms genie-BLMMSE,
indicating that the MAP objective is superior over the linear MMSE estimator
based on the Bussgang decomposition, especially in the medium to high SNR
regime. This observation is in accordance with the finding in [8] that the
Bussgang estimator significantly deteriorates from the CME in general,
justifying the adoption of nonlinear estimation approaches. The proposed GMM-
EM is on par with the related GMM-BLMMSE in the low SNR regime; however, in
medium and high SNR, it exhibits a significant performance improvement, even
outperforming the genie-BLMMSE approach, which requires genie knowledge of the
true covariance matrix for each observation. Moreover, all other state-of-the-
art estimators are outperformed with a large gap over the whole range of SNR.
In the case of one propagation cluster, the proposed GMM-EM approach even
reaches the performance of genie-EM. In contrast, the gap is generally higher
for more propagation clusters where the channel distribution is less
structured. The Toeplitz-structured GMM is on par with the full GMM, whereas
there is a larger gap for the circulant case, imposing a stronger
approximation.

Figure 2: Average number of EM iterations until convergence for one
propagation cluster, $N=64$ antennas, and $P=16$ pilot observations.

In Fig. 2, we analyze the necessary number of EM iterations until convergence
for the setup of one propagation cluster in Fig. 1. The number of iterations
is particularly low in the low SNR regime, i.e., below ten iterations, for all
considered variants utilizing the EM algorithm. This is a convenient property
since the considered operating range of low-resolution systems is the low SNR
regime where the capacity is not severely reduced [19].

Figure 3: MSE performance for varying numbers of pilot observations for one
propagation cluster, $N=64$ antennas, and $\text{SNR}=\qty{5}{dB}$.

Fig. 3 investigates the MSE performance for a varying number of pilot
observations for $N=64$ antennas and an SNR of $\qty{5}{dB}$. Since the EM-GM-
GAMP is consistently worse than global-EM, and GMM-BLMMSE suffers from a too
high complexity in the large pilot regime, we leave out these comparisons. It
can be observed that the proposed GMM-EM approach is almost on par with the
genie-aided variant over the whole range. In conclusion, the number of pilots
can be drastically reduced when using the proposed approach compared to a
simplistic Gaussian prior without performance losses. This, in turn, results
in an increased data rate, better energy efficiency, and a lower latency of
the channel estimation.

Finally, Fig. 4 evaluates the MSE performance for varying numbers of GMM
components of the GMM-EM approach for one propagation cluster, $N=64$
antennas, and $P=16$ pilot observations. The performance drastically improves
from $K=1$ (equal to the global-EM variant) to $K=16$. Afterward, a saturation
occurs where even more GMM components only lead to marginal improvements.

Figure 4: MSE performance for different numbers of GMM components for one
propagation cluster, $N=64$ antennas, and $P=16$ pilot observations.

##  VII Conclusion and Outlook

We have presented a framework for significantly enhancing a classical channel
estimation algorithm for quantized systems based on a conditionally Gaussian
latent generative prior, i.e., a GMM. Inferring the component with the highest
responsibility for a given pilot observation allows for having a strong prior
together with a low-complexity implementation due to closed-form analytic
solutions of the involved optimization problems of the EM algorithm.
Furthermore, we outlined how structurally constrained GMM covariances reduce
the memory and complexity overhead.

In future work, we aim to extend the presented framework to higher resolution
quantization [9] and different generative priors based on conditionally
Gaussian latent models, e.g., VAEs and MFAs. Additionally, we consider
learning the generative prior from quantized pilot observations as training
data [7].

## References

  * [1] M. Koller, B. Fesl, N. Turan, and W. Utschick, “An asymptotically MSE-optimal estimator based on Gaussian mixture models,” _IEEE Trans. Signal Process._ , vol. 70, pp. 4109–4123, 2022.
  * [2] B. Fesl, N. Turan, and W. Utschick, “Low-rank structured MMSE channel estimation with mixtures of factor analyzers,” in _57th Asilomar Conf. Signals, Syst., Comput._ , 2023, pp. 375–380.
  * [3] M. Baur, B. Fesl, and W. Utschick, “Leveraging variational autoencoders for parameterized MMSE estimation,” 2024, arXiv preprint: 2307.05352.
  * [4] E. Balevi, A. Doshi, A. Jalal, A. Dimakis, and J. G. Andrews, “High dimensional channel estimation using deep generative networks,” _IEEE J. Sel. Areas Commun._ , vol. 39, no. 1, pp. 18–30, 2021.
  * [5] M. Arvinte and J. I. Tamir, “MIMO channel estimation using score-based generative models,” _IEEE Trans. Wireless Commun._ , vol. 22, no. 6, pp. 3698–3713, 2023.
  * [6] B. Fesl, M. Baur, F. Strasser, M. Joham, and W. Utschick, “Diffusion-based generative prior for low-complexity MIMO channel estimation,” 2024, arXiv preprint: 2307.05352.
  * [7] B. Fesl, N. Turan, B. Böck, and W. Utschick, “Channel estimation for quantized systems based on conditionally Gaussian latent models,” _IEEE Trans. Signal Process._ , vol. 72, pp. 1475–1490, 2024.
  * [8] B. Fesl, M. Koller, and W. Utschick, “On the mean square error optimal estimator in one-bit quantized systems,” _IEEE Trans. Signal Process._ , vol. 71, pp. 1968–1980, 2023.
  * [9] A. Mezghani, F. Antreich, and J. A. Nossek, “Multiple parameter estimation with quantized channel output,” in _Int. ITG Workshop Smart Antennas (WSA)_ , 2010, pp. 143–150.
  * [10] C. Stöckle, J. Munir, A. Mezghani, and J. A. Nossek, “Channel estimation in massive MIMO systems using 1-bit quantization,” in _IEEE 17th Int. Workshop Signal Process. Advances Wireless Commun. (SPAWC)_ , 2016.
  * [11] 3GPP, “Spatial channel model for multiple input multiple output (MIMO) simulations,” 3rd Generation Partnership Project (3GPP), Tech. Rep. 25.996 (V16.0.0), Jul. 2020.
  * [12] D. Neumann, T. Wiese, and W. Utschick, “Learning the MMSE channel estimator,” _IEEE Trans. Signal Process._ , vol. 66, no. 11, pp. 2905–2917, Jun. 2018.
  * [13] P. Borjesson and C.-E. Sundberg, “Simple approximations of the error function Q(x) for communications applications,” _IEEE Trans. Commun._ , vol. 27, no. 3, pp. 639–643, 1979.
  * [14] C. M. Bishop, _Pattern Recognition and Machine Learning (Information Science and Statistics)_. Springer New York, 2006.
  * [15] G. Jacovitti and A. Neri, “Estimation of the autocorrelation function of complex Gaussian stationary processes by amplitude clipped signals,” _IEEE Trans. Inf. Theory_ , vol. 40, no. 1, pp. 239–245, 1994.
  * [16] N. Turan, B. Fesl, M. Koller, M. Joham, and W. Utschick, “A versatile low-complexity feedback scheme for FDD systems via generative modeling,” _IEEE Trans. Wireless Commun._ , pp. 1–1, 2023.
  * [17] B. Fesl, M. Joham, S. Hu, M. Koller, N. Turan, and W. Utschick, “Channel estimation based on Gaussian mixture models with structured covariances,” in _56th Asilomar Conf. Signals, Syst., Comput._ , 2022, pp. 533–537.
  * [18] J. Mo, P. Schniter, and R. W. Heath, “Channel estimation in broadband millimeter wave MIMO systems with few-bit ADCs,” _IEEE Trans. Signal Process._ , vol. 66, no. 5, pp. 1141–1154, 2018.
  * [19] J. A. Nossek and M. T. Ivrlač, “Capacity and coding for quantized MIMO systems,” in _Proc. Int. Conf. Wireless Commun. Mobile Comput._ , 2006, p. 1387–1392.