# Predicting Oxidation and Spin States by High-Dimensional Neural Networks:
Applications to Lithium Manganese Oxide Spinels

Marco Eckhoff<EMAIL_ADDRESS>Universität Göttingen,
Institut für Physikalische Chemie, Theoretische Chemie, Tammannstraße 6, 37077
Göttingen, Germany.  Knut Nikolas Lausch  Universität Göttingen, Institut für
Physikalische Chemie, Theoretische Chemie, Tammannstraße 6, 37077 Göttingen,
Germany.  Peter E. Blöchl  Technische Universität Clausthal, Institut für
Theoretische Physik, Leibnizstraße 10, 38678 Clausthal-Zellerfeld, Germany.
Universität Göttingen, Institut für Theoretische Physik, Friedrich-Hund-Platz
1, 37077 Göttingen, Germany.  Jörg Behler<EMAIL_ADDRESS>Universität Göttingen, Institut für Physikalische Chemie, Theoretische Chemie,
Tammannstraße 6, 37077 Göttingen, Germany.  Universität Göttingen,
International Center for Advanced Studies of Energy Conversion (ICASEC),
Tammannstraße 6, 37077 Göttingen, Germany.

###### Abstract

Lithium ion batteries often contain transition metal oxides like LixMn2O4
($0\leq x\leq 2$). Depending on the Li content different ratios of
Mn${}^{\text{III}}$ to Mn${}^{\text{IV}}$ ions are present. In combination
with electron hopping the Jahn-Teller distortions of the Mn${}^{\text{III}}$O6
octahedra can give rise to complex phenomena like structural transitions and
conductance. While for small model systems oxidation and spin states can be
determined using density functional theory (DFT), the investigation of
dynamical phenomena by DFT is too demanding. Previously, we have shown that a
high-dimensional neural network potential can extend molecular dynamics (MD)
simulations of LixMn2O4 to nanosecond time scales, but these simulations did
not provide information about the electronic structure. Here we extend the use
of neural networks to the prediction of atomic oxidation and spin states. The
resulting high-dimensional neural network is able to predict the spins of the
Mn ions with an error of only 0.03 $\hbar$. We find that the Mn
e${}_{\text{g}}$ electrons are correctly conserved and that the number of
Jahn-Teller distorted Mn${}^{\text{III}}$O6 octahedra is predicted precisely
for different Li loadings. A charge ordering transition is observed between
280 and 300 K, which matches resistivity measurements. Moreover, the
activation energy of the electron hopping conduction above the phase
transition is predicted to be 0.18 eV deviating only 0.02 eV from experiment.
This work demonstrates that machine learning is able to provide an accurate
representation of both, the geometric and the electronic structure dynamics of
LixMn2O4, on time and length scales that are not accessible by ab initio MD.

Machine Learning, High-Dimensional Neural Networks, PBE0r Local Hybrid Density
Functional, Oxidation States, Spin States, Lithium Manganese Oxide Spinel,
Charge Ordering, Phase Transition, Jahn-Teller Effect, Electron Hopping

##  I Introduction

Battery performance plays a key role in the technological progress of portable
electronic devices.[1, 2] Computer simulations have become an efficient tool
for a target-oriented and economic development of battery materials with
improved properties[3, 4] since they allow to reach a fundamental
understanding of underlying atomistic processes enabling systematic advances.
However, in particular lithium ion batteries pose a challenge for atomistic
simulations as they often contain complex oxides as positive electrode
materials, in which the transition metal ions can have different oxidation
states.[5, 6]

An important example is the commercially used lithium manganese oxide spinel
LixMn2O4 with $0\leq x\leq 2$.[7, 8] LixMn2O4 has a complex electronic
structure including coexisting Jahn-Teller[9] distorted MnIIIO6 octahedra and
undistorted MnIVO6 octahedra.[10, 11, 12, 13] The additional eg electron of
MnIII (t${}_{\mathrm{2g}}^{3}$e${}_{\mathrm{g}}^{1}$) leads to increased Mn-O
distances compared to MnIV (t${}_{\mathrm{2g}}^{3}$e${}_{\mathrm{g}}^{0}$),
i.e., the electronic structure has a large impact on the local geometry.
During charge and discharge of a battery the Li content varies changing the
ratio of MnIII and MnIV ions. In turn, this ratio of the oxidation states
affects the charge ordering of LixMn2O4, i.e., the distribution of MnIII and
MnIV ions, as well as the electrical conductivity, which is based on electron
hopping of the eg electrons among the different Mn species.[14, 15, 16, 17]
Thermally induced fluctuations of the Jahn-Teller distortions in combination
with electron hopping give rise to complex phenomena like structural
transitions, charge ordering transitions, and conductance. Because electrical
conduction is an indispensable property of a positive electrode material, a
detailed understanding of electrical transport in LixMn2O4 is of central
interest.

The disorder of MnIII and MnIV ions resulting from electron hopping as well as
the fluctuations in the spatial orientations of the Jahn-Teller
distortions[12, 18] on average lead to a cubic LixMn2O4 spinel structure
(space group Fd$\overline{3}$m) for temperatures above $\sim$ 290 K in the
composition range $0\leq x\leq 1$.[19, 13] Below about 290 K the LiMn2O4
spinel is orthorhombic (space group Fddd) due to an increased order in the
MnIII/MnIV arrangement and the spatial orientations of the Jahn-Teller
distortions.[20, 13, 12, 21, 18] For $1<x<2$, a tetragonal crystal structure
(space group I41/amd) with $x=2$ coexists with the cubic form.[22, 23]

Using X-ray absorption spectroscopy[24, 25, 26] or electron energy loss
spectroscopy[27, 28, 29] the eg occupancy or valency of the Mn ions can be
determined. However, the limited time and space resolution does not allow for
an investigation of the dynamics at the atomic scale. On the other hand,
oxidation and spin states are readily available in spin-polarized density
functional theory (DFT) calculations. Unfortunately, ab initio molecular
dynamics simulations are computationally too demanding to investigate changes
in the eg electron mobility with temperature, especially during the transition
from the orthorhombic to the cubic crystal structure, because simulations on
nanosecond time scales would be required.[18]

A solution could be offered by machine learning methods. The construction of
atomistic potentials representing the potential energy surfaces of materials
enabling large-scale molecular dynamics (MD) simulations has received a lot of
attention in recent years.[30, 31, 32] High-dimensional neural network
potentials (HDNNP),[33, 34, 35] which are a frequently used example for
machine learning potentials applicable for various materials,[36, 37, 38, 39,
40, 41, 42, 18] enable molecular dynamics simulations up to nanosecond time
scales but they only provide energies and forces while electronic structure
information like oxidation states are not directly accessible. However, this
information is required for a comprehensive understanding of the full
complexity of LixMn2O4.

In comparison to the high interest in the construction of machine learning
potentials, the prediction of further atomic properties related to the
electronic structure by machine learning is not as common, although also this
field is rapidly expanding. Examples are the prediction of atomic charges,[43,
44] electrostatic multipole moments,[45] atomization energies,[46]
polarizabilities, frontier orbital eigenvalues, ionization potentials,
electron affinities, excitation energies,[47, 48] redox potentials,[49]
nuclear chemical shifts, atomic core level excitations,[50] or even quantum
mechanical wavefunctions.[51] Moreover, the prediction of spin-state
splittings and metal-ligand bond distances was successfully performed for
transition metal complexes in their equilibrium configurations.[52, 53, 54,
55]

Here we introduce a high-dimensional neural network (HDNN) for the prediction
of oxidation and spin states solely based on the atomic structure, which is
applicable in the full configuration space of the system’s electronic ground
state. The HDNN is trained on atomic spins obtained from spin-polarized local
hybrid PBE0r-D3[56, 6] DFT calculations. This functional was extensively
benchmarked for LixMnyOz systems and related materials in a previous study
showing that the atomic spins obtained for LixMn2O4 reveal coexisting high-
spin MnIII and MnIV ions in agreement with experiment.[6] The HDNN is able to
complement the energies and forces in HDNNP-driven MD simulations by
electronic structure information to derive oxidation states, spin
configurations, and electron hopping. Using the HDNN we can study the eg
electron mobility of LixMn2O4 which has a major impact on the electrical
conductivity.[16, 17] Due to its relevance for battery applications, a
detailed understanding of the underlying eg electron hopping is of high
interest. Furthermore, oxidation and spin states are of fundamental importance
for reduction-oxidation reactions such as corrosion and electrocatalysis[57]
and also for spin crossovers, for example, in the binding process of oxygen to
hemoglobin.[58]

In this study we focus on the relation between the structure of the MnO6
octahedra and the spins of the Mn central ions, which we represent by a HDNN.
Using the HDNNP developed in a separate study[18] we can perform simulations
of LixMn2O4 below and above the transition temperature from the orthorhombic
to the cubic crystal structure. The HDNN provides the electronic spin
information, which is not considered in the HDNNP yielding just the energy and
forces. Using both, the HDNNP and the HDNN, we can study the transition in
detail structurally and electronically at the atomic scale. The time evolution
of the Jahn-Teller distortions of the MnIIIO6 octahedra as well as the
electron hopping of the eg electrons are studied from femto- to nanosecond
timescales. Furthermore, we investigate the charge order of LixMn2O4 at
different temperatures revealing its complex electronic structure and the
impact of finite size effects.

We start by describing the HDNN method to predict oxidation and spin states.
After providing the computational details, the reference data set of the HDNN
is analyzed with respect to the structure-spin relation, and the quality of
the HDNN approach based on the reference data is discussed. Subsequently, the
predictions are tested for different Li loadings. The discussion of the phase
transition, whose structural details were explained already in a separate
work,[18] is extended by new insights into the electronic information.
Especially, the charge ordering transition is investigated in detail. Finally,
we address the electron hopping as a function of the temperature revealing the
role of the phase transition on conductivity.

##  II Methods

The neural network topology we will use here is closely related to the high-
dimensional neural network potential approach of Behler and Parrinello.[33,
34, 35] However, instead of developing a method to calculate a global property
like the total potential energy of a given atomic structure
$\mathbf{R}=\\{\mathbf{R}_{n}^{m}\\}$, we aim to predict the atomic spin
$S_{n}^{m}$ for each atom $n$ of element $m$ in the system, i.e., the goal is
to obtain the functional relation $S_{n}^{m}(\mathbf{R})$, which is similar to
our previous work on the construction of environment-dependent atomic charges
for including long-range electrostatic energies in HDNNPs.[43, 59]

Like in case of HDNNPs, we cannot simply use Cartesian coordinates as input
for the neural network, because the derived spins must be invariant with
respect to translation and rotation of the system as well as to the
permutation of chemically equivalent atoms. Therefore, the Cartesian
coordinates are transformed to atom-centered symmetry functions (ACSFs)
$G$,[34] which meet all these requirements and describe the local chemical
environments of each atom inside a cutoff radius $R_{\mathrm{c}}$. This cutoff
has to be chosen sufficiently large to include all neighboring atoms, which
are relevant for the value of the spin of the central atom. Due to their many-
body nature, the number of ACSFs per atom is independent of the actual number
of atoms inside the cutoff sphere, which is required for the use as input of a
neural network with a fixed architecture.

For each atom the vector of ACSFs characterizing the geometric environment is
then used as input for an atomic feed-forward neural network. The great
advantage of neural networks is their ability to represent multidimensional
real-valued functions with in principle arbitrary precision.[60] Therefore,
they are well-suited as functional form for $S_{n}^{m}(\mathbf{G}_{n}^{m})$.
We note that the atomic spin only depends on the symmetry function vector of
the corresponding atom which in turn depends on all Cartesian coordinates
inside the cutoff sphere.

Figure 1: Schematic structure of an atomic neural network containing two
hidden layers that can be used to predict the atomic spin $S_{n}^{m}$ as a
function of the local geometric environment described by a vector of symmetry
functions $\mathbf{G}_{n}^{m}$. For clarity the subscript $n$ of the symmetry
functions and atomic spin as well as the superscript $m$ of all named
properties are not shown. This Figure is a visualization of Equation 1.

The schematic structure of an atomic neural network for the calculation of
atomic spins is given in Figure 1. The input layer of an atomic neural network
consists of the vector of $n_{G}$ symmetry functions $\mathbf{G}_{n}^{m}$.
This input is processed via hidden layers – in the present case we use two
hidden layers with $n_{1}$ and $n_{2}$ neurons, respectively – resulting in a
single neuron in the output layer. This output is the atomic spin of the given
atom. All neurons are connected to all neurons of the neighboring layers by
weight parameters $a^{\rho\sigma}_{\mu\nu}$, where $\mu$ and $\nu$ specify the
two connected neurons in layers $\rho$ and $\sigma$. Additionally, all neurons
in the hidden layers and the output neuron are connected to a bias neuron by a
bias weight $b_{\nu}^{\sigma}$, where $\nu$ and $\sigma$ specify the number
and layer of the target neuron. Together, the connecting weights and the bias
weights are the fitting parameters of the neural network.

For the evaluation of each neuron in the hidden layers one proceeds from the
left to the right through the network. First, the values of the neurons in the
previous layer are multiplied by the corresponding connecting weights and
combined, and the bias weight is added. To enable the representation of
nonlinear functions, to this linear combination an activation function $f$ is
applied – in our case a hyperbolic tangent – and the value is propagated to
the next layer. For the output neuron, which finally yields the atomic spin of
the given atom, the procedure is the same with the exception that the
nonlinear activation function is replaced by a linear function. As a
consequence, the atomic spin is given by the functional form of the feed-
forward neural network,

$\displaystyle\begin{split}S^{m}_{n}=&f_{1}^{3}\left\\{b_{1}^{3}+\sum_{k=1}^{n_{2}}a_{k1}^{23}\cdot
f_{k}^{2}\left[b_{k}^{2}+\sum_{j=1}^{n_{1}}a_{jk}^{12}\cdot
f_{j}^{1}\Bigg{(}b_{j}^{1}\right.\right.\\\
&\left.\left.+\sum_{i=1}^{n_{G}}a_{ij}^{01}\cdot
G_{n,i}^{m}\Bigg{)}\right]\right\\}\hbar\ .\end{split}$ (1)

The superscript $m$ is not shown for most quantities although most of them,
like number of hidden layers, neurons per layer, and also the numerical values
of the weights, can be different for each element $m$.

Figure 2: Structure of the HDNN for atomic spin prediction. The atomic
coordinates $\\{\mathbf{R}_{n}^{m}\\}$ are transformed to symmetry functions
$\\{\mathbf{G}_{n}^{m}\\}$. These are processed for each atom $n$ of element
$m$ in the atomic neural network corresponding to the element yielding the
atomic spin $S_{n}^{m}$ of the atom. The colors represent different elements.

For each element an individual atomic neural network is constructed and
replicated as many times as there are atoms of the respective element present
in the system. The entire structure of the HDNN for atomic spin prediction for
a system of $M$ elements and $N_{m}$ atoms, $1<m<M$, per element is summarized
in Figure 2.

The crucial step is the determination of the weights
$\\{a^{\rho\sigma}_{\mu\nu}\\}$ and $\\{b_{\nu}^{\sigma}\\}$ of the atomic
neural networks, for which we need a reference data set containing the atomic
spin values of a representative set of atomic structures. This data set can be
used to iteratively optimize the weights until the neural network can
reproduce the atomic spins in the configuration space spanned by the trained
data with the desired accuracy. In order to check the predictive power of the
HDNN and to reveal possible remaining errors, predictions for test data, which
are not involved in the training process, are compared with known reference
values. Only if the error of the test and the training set are comparable the
predicted spins can be trusted. Additional validation steps can be applied in
analogy to HDNNPs.[61]

A fundamental requirement for the method to work is that the atomic structure
uniquely defines the atomic spins, which is fulfilled if the correct ground
state electronic structure is consistently used in all reference data. The
larger the structural response with respect to changes in the atomic spins the
better the method will work. Because the structure does not change if the
signs of all atomic spins are inverted, we use the absolute value of the
atomic spins in the training process. While in this way the information about
the direction of the atomic spin is lost, all required information to
determine the oxidation states and spin configurations are still available.

##  III Computational Details

The DFT reference calculations were carried out using the Car-Parrinello
Projector Augmented-Wave (CP-PAW) code (version from September 28, 2016)[62,
63] employing the local hybrid PBE0r functional.[56, 6] All settings regarding
the DFT calculations were the same as in our previous work.[18] D3 dispersion
corrections were applied using the DFT-D3 software (version from June 14,
2016).[64, 65] Becke-Johnson damping was used together with the settings given
for the HSE06 functional according to our prior benchmark.[6] The atomic spins
were determined by a projection of the spin density onto the one-center
expansions of the partial waves. The sum of the contributions was calculated
inside a cutoff radius of 1.2 times the covalent radius of the respective
atom.

The HDNN for the atomic spins as well as the HDNNP for the potential energy
surface of LixMn2O4 were constructed using the RuNNer code (version from
December 4, 2018 and August 22, 2019, respectively).[61, 35, 66] The details
of the HDNNP are described in our previous study.[18] The atomic neural
networks of the HDNN for the redox-active element Mn consists of an input
layer with 63 neurons representing the ACSFs, three hidden layers with 20, 15,
and 10 neurons, and an output layer with one neuron. For the elements Li and O
the atomic neural networks contain 15 input neurons, three hidden layers with
4, 3, and 2 neurons, and the output neuron. The symmetry functions, which are
employed for the different elements, are given in the Supplementary Material.
The reference data was split into a training data set containing 90% of the
data and a test data set containing the remaining structures. The details of
the training process[67, 68] including the description of the adapted weight
initialization method[69, 18] are described in the Supplementary Material. The
spin predictions using the final HDNN were performed by the RuNNer code.

Geometry optimizations and MD simulations driven by the HDNNP were performed
using the Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS)
(version from August 7, 2019).[70, 71] In order to use HDNNPs, LAMMPS was
built with the neural network potential package (n2p2) extension (version from
December 9, 2019).[72] The simulations were run in the isothermal-isobaric
($NpT$) ensemble applying the Nosé-Hoover[73, 74] thermostat and barostat at a
pressure of $p=1\,\mathrm{bar}$. The damping constants were set to 0.05 ps and
0.5 ps, respectively, and the simulations were done with a time step of 0.5
fs. The angles of the simulation cell vectors were fixed at 90 ∘. The
simulations were analyzed after an initial equilibration period of 1 ns.

Basin-hopping Monte Carlo[75] simulations were carried out by a self-written
script using LAMMPS and n2p2 for the geometry optimizations employing the
conjugate gradient algorithm. The Monte Carlo displacements were replaced by
short MD runs at high temperatures in order to change the MnIII/MnIV
distribution. This molecular dynamics based basin-hopping Monte Carlo (MDBHMC)
approach enables an efficient sampling of various spatial distributions of the
MnIIIO6 and MnIVO6 octahedra in LixMn2O4. In the MDBHMC simulations of the
LiMn2O4 unit cell 1000 Monte Carlo steps were performed with $NpT$ MD
simulation lengths of 2000 time steps each. The time step in these MD
simulations was set to 1 fs and a temperature of 400 K has been chosen. The
configurations were classified by their potential energy. If a configuration
was already visited during the simulation, the MD temperature was increased by
10 K in the next step up to a maximum value of 550 K to escape more easily
from this region of the configuration space. The Monte Carlo temperature was
set to 200 K and was increased by 100 K if the new structure after the MD
simulation was not accepted, so that the simulation does not get stuck in a
low-lying minimum but still samples preferentially the configuration space at
low temperatures.

##  IV Results and Discussion

###  IV.1 Jahn-Teller Effect

Our aim is to construct a HDNN, which describes the relation
$S_{n}^{m}(\mathbf{R})$, i.e., the atomic spins must be uniquely defined by
the geometric environments of the atoms involving the positions of many atoms.
Therefore, we first have to investigate the relation between the atomic spins
and the structures in the reference data set. For this purpose, we use the
same reference structures that are underlying the HDNNP,[18] i.e., 15228
LixMn2O4 bulk structures, with $0\leq x\leq 2$. As our previous study has
shown, these structures sample the configuration space of LixMn2O4 accurately
up to a temperature of about 500 K. Instead of the atomic energies and forces
used for the HDNNP, now the reference data set contains the absolute values of
the atomic spins calculated as described in Section III as the target property
of the HDNN.

Figure 3: DFT atomic spin values of the Mn ions as a function of the Jahn-
Teller modes (a) $Q_{1}$ (Equation 2) and (b) $\sqrt{Q_{2}^{2}+Q_{3}^{2}}$
(Equation 3) of the corresponding MnO6 octahedra.

In order to show that the spins of the Mn ions in these structures are a
function of the atomic environment, they are plotted in Figure 3 as a function
of the breathing mode and the Jahn-Teller modes of the corresponding MnO6
octahedra

$\displaystyle Q_{1}$
$\displaystyle=\tfrac{1}{\sqrt{3}}\left(d_{x}+d_{y}+d_{z}-3\overline{d}\right)\
,$ (2) $\displaystyle\sqrt{Q_{2}^{2}+Q_{3}^{2}}$
$\displaystyle=\sqrt{\tfrac{1}{2}\left(d_{x}-d_{y}\right)^{2}+\tfrac{1}{6}(2d_{z}-d_{x}-d_{y})^{2}}$
(3)

with the distances $d_{x}$, $d_{y}$, and $d_{z}$ of the opposite O-O corners
and the mean opposite O-O distance $\overline{d}$.[76] Here, $Q_{1}$ is a
measure for the size of the octahedra, while $\sqrt{Q_{2}^{2}+Q_{3}^{2}}$
refers to the Jahn-Teller distortion.

Figures 3 (a) and (b) clearly show two main groups of MnO6 octahedra and some
intermediate data points, which is not surprising as the structures have been
obtained in MD simulations containing also continuous transitions between
octahedra with and without Jahn-Teller distortions due to electron
hopping.[18] One group is located between spin values of about 1.5 to 1.7
$\hbar$ with a mean spin value of 1.58 $\hbar$ corresponding to MnIV ions and
the other group between about 1.9 to 2.1 $\hbar$ with a mean of 2.01 $\hbar$
corresponding to high-spin MnIII ions. The spins values do not exactly match
1.5 and 2.0 $\hbar$ due to the determination by a projection of the spin
density inside a cutoff sphere of predefined size. This might miss
contributions of the atom which are outside the cutoff or include
contributions of neighboring atoms. Still, these values are qualitatively
accurate enough to determine the different groups and their corresponding
oxidation and spin states. Of the different structures included in the
reference set, all Mn ions in $\uplambda$-Mn2O4 correctly belong to the MnIV
group and all Mn ions in Li2Mn2O4 to the MnIII group as expected. For LiMn2O4,
which contains MnIII and MnIV in a ratio of one-to-one, both spins can be
identified.

The MnIIIO6 octahedra yield larger $Q_{1}$ values than the MnIVO6 octahedra
(Figure 3 (a)), i.e., the MnIIIO6 octahedra are above average in size. This is
in agreement with the additional occupation of the antibonding eg orbital with
one electron, which widens the Mn-O bonds. For about
$0.0\,\mathrm{\AA}<Q_{1}<0.2\,\mathrm{\AA}$ both types of octahedra are
present, i.e., it is not possible to reliably distinguish MnIIIO6 and MnIVO6
octahedra by their size. In this interval also the outliers of the two groups
are found. Since they have intermediate values between the two groups, they
can be classified as transition configurations.

Higher values of $\sqrt{Q_{2}^{2}+Q_{3}^{2}}$ indicate a stronger Jahn-Teller
distortion of the octahedra, and for an ideal octrahedron
$\sqrt{Q_{2}^{2}+Q_{3}^{2}}$ equals zero. Figure 3 (b) shows that the MnIVO6
octahedra are less distorted than the MnIIIO6 octahedra. The reason for small
distortions of the MnIVO6 octahedra are thermal fluctuations included in the
reference data set sampled at finite temperatures. On the other hand, the
MnIIIO6 octahedra show clear Jahn-Teller distortions caused by the single eg
electron leading to the higher values of $\sqrt{Q_{2}^{2}+Q_{3}^{2}}$. Still,
like in Figure 3 (a) in a certain interval the structural descriptor
$\sqrt{Q_{2}^{2}+Q_{3}^{2}}$ is not sufficient for an unique identification,
and the relatively low number of structures in the transition region is also
hardly distinguishable by the combination of $Q_{1}$ and
$\sqrt{Q_{2}^{2}+Q_{3}^{2}}$.

Still, we conclude that the two simple structural descriptors $Q_{1}$ and
$\sqrt{Q_{2}^{2}+Q_{3}^{2}}$ allow to distinguish the chemical environments of
MnIII and MnIV in the vast majority of structures. Therefore, in principle it
should be possible to predict the spins of these species using the HDNN. The
employed ACSFs are more sophisticated structural descriptors. They do not only
consider the opposite corners of the MnO6 octahedra, but as many-body
functions they describe for any atom of interest the positions of all
neighboring atoms including Li, O, and Mn within the atomic environment of
radius $12\,a_{0}$. Therefore, they can be expected to resolve the structural
differences of the octahedra in more detail than $Q_{1}$ and
$\sqrt{Q_{2}^{2}+Q_{3}^{2}}$. Nevertheless, we expect that neither simple
descriptors nor ACSFs will allow for the unique assignment of the oxidation
states of all intermediate transition structures, as they represent a
continuum of possible distortions and thus in principle cannot be physically
assigned to one or the other species.

###  IV.2 High-dimensional Neural Network Spin Prediction

Figure 4: Signed errors of the atomic spins for the (a) training set and (b)
test set of the HDNN as a function of the respective DFT data. The color of
the data points refers to the Li content of the structures (blue:
$\tfrac{1}{8}\leq x\leq\tfrac{7}{8}$). Data points in between the two
horizontal black lines have an absolute error smaller than 0.15 $\hbar$.

Having shown that the different spins of the Mn ions are in principle
distinguishable based on the shape of the octahedra, we now aim to represent
these spins by the HDNN. Starting from the available reference set, we used
the atomic spins of 13669 randomly chosen reference structures for training
and the remaining 1559 structures for testing. Therefore, in total 682431
atomic environments are included in the training set and 77911 in the test
set. More detailed information about the reference data can be found in the
Supplementary Material of a separate work.[18]

The atomic spins cover a range from 0 to 0.117 $\hbar$ for the Li and O ions,
i.e., only LiI and O-II ions are present, and from 1.454 to 2.144 $\hbar$ for
the Mn ions, i.e., MnIV and high-spin MnIII ions. The small spins on the Li
and O ions are a numerical result of the method we employ to determine the
reference spins in the DFT calculations, which are not uniquely defined
similar to atomic partial charges. The difference between the DFT spin density
projections and the expected zero spin values, which is up to 0.1 $\hbar$, can
thus be considered as the uncertainty of our target property. However, the
spin difference of 0.5 $\hbar$ for the two oxidation states of the Mn ions is
a physically reasonable difference and suitable as target property for the
HDNN.

The resulting root mean squared error (RMSE) of the spin for the Li and O ions
is 0.010 $\hbar$ for the training structures and 0.010 $\hbar$ for the test
structures indicating that the HDNN representation does not add significant
errors in addition to the DFT uncertainty. For the Mn ions the corresponding
values are 0.028 and 0.030 $\hbar$ corresponding to a small relative error
compared to the difference of 0.5 $\hbar$ between the two oxidation states.
The maximum spin deviation for Mn is 0.408 $\hbar$ in the training set as well
as in the test set. 0.90% of Mn ions in the training data and 1.09% of Mn ions
in the test data have an error greater than 0.15 $\hbar$. The deviations
between the HDNN and DFT atomic spins are shown in Figure 4 (a) and (b) for
both data sets as a function of the DFT atomic spin.

###  IV.3 Lithium Intercalation

The most important process in battery applications of LixMn2O4 is Li
(de)intercalation, which occurs in charging and discharging. The Li content
determines the ratio of MnIII and MnIV ions, and the number of MnIII ions must
be equal to the number of Li ions to ensure overall charge neutrality.
However, the neural network training process did not include any explicit
information about the electrons and no electron or charge conservation
condition has been imposed neither for the HDNNP nor for the present HDNN.
Therefore, in case of the HDNNP all the information must be included
indirectly in the potential energy surface, which in the simulations must only
yield configurations with a correct MnIII to Li ratio. Consequently, the
correct prediction of this ratio is a stringent test for our approach.

Figure 5: Schematic drawing of the LixMn2O4 unit cell with various Li contents
(a)-(i) $0\leq x\leq 1$. Li ions are shown in green, MnIII ions in yellow, and
MnIV ions in purple. Oxygen ions are not shown for clarity. The size of the
ions represents the $z$ position from large (top layer) to small (bottom
layer). The orientation of the Jahn-Teller elongated MnIIIO6 octahedra is
depicted by the orientation of the yellow diamond. Jahn-Teller elongated
MnIIIO6 octahedra oriented in the $z$ direction are shown as yellow squares.
Black lines highlight the superstructure of corner-sharing
$\left(\text{MnO}_{6}\right)_{4}$ tetrahedra.

Figures 5 (a) to (i) show LixMn2O4 unit cell structures with Li contents in
the range $0\leq x\leq 1$ in steps of $\tfrac{1}{8}$ optimized employing the
HDNNP. The spinel structure is based on MnO6 octahedra, which build a
superstructure of $\left(\text{MnO}_{6}\right)_{4}$ tetrahedra. The
$\left(\text{MnO}_{6}\right)_{4}$ tetrahedra share corners forming the
$\uplambda$-Mn2O4 host lattice. For $0<x\leq 1$, Li is placed in tetrahedral
sites, which form a diamond structure. A three-dimensional LiMn2O4
representation is given in the Supplementary Material.[77]

Indeed, for all Li contents the HDNN reveals that the number of MnIII ions is
equal to the number of Li ions. The structures with low Li content show that
the Li ions are located in proximity to the MnIII ions (note the $z$ direction
as well as the periodic boundary conditions). In the range
$0<x\leq\tfrac{1}{2}$ at most one Mn ion per $\left(\text{MnO}_{6}\right)_{4}$
tetrahedral superstructure formed by four MnO6 octahedra is in the MnIII state
in these minimum energy configurations. At higher Li contents at least one and
at most two MnIII ions are found per $\left(\text{MnO}_{6}\right)_{4}$
tetrahedron.

The Jahn-Teller effect leads to elongations of all the MnIIIO6 octahedra which
are aligned in the same orientation. This agrees with a previous theoretical
study applying GGA$+U$ using the PBE exchange-correlation functional.[78] The
reported optimization of a single unit cell led as well to a tetragonal
structure of LiMn2O4 with the conclusion that the Jahn-Teller distortions are
aligned. In contrast, high-resolution diffraction experiments of the LiMn2O4
spinel report a very complex ordering of the Jahn-Teller distorted MnO6
octahedra below the orthorhombic to cubic transition temperature.[10, 11, 12,
13] The orthorhombic structure exhibits five different crystallographic Mn
sites. Three of these sites are predominantly attributed to MnIII and two to
MnIV. The Jahn-Teller distortion of each MnIII site is in a different
direction. However, this structure does not hold the 1:1 ratio of MnIII and
MnIV. As a consequence, eight MnIII sites are occupied by MnIV ions.
Representing the resulting fractional occupations in this structure would
require the use of a much larger system than the $3\times 3\times 1$
supercell. This is not yet feasible, and our simulations employing cell sizes
up to $3\times 3\times 3$ supercells predict an aligned state as minimum.

Figure 6: Reorientation observed during a nudged elastic band calculation in a
cubic Li0.875Mn2O4 unit cell from (a) the initial minimum to (b) the final
minimum. The symbols are explained in the caption of Figure 5.

Next, we investigate the charge reordering during lithium diffusion along a
minimum energy pathway obtained from a nudged elastic band calculation[79, 80]
of a cubic Li0.875Mn2O4 cell in a separate work (Figure 6).[18] While the Li
ion changes from one tetrahedral site to a neighboring one, simultaneously an
eg electron hops from a close-by MnIIIO6 octahedron to another MnIVO6
octahedron. The diffusion of Li ions therefore can result in changes of the
MnIII/MnIV distribution.

Figure 7: $3\times 3\times 1$ “layer” of a $3\times 3\times 3$ Li0.5Mn2O4
supercell, in which one half of the tetrahedral Li sites is fully occupied by
Li ions, while the other half of the tetrahedral Li sites is empty. The other
two “layers” of the $3\times 3\times 3$ supercell are identical but not shown
for clarity. The symbols are explained in the caption of Figure 5.

In order to investigate if the Li-MnIII relation observed for the rather small
unit cell also holds for larger systems with inhomogeneous Li distributions we
constructed a $3\times 3\times 3$ supercell in which all tetrahedral Li sites
are empty ($x=0$) for one half of the supercell and fully occupied by Li
($x=1$) for the other half. This initial structure has then been optimized to
the closest local minimum employing the HDNNP, which prevents a global
redistribution of the Li ions in the system. The result shown in Figure 7
reveals that the final structure matches the expectation that MnIII ions are
only present in the $x=1$ half and the number of MnIII is equal to the number
of Li ions. We conclude that also in larger systems the Li and MnIII ions are
found in close vicinity.

###  IV.4 Orthorhombic to Cubic Transition

Employing the HDNNP we are able to study the atomic details of the
orthorhombic to cubic phase transition on nanosecond timescales. In addition,
the use of the HDNN further allows to investigate the underlying electronic
dynamics of the eg electrons in form of the fluctuations of the atomic spins
of the Mn ions. In this section, we will verify that HDNNP and HDNN yield a
consistent description of the dynamics of the MnIII/MnIV distribution and
prove the conservation of eg electrons during the simulations. This will
identify critical points for the prediction which are in accordance with the
critical structures observed in Section IV.1.

The orthorhombic structure is a consequence of the Jahn-Teller distorted
MnIIIO6 octahedra, which break the cubic symmetry. The transition to a cubic
phase is caused by disorder and reorientations of the MnIIIO6 octahedra.[18]
The individual MnIII octahedra are still Jahn-Teller distorted but if the
dynamics are faster than the observation time or different orientations exist
in the sample, the average result is a cubic structure. For unit cell
simulations, no spatial averaging is given due to the restricted configuration
space but due to fluctuations of the lattice constants, the time average is
still cubic above the transition temperature.[18]

Figure 8: Structural and electronic dynamics of the LiMn2O4 unit cell at 300
K. Panel (a) and (c) plot the distances of the opposite corners $d_{x}$,
$d_{y}$, and $d_{z}$ in two adjacent MnO6 octahedra and (b) and (d) display
the spins $S$ of the corresponding Mn ions of (a) and (c), respectively.

The underlying atomic processes are the dynamic Jahn-Teller effect and
electron hopping, which are shown in Figures 8 (a) and (c) by the opposite O-O
distances $d_{x}$, $d_{y}$, and $d_{z}$ in two adjacent MnO6 octahedra and
Figures 8 (b) and (d) by the spin $S$ of the corresponding Mn ions for a
simulation of the LiMn2O4 unit cell at 300 K. The eg electrons hops are always
to any of the adjacent MnIV sites. We observe an eg electron hop from the Mn
ion in Figure 8 (b) to the adjacent one in Figure 8 (d) at $t\approx 53$ ps.
At the same time also the Jahn-Teller distortion of the Mn ion in Figure 8 (a)
vanishes and a Jahn-Teller distortion of the Mn ion in Figure 8 (c) reemerges.
The orientation of the Jahn-Teller distortion on the new site is the same as
before on the old site. This is expected because the global orientation
defines the preferred orientation which does not change due to a single
electron hop.

At $t\approx 57$ ps, a spike is observed in Figure 8 (d). At first glance,
this seems to be a prediction error but at the same time also the Jahn-Teller
distortion of the MnIIIO6 vanishes for a short time interval of about 50 fs
(not resolved in Figure 8 (c)). Moreover, there are two spikes at $t\approx
14$ ps and $t\approx 39$ ps in Figure 8 (d). These do not reach a spin value
of 2.0 $\hbar$, i.e., they are still a transition state between MnIV and
MnIII. In Figure 8 (c) the distortion of the MnIVO6 octahedron is also little
larger than usual at these times. This is clearly no Jahn-Teller distorted
MnIIIO6 octahedron but it is difficult to classify if the distortion is a
result of large thermal fluctuation or if an electron transfer was about to
happen.

A change of the orientation of the Jahn-Teller distorted MnIIIO6 octahedron is
observed in Figure 8 (c) at $t\approx 55$ ps. At this time, also the global
orientation of the entire system changes from $z$ to $y$ direction. In Figure
8 (d) no change is observable at this time, i.e., the HDNN predicts the
correct oxidation state during this reorientation. Short reorientations of
only the individual MnIIIO6 octahedra happen, for example, for the MnIIIO6
octahedron in Figure 8 (a) at $t\approx 5$ ps to the $x$ direction and for the
MnIIIO6 octahedron in Figure 8 (c) at $t\approx 62$ ps to the $z$ direction.
In both cases, the orientation changes again to the global orientation after
about 0.25 ps. The oxidation state prediction of the HDNN is in both cases
correct, i.e., stays at the same spin value, during these transitions.

These results show that HDNNP and HDNN yield a consistent description of the
dynamics. There might be a few prediction discrepancies for transition
structures but these are in general very difficult to classify as they are in
between proper MnIII and MnIV states. The agreement is good even for very
short processes.

Figure 9: Electronic dynamics of the LiMn2O4 unit cell at (a) 200, (b) 300,
and (c) 400 K represented by the sum of the Mn spins
$\sum_{i=1}^{N_{\mathrm{Mn}}}S_{i}$ for the $N_{\mathrm{Mn}}=16$ Mn ions.

Employing the HDNN we can also calculate the sum of the absolute spin values
of all Mn ions in the unit cell at different temperatures (Figures 9 (a) to
(c)). From ab initio molecular dynamics simulations (see Supplementary
Material) we know that this property fluctuates around a mean value of about
$28.8\,\hbar$ due to the conservation of eg electrons. The value results from
eight MnIII ions with a spin value of about 1.58 $\hbar$ and eight MnIV ions
with a spin value of about 2.01 $\hbar$, which are included in the unit cell
of LiMn2O4. In the ab initio molecular dynamics simulations the number of
electrons is constant by construction. However, there is no condition in case
of the HDNN. Still, we observe in a simulation at 200 K (Figure 9 (a)) that
the conservation of the eg electrons is excellent.

At 300 and 400 K, where electron hopping and a dynamic Jahn-Teller effect is
present, the fluctuations get a little larger but still remain close to the
reasonable value of $28.8\,\hbar$. There are only a few spikes in Figure 9 (b)
which correspond to one eg electron too less or too much. One of these spikes
is at $t\approx 55$ ps and a very short one at $t\approx 53$ ps. At these
times the global orientation of the structure changes and the electron hop
between the two Mn ions shown in Figure 8 happens, respectively.
Reorientations of octahedra and electron hopping can cause these spikes
because the transition structures include vanishing and reemerging Jahn-Teller
distortions, which are very hard to classify. In combination with thermal
distortions this can lead to less accurate predictions during the transition
in that one Jahn-Teller distortion might already be classified as vanished but
the arising one is not yet classified as a Jahn-Teller distortion or vice
versa. Another source of error can be large thermal fluctuations which lead to
classification uncertainties. At 400 K (Figure 9 (c)) the number of these
spikes increases but the dominant value is still in agreement with the
expectation.

###  IV.5 Charge Ordering

Up to now, we used the HDNN to complement the investigation of the structural
transitions by information about the atomic oxidation states and to confirm
that the MD simulations correspond to a reasonable electronic structure. Now
we go a step further and study electronic processes during the phase
transition, which would not be possible with conventional machine learning
potentials.

The charge and Jahn-Teller order are responsible for the orthorhombic crystal
structure below about 290 K. In a separate study[18] we found a peak in the
heat capacity at about 290 K which we assumed to be a charge ordering
transition. Employing the HDNN we can verify this assignment.

A phase transition is identified by a change of an order parameter. To
describe the charge ordering transition we can employ

$\displaystyle
C=\dfrac{1}{N}\sum_{t=0}^{t_{N}}\dfrac{\mathbf{S}_{\mathrm{min}}\cdot\mathbf{S}(t)}{|\mathbf{S}_{\mathrm{min}}|\cdot|\mathbf{S}(t)|}\
,$ (4)

where $\mathbf{S}_{\mathrm{min}}$ and $\mathbf{S}(t)$ are vectors of all Mn
spins in the simulation cell in the minimum configuration and at time $t$,
respectively. Both $\mathbf{S}_{\mathrm{min}}$ and $\mathbf{S}(t)$ correspond
to the difference between the actual spin values and the mean spin value of
all Mn at all times, which is about 1.8 $\hbar$. For normalization we divide
by the length of both vectors. The sum includes all time steps from $t=0$ to
$t_{N}$ at step $N$. For $C=1$ the minimum charge order does not change with
time, for $C=0$ the charge order is not correlated to the minimum charge
order, i.e., the charge order fluctuates fast. The transition from $C=1$ to 0
as a function of the temperature therefore corresponds to the charge ordering
phase transition.

Figure 10: Charge order parameter $C$ obtained from 20 ns $NpT$ MD simulations
at various temperatures of LiMn2O4 using different cell sizes. The black lines
at $C=0$ and $C=1$ represent the ideal charge ordered and disordered phases.

Figure 10 shows the order parameter $C$ as a function of the temperature $T$
for the $2\times 2\times 2$ and $3\times 3\times 3$ supercell of LiMn2O4. The
data were obtained from 20 ns $NpT$ MD simulations. The results of the
$3\times 3\times 3$ supercell show a clear transition between roughly 280 and
300 K. This is in excellent agreement with resistivity experiments which yield
a transition temperature of about 290 K.[16] For the $2\times 2\times 2$
supercell the transition region is between about 260 and 340 K. This
continuous change of the order parameter over a wider temperature range is
expected because of finite size effects.[81, 82] The transition temperature is
shifted to higher temperatures with decreasing system size. In conclusion, our
assumption that the peak in the heat capacity at about 290 K is associated to
a charge ordering transition is verified by the HDNN analysis.

The MnIII/MnIV distribution as well as the orientations of the Jahn-Teller
distortions are disordered in the high temperature cubic structure. However,
one still can find some frequently occurring patterns. These can also be found
in simulations of a single unit cell, but of course in this case the order is
the same for the entire structure. As the charge order of the low temperature
phase in our simulations can also be constructed using the single unit cell we
will concentrate here on the low energy charge orders of the unit cell.

Figure 11: Minimum energy configurations and orientations of the MnIII ions in
the LiMn2O4 unit cell. The potential energy increases in the order from (a) to
(c). The symbols are explained in the caption of Figure 5.

To investigate the energetically most stable MnIII/MnIV distributions in more
detail MDBHMC simulations were performed for a single unit cell of LiMn2O4.
Figures 11 (a) to (c) shows the three energetically most stable structures we
found. The corresponding three charge orders agree with those found in $NpT$
MD simulations from our separate work.[18] The low temperature charge order of
the $2\times 2\times 2$ and $3\times 3\times 3$ supercells are based on the
order in Figure 11 (a). The energy difference between Figure 11 (a) with
lattice constants 8.664, 8.066, and 8.066 $\mathrm{\AA}$ and Figure 11 (b)
with lattice constants 8.631, 8.113, and 8.067 $\mathrm{\AA}$ is 1.4 meV
atom-1. The difference from Figure 11 (a) to Figure 11 (c) with lattice
constants 8.630, 8.089, and 8.089 $\mathrm{\AA}$ is 1.9 meV atom-1. These tiny
energy differences are smaller than the RMSE of the HDNNP for the test set,
which is 2.2 meV atom-1, and also below the accuracy of the underlying DFT
calculations. Therefore, the different charge orders are effectively
degenerate in energy. In PBE0r DFT calculations for a cubic LiMn2O4 cell the
charge order of Figure 11 (b) is observed as global minimum.[6]

The distribution of MnIII and MnIV ions as well as the orientations of the
Jahn-Teller distortions can change during a MD simulation. Therefore, the
different minimum charge orders can be found oriented in different spatial
directions. The charge order of Figure 11 (b) is the same as in Figure 5 (i),
which is not immediately obvious. In Figure 5 (i) the Mn ions are ordered in
an alternating pattern of helically aligned MnIII and MnIV ions with the Jahn-
Teller elongation in a direction perpendicular to the helix. A helix is
formed, for example, by the four MnIII ions on the right half of the LiMn2O4
unit cell in Figure 5 (i). The helix includes four MnIII ions per turn, i.e.,
the pitch is equal to the lattice constant. As the helix is aligned in $z$
direction the size of the diamonds in Figure 5 (i) increases from the lower
right MnIII ion to the upper right MnIII ion clock-wise with the turn of the
helix. In Figure 11 (b) the helices of MnIII and MnIV ions are aligned in $x$
direction, for example, the lower four MnIII ions form a helix at half height
in $z$ direction.

In Figure 11 (a) the Jahn-Teller distorted octahedra are aligned perpendicular
to MnIII planes resulting in a tetragonal cell. For a change of the global
orientation of the Jahn-Teller distortions the MnIII/MnIV order has to change
as already observed in a separate study.[18] For the charge order in Figure 11
(b) two degenerate configurations of the global orientation of the Jahn-Teller
distortions exist. If the MnIII helices are aligned, for example, in $z$
direction, the Jahn-Teller distortions can align in $x$ or $y$ direction.
Moreover, also the parallel alignment to the MnIII helices was observed.
Therefore, the global orientation of the Jahn-Teller distortions can change
without a redistribution of the MnIII/MnIV order. This explains the atomistic
dynamics in MD simulations of a single unit cell reported in a separate
work.[18]

The charge order in Figure 11 (c) seems to be an intermediate state between
the other two charge orders. With only two conversions between MnIII and MnIV
ions each, i.e. two eg electron hops, this charge order can transform to any
of the other two charge orders. The quite small energy differences among these
configurations indicate that modifications in the distributions of the MnIII
and MnIV ions result only in moderate energetic changes.

###  IV.6 Electrical Conductivity

Electrical transport properties are of central interest for every battery
material. Employing the HDNN we can study the underlying electronic processes
at finite temperatures in realistic model systems which are capable to
represent, for example, the charge ordering transition in LiMn2O4.

The electrical conductivity,

$\displaystyle\sigma=e\mu n\ ,$ (5)

is the product of the electron charge $e$ and the mobility $\mu$ and number
density $n$ of the electrons which are responsible for the charge
transport.[83] The mobility is proportional to the diffusion coefficient $D$,

$\displaystyle\mu=\dfrac{eD}{k_{\mathrm{B}}T}\ ,$ (6)

with the Boltzmann constant $k_{\mathrm{B}}$.[84] The diffusion coefficient in
turn is proportional to the electron hopping frequency $\nu$,

$\displaystyle D=\tfrac{1}{6}d^{2}\nu\ ,$ (7)

because the electron hopping occurs always to adjacent MnIV sites with the
nearest-neighbor Mn-Mn distance $d$.[85] This equation is derived for the
limit of low concentrations of diffusing particles. In our case the eg
electrons can only hop to adjacent MnIV sites, which are only half of the Mn
sites. It should be noted that the correlation between the eg electrons can
influence the diffusion, which is not included in this model.

Because electron hopping is a thermally activated process, the hopping
frequency $\nu$ is given by

$\displaystyle\nu=\nu_{0}\cdot\exp\left(-\dfrac{E_{\mathrm{a}}}{k_{\mathrm{B}}T}\right)\
,$ (8)

with the attempt frequency $\nu_{0}$ and the hopping activation energy
$E_{\mathrm{a}}$.[86] The activation energy for the eg electron migration is a
consequence of the associated Jahn-Teller distortion of the MnIIIO6 octahedra
locally trapping the electrons.[14, 15] For an eg electron hop thermal
distortions have to reduce the Jahn-Teller distortion of a MnIIIO6 octahedron
and create a Jahn-Teller like distortion on an adjacent MnIVO6 octahedron.
This enables the eg electron to tunnel via the O ion from one Mn site to
another. Subsequently, the Jahn-Teller distortion on the former MnO6
octahedron vanishes and emerges on the latter one.

Combining Equations 5 to 8 we obtain

$\displaystyle\sigma
T=A\cdot\exp\left(-\dfrac{E_{\mathrm{a}}}{k_{\mathrm{B}}T}\right)\ ,$ (9)

with the pre-exponential factor,

$\displaystyle A=\dfrac{e^{2}d^{2}\nu_{0}n}{6k_{\mathrm{B}}}\ .$ (10)

Experiments show that two thermally activated processes are relevant for the
electrical conductivity of LiMn2O4,[16, 17] which can both be described by the
Arrhenius relation in Equation 9. Both processes are governed by hopping of
small polarons, i.e., the eg electrons of the MnIII ions. The activation
energy of the thermally activated hopping conduction was measured to be 0.16
eV above the phase transition temperature and it is almost the same below the
transition temperature.[16, 10]

To predict the electrical conductivity using the HDNN we have to investigate
the dominant conduction process, which is the eg electron hopping. The hopping
frequency of the eg electron self-diffusion can be analyzed in $NpT$ MD
simulations of LiMn2O4 at various temperatures to obtain the $\nu(T)$ relation
(Equation 8).

Figure 12: Arrhenius relation of the $e_{\mathrm{g}}$ electron hopping
frequency $\nu$ obtained from 20 ns $NpT$ MD simulations at various
temperatures of LiMn2O4 using different cell sizes. The solid lines show
linear fits in different temperature regions of the $3\times 3\times 3$
supercell data.

In Figure 12 the natural logarithm of the $e_{\mathrm{g}}$ electron hopping
frequency $\nu$ due to self-diffusion is plotted as a function of the inverse
temperature $T^{-1}$. Electron hops are counted if the spin value of a MnIII
ion decreases below 1.65 $\hbar$ or the spin value of a MnIV ion increases
above 1.95 $\hbar$ to exclude counting of attempted incomplete transitions
which result in small spikes in $S(t)$. The absolute number of electron hops
has been divided by the number of eg electrons and by the total simulation
time to obtain the hopping frequency $\nu$. The 20 ns MD simulations were
analyzed after each picosecond.

The results obtained for the $3\times 3\times 3$ supercell clearly show two
different activated processes with a transition between 280 and 290 K
($3.45\lesssim T^{-1}/10^{-3}\,\mathrm{K}^{-1}\lesssim 3.57$). This agrees
very well with our investigation in Section IV.5 and the experimentally
observed transition at about 290 K.[16] The activation energy of the low
temperature process (250 to 270 K according to the order parameter in Figure
10) is 0.19 eV, while it is 0.18 eV for the high temperature process (310 to
400 K). This means that the transformation from the charge ordered to the
charge disordered phase reduces slightly the effective energy barrier for eg
electron diffusion among the Mn sites.

The resulting value of $\nu_{0}$ obtained from the 20 ns $3\times 3\times 3$
supercell data of the hopping frequency $\nu$ for the high temperature process
is $\nu_{0}=4\cdot 10^{13}\,\mathrm{s}^{-1}$. This matches the experimental
expectation[17] because for strongly correlated oxides this frequency is in
the range from $10^{13}$ to $10^{14}\,\mathrm{s}^{-1}$.[87] For the low
temperature process the attempt frequency is also in the order of
$10^{13}\,\mathrm{s}^{-1}$.

To obtain the pre-exponential factor $A$ (Equation 10), we need the averaged
nearest-neighbor Mn-Mn distance in the temperature region from 310 to 400 K,
which is $d=2.93\,\mathrm{\AA}$. Moreover, $n$ equals the number of eg
electrons per LiMn2O4 unit cell, which is 8, divided by the cell volume with
an averaged cubic lattice constant of 8.27 $\mathrm{\AA}$ in the investigated
temperature regime. With these data the pre-exponential factor $A$ is
predicted to be $2\cdot 10^{5}\,\Upomega^{-1}\,\mathrm{cm}^{-1}\,\mathrm{K}$.

With respect to the experiment[16] the conductivity is overestimated by our
model by two orders of magnitude, while it should be noted that the
experimentally measured conductivity changes by about four orders of magnitude
between 250 and 400 K. The overestimation is probably the consequence of the
neglected correlation of the electrons which leads to a lower diffusion
coefficient as predicted by Equation 7. Due to occupied sites the reverse
motion is more probable than another direction because the former site of the
eg electron is always empty shortly after the hop. This reduces the diffusion
coefficient. Still, the activation energy of the hopping process can be
determined precisely.

For the smaller $2\times 2\times 2$ supercell the phase transition is not as
pronounced as for the $3\times 3\times 3$ supercell. If the hopping frequency
data of the $2\times 2\times 2$ and the $3\times 3\times 3$ supercells are
compared, it can be observed that they match above about 340 K ($T^{-1}\approx
2.94\cdot 10^{-3}\,\mathrm{K}^{-1}$) which is in agreement with the results of
the order parameter (Figure 10).

In conclusion, large simulation cells and nanosecond simulation times are
required for an accurate identification of the phase transition based on eg
electron hopping. The efficiency of the HDNN method for the prediction of
oxidation and spin states fulfills these requirements enabling investigations
of the electronic dynamics.

##  V Conclusion

In this work we have trained a HDNN to the atomic spins of the redox-active Mn
ions in LixMn2O4, which could be reproduced with an error of only about 0.03
$\hbar$ with respect to the underlying DFT reference data. Using this HDNN we
have been able to predict the oxidation and spin states of the Mn ions in MD
simulations obtained using a HDNNP potential energy surface. The results
confirmed that the MD simulations in our separate study[18] correspond to
reasonable electronic structures in which the Mn eg electrons are correctly
conserved and the number of Jahn-Teller distorted MnIIIO6 octahedra has been
predicted reliably for different Li loadings.

Further, the HDNN enables investigations of a charge ordering transition and
electron hopping, which are not possible with conventional machine learning
potentials. The charge order parameter reveals a phase transition between 280
and 300 K. This agrees well with experimental resistivity measurements
observing a transition at about 290 K.[16] The hopping frequency reveals this
transition between 280 and 290 K as well. The predicted activation barrier of
the thermally activated conduction process for the high temperature phase of
0.18 eV deviates only 0.02 eV from experimental data.[16]

In conclusion, we have demonstrated that machine learning is able to provide
an accurate representation of both, the structural and the electronic dynamics
of LixMn2O4. The investigation of the phase transition requires time and
length scales that are not accessible by first-principles calculations.
Therefore, the HDNN is a promising new method for the theoretical description,
for example, of electrical transport properties of battery materials.
Moreover, the HDNN method can be applied for molecular transition metal
complexes as well. In these systems high-spin and low-spin states can occur,
which lead to very different metal-ligand distances. Since for a given
distance either the high-spin or the low-spin state is the electronic ground
state, our method can easily identify the spin states of these systems as
well.

## Supplementary Material

See Supplementary Material for (I) symmetry functions, (II) neural network
training details, and (III) a three-dimensional LiMn2O4 representation and ab
initio molecular dynamics simulations.

###### Acknowledgements.

This project was funded by the Deutsche Forschungsgemeinschaft (DFG, German
Research Foundation) - 217133147/SFB 1073, project C03. We gratefully
acknowledge the funding of this project by computing time provided by the
Paderborn Center for Parallel Computing (PC2) and by the DFG project
INST186/1294-1 FUGG (Project No. 405832858). J.B. is grateful for a DFG
Heisenberg professorship BE3264/11-2 (Project No. 329898176).

## Data Availability

The data that support the findings of this study are available from the
corresponding author upon reasonable request.

## References

## References

  * Tarascon and Armand [2001] J.-M. Tarascon and M. Armand, “Issues and challenges facing rechargeable lithium batteries,” Nature 414, 359–367 (2001).
  * Goodenough and Kim [2010] J. B. Goodenough and Y. Kim, “Challenges for rechargeable Li batteries,” Chem. Mater. 22, 587–603 (2010).
  * Hautier  _et al._ [2010] G. Hautier, C. C. Fischer, A. Jain, T. Mueller, and G. Ceder, “Finding nature’s missing ternary oxide compounds using machine learning and density functional theory,” Chem. Mater. 22, 3762–3767 (2010).
  * Sanchez-Lengeling and Aspuru-Guzik [2018] B. Sanchez-Lengeling and A. Aspuru-Guzik, “Inverse molecular design using machine learning: Generative models for matter engineering,” Science 361, 360–365 (2018).
  * Nitta  _et al._ [2015] N. Nitta, F. Wu, J. T. Lee, and G. Yushin, “Li-ion battery materials: present and future,” Mater. Today 18, 252–264 (2015).
  * Eckhoff, Blöchl, and Behler [2020] M. Eckhoff, P. E. Blöchl, and J. Behler, “Hybrid density functional theory benchmark study on lithium manganese oxides,” Phys. Rev. B 101, 205113 (2020).
  * Thackeray  _et al._ [1983] M. M. Thackeray, W. I. F. David, P. G. Bruce, and J. B. Goodenough, “Lithium insertion into manganese spinels,” Mat. Res. Bull. 18, 461–472 (1983).
  * Thackeray [1997] M. M. Thackeray, “Manganese oxides for lithium batteries,” Prog. Solid State Chem. 25, 1–71 (1997).
  * Jahn and Teller [1937] H. A. Jahn and E. Teller, “Stability of polyatomic molecules in degenerate electronic states,” Proc. Royal Soc. Lond. A 161, 220–235 (1937).
  * Rodríguez-Carvajal  _et al._ [1998] J. Rodríguez-Carvajal, G. Rousse, C. Masquelier, and M. Hervieu, “Electronic crystallization in a lithium battery material: Columnar ordering of electrons and holes in the spinel LiMn2O4,” Phys. Rev. Lett. 81, 4660–4663 (1998).
  * Massarotti  _et al._ [1999] V. Massarotti, D. Capsoni, M. Bini, P. Scardi, M. Leoni, V. Baron, and H. Berg, “LiMn2O4 low-temperature phase: synchrotron and neutron diffraction study,” J. Appl. Cryst. 32, 1186–1189 (1999).
  * Piszora [2004] P. Piszora, “Temperature dependence of the order and distribution of Mn3+ and Mn4+ cations in orthorhombic LiMn2O4,” J. Alloy. Comp. 382, 112–118 (2004).
  * Akimoto  _et al._ [2004] J. Akimoto, Y. Takahashi, N. Kijima, and Y. Gotoh, “Single-crystal X-ray structure analysis of the low temperature form of LiMn2O4,” Solid State Ion. 172, 491–494 (2004).
  * Schütte, Colsmann, and Reuter [1979] L. Schütte, G. Colsmann, and B. Reuter, “Kristallographische, elektronische und magnetische Eigenschaften des Spinells Li[Mn2]O4,” J. Solid State Chem. 27, 227–231 (1979).
  * Goodenough, Manthiram, and Wnetrzewski [1993] J. B. Goodenough, A. Manthiram, and B. Wnetrzewski, “Electrodes for lithium batteries,” J. Power Sources 43, 269–275 (1993).
  * Shimakawa, Numata, and Tabuchi [1997] Y. Shimakawa, T. Numata, and J. Tabuchi, “Verwey-type transition and magnetic properties of the LiMn2O4 spinels,” J. Solid State Chem. 131, 138–143 (1997).
  * Iguchi, Nakamura, and Aoki [1998] E. Iguchi, N. Nakamura, and A. Aoki, “Electrical transport properties in LiMn2O4,” Philos. Mag. B 78, 65–77 (1998).
  * Eckhoff  _et al._ [2020] M. Eckhoff, F. Schönewald, M. Risch, C. A. Volkert, P. E. Blöchl, and J. Behler, “Closing the gap between theory and experiment for lithium manganese oxide spinels using a high-dimensional neural network potential,” (2020), arXiv:2007.00327 [cond-mat.mtrl-sci] .
  * Takahashi  _et al._ [2003] Y. Takahashi, J. Akimoto, Y. Gotoh, K. Dokko, M. Nishizawa, and I. Uchida, “Structure and electron density analysis of lithium manganese oxides by single-crystal X-ray diffraction,” J. Phys. Soc. Jpn. 72, 1483–1490 (2003).
  * Akimoto  _et al._ [2000] J. Akimoto, Y. Takahashi, Y. Gotoh, and S. Mizuta, “Single crystal X-ray diffraction study of the spinel-type LiMn2O4,” Chem. Mater. 12, 3246–3248 (2000).
  * Ouyang, Shi, and Lei [2009] C. Y. Ouyang, S. Q. Shi, and M. S. Lei, “Jahn-Teller distortion and electronic structure of LiMn2O4,” J. Alloy. Comp. 474, 370–374 (2009).
  * Mosbah, Verbaere, and Tournoux [1983] A. Mosbah, A. Verbaere, and M. Tournoux, “Phases LixMnO2 rattachees au type spinelle,” Mat. Res. Bull. 18, 1375–1381 (1983).
  * Ohzuku, Kitagawa, and Hirai [1989] T. Ohzuku, M. Kitagawa, and T. Hirai, “Electrochemistry of manganese dioxide in lithium nonaqueous cell,” J. Electrochem. Soc. 136, 3169–3174 (1989).
  * de Groot [2001] F. de Groot, “High-resolution X-ray emission and X-ray absorption spectroscopy,” Chem. Rev. 101, 1779–1808 (2001).
  * de Groot [2005] F. de Groot, “Multiplet effects in X-ray spectroscopy,” Coord. Chem. Rev. 249, 31–63 (2005).
  * de Groot and Kotani [2008] F. de Groot and A. Kotani,  _Core Level Spectroscopy of Solids_ (CRC Press, Boca Raton, 2008).
  * Varela  _et al._ [2009] M. Varela, M. P. Oxley, W. Luo, J. Tao, M. Watanabe, A. R. Lupini, S. T. Pantelides, and S. J. Pennycook, “Atomic-resolution imaging of oxidation states in manganites,” Phys. Rev. B 79, 085117 (2009).
  * Zhang  _et al._ [2010] S. Zhang, K. J. T. Livi, A.-C. Gaillot, A. T. Stone, and D. R. Veblen, “Determination of manganese valence states in (Mn3+, Mn4+) minerals by electron energy-loss spectroscopy,” Am. Mineral. 95, 1741–1746 (2010).
  * Schönewald  _et al._ [2020] F. Schönewald, M. Eckhoff, M. Baumung, M. Risch, P. E. Blöchl, J. Behler, and C. A. Volkert, “A criticial view on eg occupancy as a descriptor for oxygen evolution catalytic activity in LiMn2O4 nanoparticles,” (2020), arXiv:2007.04217 [cond-mat.mtrl-sci] .
  * Behler [2016] J. Behler, “Perspective: Machine learning potentials for atomistic simulations,” J. Chem. Phys. 145, 170901 (2016).
  * Deringer, Caro, and Csanyi [2019] V. L. Deringer, M. A. Caro, and G. Csanyi, “Machine learning interatomic potentials as emerging tools for materials science,” Adv. Mater. 31, 1902765 (2019).
  * Noé  _et al._ [2020] F. Noé, A. Tkatchenko, K.-R. Müller, and C. Clementi, “Machine learning for molecular simulation,” Ann. Rev. Phys. Chem. 71, 361–390 (2020).
  * Behler and Parrinello [2007] J. Behler and M. Parrinello, “Generalized neural-network representation of high-dimensional potential-energy surfaces,” Phys. Rev. Lett. 98, 146401 (2007).
  * Behler [2011] J. Behler, “Atom-centered symmetry functions for constructing high-dimensional neural network potentials,” J. Chem. Phys. 134, 074106 (2011).
  * Behler [2017] J. Behler, “First principles neural network potentials for reactive simulations of large molecular and condensed systems,” Angew. Chem. Int. Ed. 56, 12828–12840 (2017).
  * Behler  _et al._ [2008] J. Behler, R. Martoňák, D. Donadio, and M. Parrinello, “Metadynamics simulations of the high-pressure phases of silicon employing a high-dimensional neural network potential,” Phys. Rev. Lett. 100, 185501 (2008).
  * Artrith, Hiller, and Behler [2013] N. Artrith, B. Hiller, and J. Behler, “Neural network potentials for metals and oxides - first applications to copper clusters at zinc oxide,” Phys. Status Solidi B 250, 1191–1203 (2013).
  * Gastegger and Marquetand [2015] M. Gastegger and P. Marquetand, “High-dimensional neural network potentials for organic reactions and an improved training algorithm,” J. Chem. Theory Comput. 11, 2187–2198 (2015).
  * Morawietz  _et al._ [2016] T. Morawietz, A. Singraber, C. Dellago, and J. Behler, “How van der Waals interactions determine the unique properties of water,” PNAS 113, 8368–8373 (2016).
  * Hellström and Behler [2016] M. Hellström and J. Behler, “Concentration-dependent proton transfer mechanisms in aqueous NaOH solutions: From acceptor-driven to donor-driven and back,” J. Phys. Chem. Lett. 7, 3302–3306 (2016).
  * Natarajan and Behler [2016] S. K. Natarajan and J. Behler, “Neural network molecular dynamics simulations of solid-liquid interfaces: water at low-index copper surfaces,” Phys. Chem. Chem. Phys. 18, 28704–28725 (2016).
  * Eckhoff and Behler [2019] M. Eckhoff and J. Behler, “From molecular fragments to the bulk: Development of a neural network potential for MOF-5,” J. Chem. Theory Comput. 15, 3793–3809 (2019).
  * Artrith, Morawietz, and Behler [2011] N. Artrith, T. Morawietz, and J. Behler, “High-dimensional neural-network potentials for multicomponent systems: Applications to zinc oxide,” Phys. Rev. B 83, 153101 (2011).
  * Ghasemi  _et al._ [2015] S. A. Ghasemi, A. Hofstetter, S. Saha, and S. Goedecker, “Interatomic potentials for ionic systems with density functional accuracy based on charge densities obtained by a neural network,” Phys. Rev. B 92, 045131 (2015).
  * Houlding, Liem, and Popelier [2007] S. Houlding, S. Y. Liem, and P. L. A. Popelier, “A polarizable high-rank quantum topological electrostatic potential developed using neural networks: Molecular dynamics simulations on the hydrogen fluoride dimer,” Int. J. Quantum Chem. 107, 2817–2827 (2007).
  * Rupp  _et al._ [2012] M. Rupp, A. Tkatchenko, K.-R. Müller, and O. A. von Lilienfeld, “Fast and accurate modeling of molecular atomization energies with machine learning,” Phys. Rev. Lett. 108, 058301 (2012).
  * Montavon  _et al._ [2013] G. Montavon, M. Rupp, V. Gobre, A. Vazquez-Mayagoitia, K. Hansen, A. Tkatchenko, K.-R. Mueller, and O. A. von Lilienfeld, “Machine learning of molecular electronic properties in chemical compound space,” New J. Phys. 15, 095003 (2013).
  * Li  _et al._ [2018] Z. Li, N. Omidvar, W. S. Chin, E. Robb, A. Morris, L. Achenie, and H. Xin, “Machine-learning energy gaps of porphyrins with molecular graph representations,” J. Phys. Chem. A 122, 4571–4578 (2018).
  * Janet  _et al._ [2020] J. P. Janet, S. Ramesh, C. Duan, and H. J. Kulik, “Accurate multiobjective design in a space of millions of transition metal complexes with neural-network-driven efficient global optimization,” ACS Cent. Sci. 6, 513–524 (2020).
  * Rupp, Ramakrishnan, and von Lilienfeld [2015] M. Rupp, R. Ramakrishnan, and O. A. von Lilienfeld, “Machine learning for quantum mechanical properties of atoms in molecules,” J. Phys. Chem. Lett. 6, 3309 (2015).
  * Schütt  _et al._ [2019] K. Schütt, M. Gastegger, A. Tkatchenko, K.-R. Müller, and R. Maurer, “Unifying machine learning and quantum chemistry with a deep neural network for molecular wavefunctions,” Nature Comm. 10, 5024 (2019).
  * Janet and Kulik [2017a] J. P. Janet and H. J. Kulik, “Predicting electronic structure properties of transition metal complexes with neural networks,” Chem. Sci. 8, 5137–5152 (2017a).
  * Janet and Kulik [2017b] J. P. Janet and H. J. Kulik, “Resolving transition metal chemical space: Feature selection for machine learning and structure–property relationships,” J. Phys. Chem. A 121, 8939–8954 (2017b).
  * Janet, Chan, and Kulik [2018] J. P. Janet, L. Chan, and H. J. Kulik, “Accelerating chemical discovery with machine learning: Simulated evolution of spin crossover complexes with an artificial neural network,” J. Phys. Chem. Lett. 9, 1064–1071 (2018).
  * Taylor  _et al._ [2020] M. G. Taylor, T. Yang, S. Lin, A. Nandy, J. P. Janet, C. Duan, and H. J. Kulik, “Seeing is believing: Experimental spin states from machine learning model structure predictions,” J. Phys. Chem. A 124, 3286–3299 (2020).
  * Sotoudeh  _et al._ [2017] M. Sotoudeh, S. Rajpurohit, P. Blöchl, D. Mierwaldt, J. Norpoth, V. Roddatis, S. Mildner, B. Kressdorf, B. Ifland, and C. Jooss, “Electronic structure of Pr1-xCaxMnO3,” Phys. Rev. B 95, 235150 (2017).
  * McCafferty [2010] E. McCafferty,  _Introduction to Corrosion Science_ (Springer, New York, 2010).
  * Perutz  _et al._ [1998] M. F. Perutz, A. J. Wilkinson, M. Paoli, and G. G. Dodson, “The stereochemical mechanism of the cooperative effects in hemoglobin revisited,” Annu. Rev. Biophys. Biomol. Struct. 27, 1–34 (1998).
  * Morawietz, Sharma, and Behler [2012] T. Morawietz, V. Sharma, and J. Behler, “A neural network potential-energy surface for the water dimer based on environment-dependent atomic energies and charges,” J. Chem. Phys. 136, 064103 (2012).
  * Hornik, Stinchcombe, and White [1989] K. Hornik, M. Stinchcombe, and H. White, “Multilayer feedforward networks are universal approximators,” Neural Netw. 2, 359–366 (1989).
  * Behler [2015] J. Behler, “Constructing high-dimensional neural network potentials: A tutorial review,” Int. J. Quantum Chem. 115, 1032–1050 (2015).
  * Blöchl [1994] P. E. Blöchl, “Projector augmented-wave method,” Phys. Rev. B 50, 17953–17979 (1994).
  * Blöchl [2016] P. E. Blöchl, “CP-PAW,” https://www2.pt.tu-clausthal.de/paw/ (September 28, 2016).
  * Grimme  _et al._ [2010] S. Grimme, J. Antony, S. Ehrlich, and H. Krieg, “A consistent and accurate ab initio parametrization of density functional dispersion correction (DFT-D) for the 94 elements H-Pu,” J. Chem. Phys. 132, 154104 (2010).
  * Grimme, Ehrlich, and Goerigk [2011] S. Grimme, S. Ehrlich, and L. Goerigk, “Effect of the damping function in dispersion corrected density functional theory,” J. Comput. Chem. 32, 1456–1465 (2011).
  * Behler [2019] J. Behler, “RuNNer,” http://gitlab.com/TheochemGoettingen/RuNNer (August 22, 2019).
  * Kalman [1960] R. E. Kalman, “A new approach to linear filtering and prediction problems,” J. Basic Eng. 82, 35–45 (1960).
  * Blank  _et al._ [1995] T. B. Blank, S. D. Brown, A. W. Calhoun, and D. J. Doren, “Neural network models of potential energy surfaces,” J. Chem. Phys. 103, 4129–4137 (1995).
  * Glorot and Bengio [2010] X. Glorot and Y. Bengio, “Understanding the difficulty of training deep feedforward neural networks,” J. Mach. Learn. Res. 9, 249–256 (2010).
  * Plimpton [1995] S. Plimpton, “Fast parallel algorithms for short-range molecular dynamics,” J. Comput. Phys. 117, 1–19 (1995).
  * LAM [2019] “LAMMPS – Large-scale atomic/molecular massively parallel simulator,” http://lammps.sandia.gov (August 7, 2019).
  * Singraber [2019] A. Singraber, “n2p2 – A neural network potential package,” https://github.com/CompPhysVienna/n2p2 (December 9, 2019).
  * Nosé [1984] S. Nosé, “A molecular dynamics method for simulations in the canonical ensemble,” Mol. Phys. 52, 255–268 (1984).
  * Hoover [1985] W. G. Hoover, “Canonical dynamics: Equilibrium phase-space distributions,” Phys. Rev. A 31, 1695–1697 (1985).
  * Wales and Doye [1997] D. J. Wales and J. P. K. Doye, “Global optimization by basin-hopping and the lowest energy structures of Lennard-Jones clusters containing up to 110 atoms,” J. Phys. Chem. A 101, 5111–5116 (1997).
  * Kanamori [1960] J. Kanamori, “Crystal distortion in magnetic compounds,” J. Appl. Phys. 31, S14–S23 (1960).
  * Stukowski [2010] A. Stukowski, “Visualization and analysis of atomistic simulation data with OVITO – the open visualization tool,” Model. Simul. Mater. Sci. Eng. 18, 015012 (2010).
  * Sun  _et al._ [2020] Y. Sun, H. Liao, J. Wang, B. Chen, S. Sun, S. J. H. Ong, S. Xi, C. Diao, Y. Du, J.-O. Wang, M. B. H. Breese, S. Li, H. Zhang, and Z. J. Xu, “Covalency competition dominates the water oxidation structure-activity relationship on spinel oxides,” Nat. Catal. 3, 554–563 (2020).
  * Mills and Jónsson [1994] G. Mills and H. Jónsson, “Quantum and thermal effects in H2 dissociative adsorption: Evaluation of free energy barriers in multidimensional quantum systems,” Phys. Rev. Lett. 72, 1124–1127 (1994).
  * Jónsson, Mills, and Jacobsen [1998] H. Jónsson, G. Mills, and K. W. Jacobsen,  _Classical and Quantum Dynamics in Condensed Phase Simulations_ (World Scientific, Singapore, 1998).
  * Challa, Landau, and Binder [1986] M. S. S. Challa, D. P. Landau, and K. Binder, “Finite-size effects at temperature-driven first-order transitions,” Phys. Rev. B 34, 1841–1852 (1986).
  * Binder [1987] K. Binder, “Theory of first-order phase transitions,” Rep. Prog. Phys. 50, 783–859 (1987).
  * Ashcroft and Mermin [1976] N. W. Ashcroft and N. D. Mermin,  _Solid State Physics_ (Saunders College, New York, 1976).
  * Bisquert [2008] J. Bisquert, “Interpretation of electron diffusion coefficient in organic and inorganic semiconductors with broad distributions of states,” Phys. Chem. Chem. Phys. 10, 3175–3194 (2008).
  * Bisquert [2004] J. Bisquert, “Chemical diffusion coefficient of electrons in nanostructured semiconductor electrodes and dye-sensitized solar cells,” J. Phys. Chem. B 108, 2323–2332 (2004).
  * Mehrer [2007] H. Mehrer,  _Diffusion in Solids_ (Springer, Berlin, 2007).
  * Phillips [1989] J. C. Phillips,  _Physics of High- $T_{\mathrm{c}}$ Superconductors_ (Academic Press, San Diego, 1989).


# A note on ”The Klein-Gordon oscillator in (1+2)-dimensions Gürses space-time
backgrounds (Ann. Phys. (N. Y.) 404 (2019) 1)”

Omar Mustafa<EMAIL_ADDRESS>Department of Physics, Eastern
Mediterranean University, G. Magusa, north Cyprus, Mersin 10 - Turkey.
Faizuddin Ahmed<EMAIL_ADDRESS>Department of Physics, University
of Science and Technology Meghalaya, Ri-Bhoi, Meghalaya-793101, India.

###### Abstract

Abstract: We revisit and discuss the KG-oscillator in the $(1+2)$-dimensional
Gürses space-time studied in Ref. Re1  (F. Ahmed, Ann. Phys. (N. Y.) 404
(2019) 1). The modified oscillator frequency i.e.,
$\tilde{\omega}^{2}=(\Omega^{2}\,E^{2}+\eta^{2})$ appeared in the eigenvalue
equation is an energy-dependent parameter, and consequently, the results are
accurately reported here. Moreover, we show some interesting spectroscopic
features indulged within the very nature of Gürses space-time for the KG-
Gürses oscillators.

PACS numbers: 05.45.-a, 03.50.Kk, 03.65.-w

Keywords: Klein-Gordon oscillators, Gürses space-time, vorticity-energy
correlations.

In a recent paper, one of us Re1  has investigated the Klein-Gordon (KG)
particles in the $(1+2)$-dimensions Gürses space-time background Re2
described by the metric

$ds^{2}=-dt^{2}+dr^{2}-2\,\Omega\,r^{2}\,dt\,d\theta+r^{2}\,\left(1-\Omega^{2}\,r^{2}\right)\,d\theta^{2}=g_{\mu\nu}\,dx^{\mu}\,dx^{\nu};\text{
}\mu,\nu=0,1,2.$ (1)

Under such space-time metric, the covariant and contravariant metric tensors,
respectively, read

$g_{\mu\nu}=\left(\begin{tabular}[]{ccc}$-1\vskip 3.0pt plus 1.0pt minus
1.0pt$&$0$&$-\Omega r^{2}$\\\ $0$&$1\vskip 3.0pt plus 1.0pt minus
1.0pt$&$0$\\\ $-\Omega
r^{2}$&$\,0$&$\,r^{2}\left(1-\Omega^{2}r^{2}\right)$\end{tabular}\right)\Longleftrightarrow
g^{\mu\nu}=\left(\begin{tabular}[]{ccc}$\left(\Omega^{2}r^{2}-1\right)$&$0\vskip
3.0pt plus 1.0pt minus 1.0pt$&$-\Omega$\\\ $0$&$1\vskip 3.0pt plus 1.0pt minus
1.0pt$&$0$\\\ $-\Omega$&$\,0$&$\,1/r^{2}$\end{tabular}\right)\text{ };\text{ \
}g=\det\left(g_{\mu\nu}\right)=-r^{2}.$ (2)

Then the corresponding KG-equation is given by

$\frac{1}{\sqrt{-g}}\partial_{\mu}\left(\sqrt{-g}g^{\mu\nu}\partial_{\nu}\Psi\right)=m^{2}\Psi.$
(3)

At this point, the KG-oscillator is introduced through Moshinsky and
Szczepaniak recipe Re3  so that the momentum operator takes its non-minimal
coupling form

$p_{\mu}\longrightarrow
p_{\mu}+i\mathcal{F}_{\mu}\Longrightarrow\partial_{\mu}\longrightarrow\partial_{\mu}+\mathcal{F}_{\mu},$
(4)

where
$\mathcal{F}_{\mu}=\left(0,\mathcal{F}_{r},0\right)=\left(0,m\,\omega\,r,0\right)$
as proposed in Re3 ; Re4 ; Re5 ; Re6 ; Re8 . We shall, however, define
$\mathcal{F}_{\mu}=\left(0,\eta\,r,0\right)$ to avoid the confusion between
$m^{2}\equiv\left(mc^{2}\right)^{2}$ in (3) and $m$ in
$\mathcal{F}_{\mu}=\left(0,m\,\omega r,0\right)$ which identifies the non-
relativistic mass of the quantum particle (as in Schrödinger equation). To
retrieve Ahmed’s results, we substitute $\eta=m\,\omega$, therefore. Under
such settings, one rewrites (3) as

$\frac{1}{\sqrt{-g}}\left(\partial_{\mu}+\mathcal{F}_{\mu}\right)\left[\sqrt{-g}\,g^{\mu\nu}\left(\partial_{\nu}-\mathcal{F}_{\nu}\right)\Psi\right]=m^{2}\Psi,$
(5)

to obtain

$\left\\{-\partial_{t}^{2}+\left(\Omega\,r\,\partial_{t}-\frac{1}{r}\partial_{\theta}\right)^{2}+\partial_{r}^{2}+\frac{1}{r}\partial_{r}-\eta^{2}r^{2}-2\eta-m^{2}\right\\}\Psi=0.$
(6)

Then the substitution of
$\Psi=\Psi\left(t,r,\theta\right)=\exp\left(-i\,E\,t+i\,\ell\,\theta\right)\,\psi\left(r\right)$
(where $\ell$ denotes the magnetic quantum number $\ell=0,\pm 1,\pm 2,\cdots$)
would result in

$\psi^{{}^{\prime\prime}}\left(r\right)+\frac{1}{r}\psi^{{}^{\prime}}\left(r\right)+\left(\lambda-\tilde{\omega}^{2}r^{2}-\frac{\ell^{2}}{r^{2}}\right)\psi\left(r\right)=0,$
(7)

where

$\lambda=E^{2}-2\,\Omega\,\ell\,E-m^{2}-2\,\eta,\text{
}\tilde{\omega}=\sqrt{\Omega^{2}\,E^{2}+\eta^{2}}.$ (8)

We observe that equation (7) is in the form of the two-dimensional radial
Schrödinger oscillator and admits the textbook solution

$\psi\left(r\right)\sim
r^{\left|\ell\right|}\exp\left(-\frac{\left|\tilde{\omega}\right|r^{2}}{2}\right)L_{n_{r}}^{\left|\ell\right|}\left(\left|\tilde{\omega}\right|r^{2}\right)$
(9)

as the eigenfunctions, and

$\lambda=2|\tilde{\omega}|\left(2n_{r}+|\ell|+1\right)$ (10)

as the eigenvalues, where $n_{r}=0,1,2,\cdots$ is the radial quantum number.
Consequently, using (8) and (10), one would obtain the following eigenvalue
polynomial equation

$E^{2}_{n_{r},\ell}-2\,\Omega\,\ell\,E_{n_{r},\ell}-C_{n_{r},\ell}^{2}=0,$
(11)

where

$\text{ \
}C_{n_{r},\ell}^{2}=2|\tilde{\omega}|\,\left(2n_{r}+|\ell|+1\right)+m^{2}+2\,\eta$
(12)

which is exactly the same expression as equation (20), along with (22),
reported in Ref. Re1  using the parametric Nikiforov-Uvarov method.

However, the following discussion was missing in Ref. Re1  which we feel
necessary to explain. The modified oscillator frequency $\tilde{\omega}$,
given by (8) with
$|\tilde{\omega}|=|\Omega\,E|\sqrt{1+\eta^{2}/\Omega^{2}E^{2}}$, is actually
an energy-dependent parameter, and $C_{n_{r},\ell}$’s of (12) are given by

$C_{n_{r},\ell}^{2}=C_{n_{r},\ell}^{2}\left(E_{n_{r},\ell}\right)=2\,\sqrt{\Omega^{2}\,E^{2}_{n_{r},\ell}+\eta^{2}}\left(2\,n_{r}+|\ell|+1\right)+m^{2}+2\,\eta.$
(13)

Therefore, for any $n_{r}\geq 0$ the results starting from equations (25) to
(30) in Ref. Re1  seem insecure (since they are in an incomplete form). For
example, $E_{0,\ell}$ given in equation (25) in the original work is
explicitly given by

$E_{0,\ell}=\Omega\,\ell+\sqrt{\Omega^{2}\,\ell^{2}+2\,\sqrt{\Omega^{2}\,E_{0,\ell}^{2}+\eta^{2}}\left(|\ell|+1\right)+m^{2}+2\,\eta}$
(14)

for $n_{r}=0$ and so on. Actually, the eigenvalue equation (11) is a
polynomial of $E_{n_{r},\ell}$ given by

$\displaystyle
E^{2}_{n_{r},\ell}-2\,\Omega\,\ell\,E_{n_{r},\ell}-m^{2}-2\,\eta=2\,\sqrt{\Omega^{2}\,E^{2}_{n_{r},\ell}+\eta^{2}}\,\left(2\,n_{r}+|\ell|+1\right)$
$\displaystyle\Rightarrow$
$\displaystyle\left(E^{2}_{n_{r},\ell}-2\,\Omega\,\ell\,E_{n_{r},\ell}-m^{2}-2\,\eta\right)^{2}=4\,\left(\Omega^{2}\,E^{2}_{n_{r},\ell}+\eta^{2}\right)\,\left(2\,n_{r}+|\ell|+1\right)^{2}$
$\displaystyle\Rightarrow$ $\displaystyle
f_{4}\,E^{4}_{n_{r},\ell}+f_{3}\,E^{3}_{n_{r},\ell}+f_{2}\,E^{2}_{n_{r},\ell}+f_{1}\,E_{n_{r},\ell}+f_{0}=0$
(15)

where $f_{i}\quad(i=0\ldots 4)$ are non-zero parameter. The real solution of
the above polynomial (15) gives us the exact energy eigenvalue expression
$E_{n_{r},\ell}$ of the KG-oscillator field associated with the modes
$\left\\{n_{r},\ell\right\\}$ in the Gürses space-time.

Nevertheless, it would be interesting to mention that the Gürses space-time
metric (1) provides the KG-oscillators (hence KG-Gürses oscillators) as a
natural Gürses space-time byproduct. Namely, setting $\eta=0$ (i.e.,
retrieving the usual momentum operator in (4)) would result that
$\tilde{\omega}^{2}=\Omega^{2}E^{2}$ in (7) and the corresponding quadratic
energy form of equation (11) now reads

$E^{2}-2\,\Omega\,\ell\,E-2\,|\Omega\,E|\,\left(2\,n_{r}+|\ell|+1\right)-m^{2}=0.$
(16)

This energy relation is by itself interesting in the sense that $E$ denotes
the energies for both KG-particles and anti-particles, i. e., $E=E_{+}=+|E|$
for particles and $E=E_{-}=-|E|$ for anti-particles. Moreover, the vorticity
parameter $\Omega$ denotes rotation of the Gürses universe so that
$\Omega=\pm|\Omega|$ to imply that $\Omega_{+}=|\Omega|$ and
$\Omega_{-}=-|\Omega|$. Under such assumptions, one would obtain

$E_{\pm}^{2}-2\,\Omega_{\pm}E_{\pm}\,\tilde{n}_{+}-m^{2}=0;\;\tilde{n}_{+}=2n_{r}+\left|\ell\right|+\ell\,+1,$
(17)

for $\left|\Omega E\right|=\Omega_{\pm}E_{\pm}$, and

$E_{\pm}^{2}+2\,\Omega_{\mp}E_{\pm}\,\tilde{n}_{-}\,-m^{2}=0;\;\tilde{n}_{-}=2n_{r}+\left|\ell\right|-\ell\,+1.$
(18)

for $\left|\Omega E\right|=-\Omega_{\mp}E_{\pm}$. Which would allow us to cast

$E_{\pm}=\Omega_{\pm}\,\tilde{n}_{+}\pm\sqrt{\Omega^{2}\tilde{n}_{+}^{2}+m^{2}}\Rightarrow\left\\{\begin{tabular}[]{l}$E_{+}=\Omega_{\pm}\,\tilde{n}_{+}+\sqrt{\Omega^{2}\tilde{n}_{+}^{2}+m^{2}}$\\\
$E_{-}=\Omega_{\pm}\,\tilde{n}_{+}-\sqrt{\Omega^{2}\tilde{n}_{+}^{2}+m^{2}}$\end{tabular}\right.,$
(19)

for $\left|\Omega E\right|=\Omega_{\pm}E_{\pm}$ and

$E_{\pm}=-\Omega_{\mp\,}\tilde{n}_{-}\,\pm\sqrt{\Omega^{2}\tilde{n}_{-}^{2}+m^{2}}\Rightarrow\left\\{\begin{tabular}[]{l}$E_{+}=-\Omega_{-\,}\tilde{n}_{-}\,+\sqrt{\Omega^{2}\tilde{n}_{-}^{2}+m^{2}}$\\\
$E_{-}=-\Omega_{+\,}\tilde{n}_{-}\,-\sqrt{\Omega^{2}\tilde{n}_{-}^{2}+m^{2}}$\end{tabular}\right..$
(20)

However, it is more interesting to express the energies in terms of
$\Omega_{+}$ and $\Omega_{-}$so that

$E_{\pm}^{\left(\Omega_{+}\right)}=\pm\left|\Omega\right|\,\tilde{n}_{\pm}\pm\sqrt{\Omega^{2}\tilde{n}_{\pm}^{2}+m^{2}},$
(21)

for positive vorticity parameter, and

$E_{\pm}^{\left(\Omega_{-}\right)}=\pm\left|\Omega\right|\,\tilde{n}_{\mp}\pm\sqrt{\Omega^{2}\tilde{n}_{\mp}^{2}+m^{2}}.$
(22)

for negative vorticity parameter. Notably, we observe that
$\tilde{n}_{\pm}\left(\ell=\pm\left|\ell\right|\right)=\tilde{n}_{\mp}\left(\ell=\mp\left|\ell\right|\right)$
which would in effect introduce the so called vorticity-energy correlations so
that
$E_{\pm}^{\left(\Omega_{+}\right)}\left(\ell=\pm\left|\ell\right|\right)=E_{\pm}^{\left(\Omega_{-}\right)}\left(\ell=\mp\left|\ell\right|\right)$.

On the other hand, the inclusion of the coupling $\xi$ of the gravitational
field with the background curvature $R=2\Omega^{2}$ would imply that

$E_{\pm}^{\left(\Omega_{+}\right)}=\pm\left|\Omega\right|\,\tilde{n}_{\pm}\pm\sqrt{\Omega^{2}\tilde{n}_{\pm}^{2}+m^{2}+2\xi\Omega^{2}},$
(23)

for the energies with positive vorticity in (21), and

$E_{\pm}^{\left(\Omega_{-}\right)}=\pm\left|\Omega\right|\,\tilde{n}_{\mp}\pm\sqrt{\Omega^{2}\tilde{n}_{\mp}^{2}+m^{2}+2\xi\Omega^{2}}.$
(24)

for the energies with negative vorticity in (22). We again observe the
vorticity-energy correlations
$E_{\pm}^{\left(\Omega_{+}\right)}\left(\ell=\pm\left|\ell\right|\right)=E_{\pm}^{\left(\Omega_{-}\right)}\left(\ell=\mp\left|\ell\right|\right)$
for the energies in (23) and (24).

Finally, the current note is intended to show some interesting spectroscopic
features indulged within the Gürses spacetime structure, especially for the
KG-Gürses equation. It is worthwhile mentioning that in the original work Re1
, the vorticity parameter $\Omega$ is chosen to be positive and forgot to
interpreted the modified oscillator frequency $\tilde{\omega}$ for which the
exact energy is now possible to obtain from the polynomial equation (15) (but
not the equation (21) in Re1  written in an incomplete form). In this note, we
have shown that since the modified oscillator frequency $\tilde{\omega}$ given
by (8) is an energy-dependent parameter, and hence, the exact energy
expression of the KG-oscillator can be obtained by solving the polynomial
equation (15) of the real solution. Consequently, the results reported in the
original paper Re1  for different modes, such as $n_{r}=0,1,2$ are insecure
and incomplete.

Data Availability Statement: No new data generated or analyzed in this study.

Conflicts of Interest: There is no conflict of interests.

## References

  * (1) F. Ahmed, Ann. Phys. (N. Y.) 404 (2019) 1.
  * (2) M. Gürses, Class. Quantum Grav.11 (1994) 2585.
  * (3) M. Moshinsky and A. Szczepaniak, J. Phys. A: Math. Gen. 22 (1989) L817.
  * (4) S. Bruce and P. Minning, Nuovo Cimento II A 106 (1993) 711\.
  * (5) V. V. Dvoeglazov, Nuovo Cimento II A 107 (1994) 1413.
  * (6) B. Mirza and M. Mohadesi, Commun. Theor. Phys. 42 (2004) 664\.
  * (7) J. Carvalho, A. M. de M. Carvalho, E. Cavalcante and C. Furtado, Eur. Phys. J. C 76 (2016) 365.