#
Developing correlation-consistent numeric atom-centered orbital basis sets for
Krypton: Applications in RPA-based correlated calculations

Sixian Yang  Key Laboratory of Quantum Information, University of Science and
Technology of China, Hefei, 230026, China  Institute of Physics, Chinese
Academy of Sciences, Beijing 100190, China  Igor Ying Zhang
<EMAIL_ADDRESS>Collaborative Innovation Center of Chemistry for
Energy Materials, Shanghai, Key Laboratory of Molecular Catalysis and
Innovative Materials, MOE Key Laboratory of Computational Physical Sciences,
Shanghai Key Laboratory of Bioactive Small Molecules, Department of Chemistry,
Fudan University, Shanghai 200433, China  Xinguo Ren<EMAIL_ADDRESS>Institute of Physics, Chinese Academy of Sciences, Beijing 100190, China

###### Abstract

Localized atomic orbitals are the preferred basis-set choice for large-scale
explicit correlated calculations, and high-quality hierarchical correlation-
consistent basis sets are a prerequisite for correlated methods to deliver
numerically reliable results. At present, Numeric Atom-centered Orbital (NAO)
basis sets with valence correlation consistency (VCC), designated as NAO-
VCC-$n$Z, are only available for light elements from hydrogen (H) to argon
(Ar) (New J. Phys. 15, 123033, (2013) ). In this work, we extend this series
by developing NAO-VCC-$n$Z basis sets for krypton (Kr), a prototypical element
in the fourth row of periodic table. We demonstrate that NAO-VCC-$n$Z basis
sets facilitate the convergence of electronic total-energy calculations using
the Random Phase Approximation (RPA), which can be used together with a two-
point extrapolation scheme to approach the complete-basis-set (CBS) limit.
Notably, the Basis Set Superposition Error (BSSE) associated with the newly
generated NAO basis sets is minimal, making them suitable for applications
where BSSE correction is either cumbersome or impractical to do. After
confirming the reliability of NAO basis sets for Kr, we proceed to calculate
the Helmholtz free energy for Kr crystal at the theoretical level of RPA plus
renormalized single excitation (rSE) correction. From this, we derive the
pressure-volume ($P$-$V$) diagram, which shows excellent agreement with the
latest experimental data. Our work demonstrates the capability of correlation-
consistent NAO basis sets for heavy elements, paving the way toward
numerically reliable correlated calculations for bulk materials.

##  I Introduction

Kohn-Sham (KS) density functional theory (DFT) [1, 2] has become an
indispensable tool for computational studies in physics, chemistry, and
materials science. Practical applications of DFT rely on useful approximations
to the exchange-correlation (XC) functional – the only unknown component in
KS-DFT. The success of DFT lies largely in the fact that good compromise
between accuracy and the computational cost can be achieved at the level of
local and semi-local approximations [2, 3, 4, 5]. Despite their success, these
functionals suffer from intrinsic deficiencies, such as severe delocalization
errors, the lack of long-range van der Waals interaction, and the incapability
of dealing with strong electron correlations, to name a few [6]. To further
improve the accuracy of first-principles approaches, recent advancements have
focused on the development of so-called “fifth-rung” functionals. These
include the Random Phase Approximation (RPA) [7, 8, 9, 10, 11] formulated
within the adiabatic-connection fluctuation-dissipation theorem framework [12,
13], and doubly hybrid functionals [14, 15, 16, 17] that integrate Density
Functional Approximations (DFAs) with explicit non-local correlated methods,
such as the Møller–Plesset perturbation theory (MP2) [18].

One significant challenge associated with these fifth-rung functionals (and
explicit correlated methods in general) is the need for a large amount of
unoccupied single-particle orbitals. In practice, the corresponding
correlation energy exhibits slow convergence with respect to the basis set
size [19, 20, 9, 21, 22]. This issue stems from the inadequate representation
of the sharp Coulomb cusp at the electron-electron coordinate coalescence,
where the first derivative of the many-electron wave function is
discontinuous. Achieving an accurate description of such a cusp using smooth
orbital expansions necessitates a large number of basis functions [23, 24].
Over the years, considerable efforts have been devoted to addressing the slow
convergence of the electron correlation energy. One prevalent approach is the
Complete Basis Set (CBS) extrapolation scheme [25, 20, 21, 26, 27, 28], which
is widely adopted by the quantum chemistry community to mitigate the Basis Set
Incompleteness Errors (BSIE). This is accomplished by utilizing hierarchical
basis sets allowing the convergence pattern of the correlation energy to be
described through straightforward analytical formula. Another approach
involves the utilization of explicit wave-function methods, such as the
$F_{12}$ method, [24, 29, 30, 31], which explicitly account for the electronic
cusp and thereby reduce BSIE at finite basis sets.

Both techniques mentioned above have proven to be successful for light
elements and molecular systems. Gaussian-Type Orbitals (GTOs) are often the
basis set choice for molecular calculations due to their nice analytical
properties that facilitate the computation of multi-center integrals. In this
regard, the correlation-consistent GTOs developed by Dunning [32], known as
“cc-pVXZ” basis sets, have gained widespread acceptance in the field of
quantum chemistry. Their systematic convergence behavior allows for reliable
extrapolation to the CBS limit. In recent years, there has been considerable
interest in applying quantum chemistry methods and fifth-rung functionals to
extended periodic systems [33, 34]. In such cases, the standard cc-pVXZ basis
sets, while effective for molecular calculations, are not optimal. One issue
is that the extended tails of GTOs can result in an overlap matrix with a
large condition number for extended periodic systems, rendering the preceding
self-consistent calculations unstable before a sufficiently large GTO basis
set can be employed for a reliable extrapolation. Secondly, the Basis Set
Superposition Error (BSSE) associated with commonly used cc-pVXZ basis sets
can be significant for correlation methods, which is difficult to correct in
periodic systems. Within the quantum chemistry community, there is ongoing
research aimed at developing GTOs better suited for correlated calculations in
periodic systems. For example, Ye et al. [35] recently developed the “GTH-cc-
pVXZ” GTOs, which are combined with the Goedecker-Teter-Hutter (GTH) type
pseudopotentials [36, 37], to facilitate solid-state calculations using MP2
and Coupled Cluster theory at the level of Single, Double, and perturbative
Triple excitations (CCSD(T)).

Alternatively, plane-wave basis sets can be used for correlated calculations
in solid-state systems [38, 39]. However, it should be noted that the
pseudopotentials designed for ground-state DFAs may not be suitable for such
correlated calculations. For example, special pseudopotentials must be
developed for high-precision $GW$ calculations in the projector augmented wave
framework [40]. Furthermore, plane-wave basis sets make it challenging to
formulate low-scaling algorithms that fully exploit Kohn’s nearsightness
principle [41]. Although $O(N^{3})$-scaling algorithm has been developed for
RPA and $GW$ utilizing the fast Fourier transform between real and reciprocal
space [42] and the efficient Fourier transform between the imaginary time and
frequency domains [43], further scaling reductions that make use of the
intrinsic locality of the electronic structure have not yet been reported.

In the case of Numeric Atom-centered Orbitals (NAOs), the situation is even
less satisfactory. In 2013, Zhang et al. developed the valence-correlation-
consistent NAO basis sets, termed as NAO-VCC-$n$Z [44], which are apt for
correlated calculations in molecular calculations. However, NAO-VCC-$n$Z is
currently only available for light elements from H to Ar. For solid-state
calculations, these basis sets require further refinement and modification to
alleviate the ill-conditioning problem [45], which is more pronounced in
solids compared to molecules when using large basis sets. Consequently, there
is an urgent need to develop correlation-consistent NAO basis sets that are 1)
applicable to heavier elements beyond Ar, and 2) suitable for solid-state
calculations.

In this work, we first construct hierarchical valence-correlation-consistent
NAO basis sets for the element krypton (Kr), largely adhering to the
prescription developed in Ref. [44]. Specifically, the basis functions are
generated by minimizing the frozen-core PBE-based RPA total energy for a
single atom, while subject to the “even-tempered” principle, akin to the
approach used in the construction of cc-pVXZ basis sets. The candidate basis
functions are all “hydrogen-like” functions [46]. We then assess the
performance of the newly generated NAO-VCC-$n$Z ($n=$2-5) basis sets for both
Kr dimer (Kr2) and Kr crystal systems. Systematic convergence behaviors of the
binding energies are observed for both the dimer and bulk systems, allowing
for reliable extrapolation to the CBS limit. More importantly, our comparison
with and without counterpoise corrections reveals that these new NAO basis
sets carry rather small BSSE, permitting reliable calculations even without
BSSE corrections. In contrast, the cc-pVXZ GTOs and FHI-aims-2009 (tier’s)
NAOs for Kr are contaminated with huge BSSE, resulting in the binding energies
that necessitate BSSE corrections for obtaining meaningful results. The
minimal and controllable BSSE in the newly generated NAO-VCC-$n$Z for Kr
offers a significant advantage, especially in situations where BSSE
corrections are cumbersome or impractical to perform. Finally, we calculate
the Helmholtz free energy of the face-centered-cubic (FCC) Kr crystal, based
on which the $P$-$V$ curve at finite temperatures can be attained. The
electronic part of the free energy is obtained from converged RPA+rSE
calculations. The final RPA+rSE based $P$-$V$ curve shows excellent agreement
with experimental measurements for compressed Kr crystal across a broad
pressure range. The experience gained from developing these correlation-
consistent NAO basis sets for Kr is highly useful for guiding similar efforts
for other heavy elements. Works on constructing correlation-consistent NAO
basis sets for other fourth-row elements are still ongoing, and will be
published elsewhere.

The paper is organized as follows. Sec. II outlines the foundamentals of NAO
basis sets, focusing on the prescription used to generate correlation-
consistent NAOs for Kr, along with the computational details of the present
work. The major results are presented in Sec. III, including the detailed
layout of the NAO-VCC-$n$Z basis sets for Kr (Sec. III.1), a systematic study
of the convergence behavior of the new basis sets for both Kr2 and the bulk
system (Sec. III.2), and the $P$-$V$ curve based on converged RPA+rSE
calculations (Sec. III.3). The paper concludes in Sec. IV.

##  II Methodology and computational details

###  II.1 Numerical atomic-centered basis sets

In this work, all calculations are performed using the all-electron Fritz
Haber Institute ab initio molecular simulations (FHI-aims) code package [46].
The primary basis-set choice of FHI-aims is NAO $\\{\phi_{i}({\bf r})\\}$,
given by numerically tabulated radial functions multiplied by real spherical
harmonics,

$\phi_{i}({\bf r})=\frac{u_{i}(r)}{r}Y_{lm}(\Omega)\,,$ (1)

where $Y_{lm}(\Omega)$ denotes real parts ($m=0,\cdots,l$) and imaginary parts
($m=-l,\cdots,-1$) of complex spherical harmonics. In FHI-aims, one has the
flexibility to use radial functions with analytic forms, such as Slater-Type
Orbitals (STOs) or GTOs. Generally, the radial part of the NAOs $u_{i}(r)$ is
obtained by solving the following radial Schrödinger equation,

$\left[-\frac{1}{2}\frac{d^{2}}{dr^{2}}+\frac{l(l+1)}{2r^{2}}+v_{i}(r)+v_{cut}(r)\right]u_{i}(r)=\epsilon_{i}u_{i}(r)\,.$
(2)

In Eq. 2, the total effective potential contains two parts: a radial potential
$v_{i}(r)$, which governs the main shape of $u_{i}(r)$, and a confining
potential $v_{cut}(r)$ which is switched on beyond a certain distance
$r_{onset}$ from the nucleus and becomes infinite at $r_{onset}+w$. Careful
tests indicate that $w=2.0$ Å is a good default value, which ensures a smooth
decay of the radial tails of $u_{i}(r)$ to strict zero at $r_{onset}+w$. The
choice of the cutoff radius $r_{onset}$ of the confining potential
$v_{cut}(r)$ has been tested by Zhang et al. [44] and it was shown that the
energy variation is at the level of meV for $r_{onset}$ changes from 4.0 to
6.0 Å. This uncertainty is negligible compared to the typical BSIE for
correlation methods [9].

In FHI-aims, NAO basis sets feature the use of so-called “minimal basis
functions”, which are generated from self-consistent, non-spin-polarized, and
spherically symmetric free-atom radial potential $v_{i}(r)=v_{at}^{free}(r)$.
The minimal basis consists of core and valence functions for all occupied
states, exhibiting excellent performance by naturally capturing wavefunction
oscillations near the nucleus, where the nuclear $Z/r$ potential dominates.
Hence the minimal basis serves as a natural foundation for constructing NAO
basis sets for correlated methods. Also, it is worth noting that the
description of the minimal basis of NAOs requires less dense radial grids for
light elements compared to GTOs such as Dunning’s cc-pVXZ basis, as
demonstrated in Ref. [47].

###  II.2 NAO-VCC-$n$Z basis sets with valence correlation consistency

Currently, there exists a set of optimized NAO basis sets, referred to as FHI-
aims-2009 or more commonly as “tier-$n$” basis sets [46], which are the
standard choice for the ground-state calculations in FHI-aims. FHI-aims-2009
basis sets are constructed by a step-wise minimization of total energy from
LDA calculations of selected dimers. Studies employing tier-$n$ basis sets for
conventional (semi)-local and hybrid DFAs have demonstrated their high
efficacy [46, 10]. However, reliable results for correlation methods that
require virtual states can only be achieved when a counterpoise (CP)
correction is applied. In 2013, a new sequence of NAO basis sets, inspired by
Dunning’s correlation-consistency strategy, was developed [47]. These NAO
basis sets, designated as ”NAO-VCC-$n$Z”, achieve systematic convergence of
the frozen-core RPA total energy to the CBS limit. More importantly, a two-
point extrapolation to the CBS limit effectively mitigates the BSSE issue,
which is important in situations where a CP correction is impractical.
Unfortunately, these NAO-VCC-$n$Z basis sets are currently limited to light
elements from H to Ar.

In this work, we extend the work of Ref. [44, 45] to construct the NAO-
VCC-$n$Z basis sets for a prototypical fourth-row main-group element – Kr. We
group our NAO basis into “one base + two subsets”, where “one base” refers to
the minimal basis described above, and the “two subsets” include: i) the
($sp$) correlation set, ii) the polarization set on top of the FHI-aims
original minimal basis. Here, ($sp$) indicates that only $s$\- and $p$-type
orbitals are used, and the polarization set includes NAOs with higher momenta,
i.e. $d,f,g$ functions, etc. The $4s^{2}$ $4p^{6}$ valence-shell configuration
is adopted for Kr in post-KS RPA calculations, while occupied $3d$ and lower
states are frozen. We employ radial functions $u_{i}(r)$ generated via
hydrogen-like radial potential $v_{i}(r)=z_{i}/r$ for all angular momenta $l$,
and each radial function yields $2l+1$ basis functions by multiplying the real
spherical harmonics (cf. Eq. 1). Here we set the principal quantum number
$n=l+1$, resulting in nodeless radial functions similar to STOs [48, 49].

The effective nuclear charge $\\{z_{i}\\}$, which determines the actual shape
of hydrogen-like basis functions, is the major parameter to be optimized. For
those nodeless radial functions, we use “even-tempered”[50] scheme to expand
the effective nuclear charges:

$z_{i}=\alpha\beta^{i-1}\qquad i=1,2,...,N_{orb}\,.$ (3)

Here the nuclear charges $\\{z_{i}\\}$ are represented by two parameters
$(\alpha,\beta)$ which are the actual parameters to be optimized. Individual
subsets of radial functions with the same angular momentum $l$ share the same
$(\alpha,\beta)$ parameters. For example, in the case of a quadruple–$\zeta$
basis one could expand all three $s$-type functions via the optimization of a
pair of parameters $(\alpha_{1},\beta_{1})$, while a new pair of parameters
$(\alpha_{2},\beta_{2})$ should be used for three $d$-type functions in the
polarization set, and so on.

Following Ref. [47], we chose the frozen-core RPA@PBE energy of a single Kr
atom as the optimization target for the parameters characterizing the basis
functions. Our optimization procedure proceeds as follows: Initially, the $sp$
correlation subset is optimized, followed by the optimization of the
polarization subset. Subsequently, the $sp$ correlation subset undergoes a re-
optimization as the final step. The Nealdar-Mead downhill simplex algorithm
[51] is used for executing the optimization of the parameters
$(\alpha,\beta)$, with a convergence criterion of $10^{-3}$ for the effective
nuclear charge $Z$. This corresponds to an approximate convergence criterion
of $10^{-6}$ eV for the RPA total energy. To search for the global minimum, we
explored various initial guesses for $(\alpha,\beta)$ within a reasonable
window, using a tabulated grid.

###  II.3 Computational details

The total energy of RPA as a post-KS computational method is given by

$E^{\text{RPA}}=E^{\text{DFA}}-E_{xc}^{\text{DFA}}+E_{x}^{\text{EX}}+E_{c}^{\text{RPA}}\,,$
(4)

whereby the XC component $E_{xc}^{\text{DFA}}$ of the total energy of a
conventional DFA, $E^{\text{DFA}}$, is replaced by the exact-exchange energy
$E_{x}^{\text{EX}}$ plus RPA correlation energy $E_{c}^{\text{RPA}}$, both
evaluated in terms of pre-determined DFA orbitals and orbital energies.
Furthermore, a “single excitation” (SE) term, and in particular its
renormalized version (termed as renormalized SE, or rSE) $E_{c}^{\text{rSE}}$
can be added to the standard RPA total energy to alleviate the underbinding
behavior of RPA for molecules and solids. This correction scheme has been
discussed in Refs. [52, 11], and subsequently in Ref. [40]. Benchmark tests
for NAO-based RPA and RPA+rSE calculations have demonstrated the excellent
performance of these computational schemes for accurately describing both
molecular and periodic systems [10, 11, 45, 53, 54].

The employment of NAO-VCC-$n$Z basis sets allows the RPA or RPA+rSE total
energy to be extrapolated to the CBS limit via the following formula [55, 56],

$E(n)=E(\infty)-C/n^{3}\,.$ (5)

In Eq. 5, $n$ is the cardinal number of the NAO-VCC-$n$Z basis sets, and $C$
is a fitting parameter. Equation 5 is a widely used two-point extrapolation
formula initially developed for correlation-consistent Gaussian basis sets
[56], but also works well for NAO-VCC-$n$Z basis sets [44]. In this work,
basis sets with $n$=4, 5 are used for the extrapolation, implying that

$E(\infty)=\frac{E(4)4^{3}-E(5)5^{3}}{4^{3}-5^{3}}\,.$ (6)

The CBS limit results obtained via Eq. 6 are denoted as CBS(4,5). Conventional
DFAs that require only occupied states converge much faster with respect to
the basis set size. In this work, a finite NAO-VCC-5Z basis set is used for
conventional DFA calculations, which proves to be adequate. Also, for periodic
systems, a finite ${\bf k}$-point mesh is used for Brillouin-zone sampling.
For FCC primitive cell, an $8\times 8\times 8$ ${\bf k}$ grid is used , which
is sufficient to achieve sub-meV precision for the RPA total energy. Moreover,
the scaled zeroth-order regular approximation (ZORA) [57] scheme is employed
to account for the relativistic effect.

For free energy calculations, one needs to take into account of the entropy
contribution from lattice vibrations. In this work, we calculate the phonon
spectra at the PBE level, and the Helmholtz free energy at a given volume $V$
and temperature $T$ invoked to determine the $P$-$V$ diagram is given as

$\displaystyle F(V,T)=$ $\displaystyle
E^{\text{RPA+rSE}}(V)+\frac{1}{2}\sum_{{\bf q},j}\hbar\omega_{{\bf q},j}(V)+$
$\displaystyle\emph{k}_{B}T\sum_{{\bf
q},j}ln\left[2sinh\left(\frac{\hbar\omega_{{\bf
q},j}(V)}{2\emph{k}_{B}T}\right)\right]$ (7)

Here the summation over ${\bf q}$ goes over the Brillouin zone and $j$ labels
different phonon modes. Furthermore, in Eq. 7 $E^{\text{RPA+rSE}}(V)$ is the
RPA+rSE total energy at 0 K, and the second and third terms correspond to the
contributions from zero-point energies and lattice vibrations at finite
temperatures. The phonon frequencies are calculated by the finite displacement
method as implemented in the PHONOPY code [58] interfaced with FHI-aims.
Convergence criteria for the total energy and force are set to $10^{-6}$ eV
and $10^{-5}$ eV/A, respectively. A displacement of 0.002 A for the finite
displacement method is adopted. In Ref. [53], the influence of the employed
functionals on the phonon spectrum and the subsequent free energy calculations
was investigated, the conclusion was that the use of PBE functional instead of
the RPA+rSE in the phonon calculations should be adequate for calculating the
$P$-$V$ curve.

##  III Results and discussions

###  III.1 NAO-VCC-$n$Z basis sets of Kr

As discussed in Sec. II.2, each of the NAO-VCC-$n$Z basis sets for Kr consists
of a subset of $sp$ correlation functions, a subset of polarization functions,
and a base of minimal basis functions. The former two subsets are optimized in
this work, and the layouts of these two subsets for $n=2,3,4,5$ are shown in
Table 1. In addition, the same minimal basis functions of (4$s$ 3$p$ 1$d$) for
Kr are included in all NAO-VCC basis sets and added up in total count in Table
1. Note that the minimal basis functions are determined using Eq. 2 and hence
their exact shape will (slightly) depend on the actual XC potential used in
the SCF calculations.

Table 1: Number of radial functions in Kr NAO-VCC-$n$Z ($n$=2,3,4,5) basis sets. The minimal basis set of (4$s$ 3$p$ 1$d$) is not shown in this table, but is added up in the total count. | (sp) correlation subset |  | polarization subset |  | total
---|---|---|---|---|---
NAO-VCC-2Z | (1$s$ 1$p$) |  | (1$d$) |  | (5$s$ 4$p$ 2$d$)
NAO-VCC-3Z | (2$s$ 2$p$) |  | (2$d$ 1$f$) |  | (6$s$ 5$p$ 3$d$ 1$f$)
NAO-VCC-4Z | (3$s$ 3$p$) |  | (3$d$ 2$f$ 1$g$) |  | (7$s$ 6$p$ 4$d$ 2$f$ 1$g$)
NAO-VCC-5Z | (4$s$ 4$p$) |  | (4$d$ 3$f$ 2$g$ 1$h$) |  | (8$s$ 7$p$ 5$d$ 3$f$ 2$g$ 1$h$)

Similar to the Dunning’s cc-pVXZ basis sets [32] and the NAO-VCC-$n$ basis
sets for light elements [47], the newly generated NAO-VCC-$n$Z basis sets for
Kr adhere to the “correlation consistency” principle. Specifically, the
increase of the correlation energy resulting from the inclusion of the $n$-th
polarization function with angular momentum $l$ is approximately equal to the
increase stemming from the inclusion of the $(n-1)$-th polarization function
with angular momentum $l+1$. As illustrated in Fig. 1, the RPA and RPA+rSE
correlation energies for the Kr atom exhibit a clear exponential convergence
as a function of the basis set index $n$. In Table. S1 of the Supplementary
Materials (SM), the parameter settings of the NAO-VCC-$n$Z basis sets for Kr
as used in FHI-aims calculations are presented.

Figure 1:  Basis set errors of the RPA@PBE and (RPA+rSE)@PBE correlation
energies for isolated Kr atom. The CBS(4,5) results are taken as the
references.

Figure 2:  Binding energy curves of Kr2 without CP correction(left panel) and
with CP correction(right panel) are shown. NAO-VCC-$n$Z basis set is
extrapolated to CBS limit CBS(4,5) and the other basis sets are using the
largest size of its family (e.g. tier2+, cc-pV5Z).

###  III.2 Binding energy curves of the Kr dimer and solid bulk

We now examine the performance of NAO-VCC-$n$Z basis sets on calculating the
binding energies of the Kr dimer – Kr2. Figure 2 presents the (RPA+rSE)@PBE
binding energy curves for Kr2, which are calculated using NAO-VCC-$n$Z with
$n=2$ to 5, as well as the extrapolated CBS(4,5) limit. The left and right
panels of Fig. 2 show the results without and with CP corrections,
respectively. Additionally, we performed RPA+rSE calculations using Gaussian
cc-pV5Z and FHI-aims-2009 NAO basis sets for comparative analysis. For the
latter, we choose tier 2 plus some additional functions (denoted here as
“tier2+”). It is the largest basis set of this type available for Kr. To
isolate the impact of basis-set errors, a very tight convergence criterion is
applied to numerical integration in both radial and angular parts. The
calculated equilibrium bond lengths and the corresponding binding energies for
various basis sets are listed in Table 2.

Table 2: Equilibrium bond lengths $r_{e}$ (in Å) and interaction energies $D_{e}$ (in meV) of RPA+rSE obtained using various different basis sets. Results obtained both with (right side) and without (left side) CP corrections are presented. Basis set | no CP correction |  | CP corrected
---|---|---|---
|  $r_{e}$ |  $D_{e}$ |  |  $r_{e}$ |  $D_{e}$
NAO-VCC-3Z | 4.31 | -6.71 |  | 4.46 | -4.96
NAO-VCC-4Z | 4.15 | -9.90 |  | 4.29 | -6.63
NAO-VCC-5Z | 4.10 | -12.37 |  | 4.21 | -8.76
NAO-VCC CBS(4,5) | 4.05 | -15.12 |  | 4.15 | -11.22
tier2+ | 3.58 | -90.92 |  | 4.08 | -15.39
cc-pV5Z | 3.69 | -121.53 |  | 4.15 | -10.79

Figure 2 shows that the newly generated NAO-VCC-$n$Z basis sets for Kr yield
systematically converging RPA+rSE binding energies for the Kr2. More
importantly, a comparison between the left and right panels of Fig. 2 reveals
that the BSSEs associated with NAO-VCC-$n$Z are rather small, corresponding to
a small fraction of the total binding energy. In stark contrast, binding
energies calculated using cc-pV5Z and the NAO tier2+ basis sets are
contaminated with substantial BSSE. Without performing CP corrections, both
cc-pV$5$Z and tier2+ significantly overestimate the binding energies by an
order of magnitude, rendering the results entirely unusable. These
discrepancies are evident in both Fig. 2 and Table 2. Further results for
other cc-pV$X$Z basis sets are presented in Fig. S1 of the SM.

It should be noted that, by introducing an augmented set of diffuse functions,
aug-cc-pV$X$Z basis sets usually outperform cc-pV$X$Z, in particular for weak
interactions [47]. However, unlike aug-cc-pV$X$Z basis sets, our NAO-VCC-$n$Z
basis sets do not incorporate an augmented set of diffuse functions. As shown
in Fig. S2 in the SM, the introduction of diffuse functions to NAO-VCC-$5$Z
leads to a CP-corrected bonding energy of -13.70 meV, which is similar to the
NAO-VCC-5Z result without CP corrections (-12.37 meV, see Table 2). This is
appealing since in many situations it is too cumbersome or impractical to
perform BSSE corrections, and in such situations the newly generated NAO-
VCC-$n$Z basis set for Kr can be directly used, expecting that the small and
controllable BSSE of this basis set compensates the effect of missing diffuse
functions.

Figure 3:  Fitting curve under second-order Birch–Murnaghan equation of state
method for the FCC Kr crystal. The experimental cohesive energy corrected by
zero-point vibration energy from accurate CCSD(T) is shown as black dashed
line, -121.8 meV.

Figure 4:  Cohesive energies of the FCC Kr crystal determined by different
functionals. The RPA and RPA+rSE results are obtained by extrapolating the
NAO-VCC-4Z and 5Z results to the complete basis set (CBS) limit, and all other
results are obtained using the NAO-VCC-5Z basis set. A FCC primitive unit cell
and $8\times 8\times 8$ k grid sampling are used in all calculations.
Experimental values are corrected with zero-point vibration energies by
accurate quantum chemical method from reference [59, 60, 61].

.

Next, we examine the performance of NAO-VCC-$n$Z basis sets in calculating
cohesive energies using RPA-type methods for bulk systems. Figure 3 presents
the cohesive energy curves for Kr FCC crystal based on RPA (left side) and
RPA+rSE (right side) methods using NAO-VCC-$n$Z basis sets with $n$=2-5. As
before, the figure also includes the extrapolated CBS(4,5) results. A
consistent convergence pattern relative to the cardinal number $n$ is again
observed. However, it’s worth noting that for bulk Kr, the results obtained
using 2Z or 3Z basis sets are still quite far from the CBS limit. Therefore,
extrapolation based on larger basis sets, such as 4Z and 5Z, is necessary.
Furthermore, a comparison between the left and right sides of Fig. 3 reveals
that the rSE correction significantly enhances the binding strength.
Concretely, RPA tends to underbind, while RPA+rSE tends to overbind the Kr
crystal, respectively.

In Fig. 4, we present the cohesive energies of the FCC Kr crystal as a
function of the lattice constant calculated using RPA@PBE and (RPA+rSE)@PBE.
For comparison, we also include results of lower-rung functionals, including
LDA, GGA–PBE, the global hybrid functional PBE0 [62], the
Heyd–Scuseria–Ernzerhof (HSE) screened hybrid functional [63], as well as PBE
complemented with vdW(TS) [64] and with the many-body dispersion (MBD) [65].
The experimental equilibrium lattice constant and cohesive energy are marked
by dashed lines in the figure. Here the reference cohesive energy has been
corrected for zero-point vibration energies as determined by accurate quantum
chemical methods [61]. Similar to other noble-gas crystals, an accurate
description of the cohesive properties of the Kr crystal is challenging for
conventional DFAs. LDA and PBE exhibit substantial overbinding and
underbinding tendencies, respectively. Hybrid functionals like PBE0 or HSE
offer no improvement over PBE, due to their inability to capture long-range
vdW interactions. Conventional GGA-PBE combined with semi-empirical vdW
corrections, like PBE+vdW(TS) and PBE+MBD, improves accuracy for weak
interactions but still underestimates the cohesive energy and overestimates
the equilibrium lattice constant.

Regarding the performance of RPA, we obtained an equilibrium lattice constant
of $5.770$ Åand a cohesive energy of -101.7 meV at the CBS(4,5) level.
Compared to the experimental data of $5.670\AA$ and -121.8 meV [60], the
standard RPA@PBE scheme exhibits a noticeable underbinding, despite its
significant improvement over conventional DFAs. This underbinding tendency of
RPA has also been observed by Klimeš et al. [40], based on plane-wave-based
RPA implementation [66]. As previously discussed, RPA@PBE typically underbinds
molecules and solids due to overestimated Pauli repulsion, which can be
largely corrected by adding single-excitation corrections [52, 11, 40]. Figure
4 shows that incorporating the rSE correction enhances the cohesive energies
of the Kr crystal, resulting in an equilibrium lattice constant of $5.633$
Åand a cohesive energy of -136.3 meV. While this is generally satisfactory, an
overestimation of cohesive energy by 14 meV is noticeable, which can be partly
attributed to the lack of BSSE correction in the bulk calculations. Table. S2
of the SM indicates that, with the CP correction, the binding energies are
reduced by about 24 meV, leading to a BSSE-corrected RPA+rSE binding energy of
-112.1 meV. It is worth pointing out that this value is close to the NAO-
VCC-$5$Z result (see Fig 4), indicating a compensation between the small BSSE
and the BSIE in the finite basis set. The fully numerically converged RPA+rSE
cohesive energy for Kr FCC crystal is most likely in between the CP-corrected
and uncorrected values, aligning well with the experimental result.

###  III.3 $P$-$V$ curve of Kr FCC crystal in the high-pressure regime

The study of matter under high pressure is an important area of research in
condensed-matter physics, as a plethora of intriguing phenomena might emerge
in this regime. Fascinating structural transition or abrupt changes in
conductivity and other electronic-structure-derived properties are often
observed in experiments. For example, it has been shown that the application
of high pressure can induce a phase transition from trilayer graphene to
hexagonal diamond structure [67], and a structural transformation from 2Hc to
2Ha in MoS2 [68]. At ultrahigh pressure, superconductivity has even been
observed in pristine 2Ha-MoS2 [69]. Similar phenomena have been reported in
the family of rare-gas solids. Specifically, a FCC-HCP phase transition has
been experimentally observed [70, 71, 72, 73, 74], accompanied by changes in
the band structure [75, 76]. More recently, it was found that Kr may exhibit
significant stacking disorders in its FCC lattice. Furthermore, the FCC and
HCP phases have been found to coexist over a wide pressure range from 2.7 to
140 GPa upon high-temperature annealing via laser heating [77]. While
considerable theoretical efforts have been devoted to the understanding of the
physical properties of rare gas systems under high pressure [78, 79, 80, 81],
a quantitative description achieving experimental accuracy remains a
significant challenge.

Figure 5:  Calculated $P$–$V$ diagram at 300 K for the FCC phase of the Kr
crystal. The results are extrapolated to the CBS limit [CBS(4,5)]. The
electronic part of the free energy is obtained using RPA+rSE, complemented
with PBE phonon free energy at 300 K (orange line). Reference data from four
experimental studies are included for comparison.

Having verified the reliability of our NAO-VCC-$n$Z basis set for Kr and the
good performance of the RPA+rSE method in describing the cohesive property of
the Kr crystal, we set out to determine the $P$-$V$ diagram theoretically. To
achieve this, we calculate the Helmholtz free energy $F(V,T)$ according to Eq
7. Subsequently, we determine the pressure at a given temperature and volume
via $P(V,T)=-\partial F(V,T)/\partial V|_{T}$, which allows us to generate the
$P$-$V$ curves at various temperatures. In Fig. 5, we present the calculated
(RPA+rSE)@PBE and PBE $P$-$V$ curves at $T=300$ K for the FCC phase of Kr
crystal. Furthermore, we include four sets of experimental data for
comparative analysis.

Figure 5 shows that our calculated RPA+rSE $P$–$V$ curve is in overall good
agreement with the experimental results of Polian et al. [82] and Jephcoat et
al. [83] extracted from energy-dispersive x-ray diffraction data, and those of
Errandonea et al. [73] and Rosa et al. [77] extracted from in-situ x-ray
diffraction and absorption measurements. Note that, there are noticeable
discrepancies among the various sets of experimental data themselves. Here, we
focus on the comparison to the experimental data from Rosa et al., as it
represents the most recent and comprehensive investigation of Kr under high-
pressure conditions, extending up to 140 GPa. Overall our theoretical $P$-$V$
curve agrees well with the experimental data of Rosa et al. in both high and
low pressure regimes. For comparison, we also include the calculated PBE
$P$-$V$ curve in Fig. 5. It is evident that the PBE method tends to
overestimate the pressure (i.e., a higher pressure is needed to compress the
Kr crystal to the same volume), a tendency that is pronounced in the low and
medium pressure range.

##  IV Conclusion

In this work, we developed valence-correlation-consistent NAO basis sets for
the Kr element. Our tests on the binding energies of Kr2 and the cohesive
energies of Kr FCC crystal demonstrate that these basis sets yield systematic
convergence patterns for RPA-type correlated calculations. This allows for a
reliable extrapolation to the CBS limit. More importantly, these newly
generated NAO-VCC-$n$Z basis sets for Kr exhibit minimal BSSE, in stark
contrast with the Gaussian cc-pVXZ basis sets and the FHI-aims-2009 NAO basis
sets. This feature makes them particularly useful in the correlated
calculations of extended materials, where a CP correction for BSSE is often
impractical. We then calculated the Helmholtz free energy and derived the
$P$-$V$ curve for the FCC phase of the Kr crystal, based on converged RPA+rSE
calculations for the electronic part of the free energy. The resulting
theoretical $P$-$V$ curve is in excellent agreement with the most recent
experiment data. The insights and methodologies developed in this study offer
valuable guidance for constructing NAO-VCC basis sets for other heavy
elements. Overall, our research investigates the capability of correlation-
consistent NAO basis sets for heavy elements, laying the groundwork for
conducting reliable correlated calculations at the RPA level and beyond in
extended materials.

###### Acknowledgements.

The work is supported by the National Key R & D Program of China (Grant No.
2022YFA1403800) and the National Natural Science Foundation of China (Grant
Nos. 12188101, 12134012, 21973015, and 22125301), as well as by the Max Planck
Partner Group for Advanced Electronic Structure Methods, and innovative
research team of high-level local universities in Shanghai and a key
laboratory program of the Education Commission of Shanghai Municipality
(ZDSYS14005).

## Additional Information

Competing interest: The authors declare no competing financial or non-
financial interests.

## Data Availability

All the input and output files of the calculations involved in this work have
been uploaded to the NOMAD repository and can be found under the link
https://dx.doi.org/10.17172/NOMAD/2023.09.01-1.

## References

  * Hohenberg and Kohn [1964] P. Hohenberg and W. Kohn, Phys. Rev. 136, B864 (1964).
  * Kohn and Sham [1965] W. Kohn and L. J. Sham, Phys. Rev. 140, A1133 (1965).
  * Perdew and Zunger [1981] J. P. Perdew and A. Zunger, Self-Interaction Correction to Density-Functional Approximation for Many-Electron Systems, Phys. Rev. B 23, 5048 (1981).
  * Perdew  _et al._ [1996a] J. P. Perdew, K. Burke, and M. Ernzerhof, Generalized Gradient Approximation Made Simple, Phys. Rev. Lett 77, 3865 (1996a).
  * Sun  _et al._ [2015] J. Sun, A. Ruzsinszky, and J. P. Perdew, Strongly constrained and appropriately normed semilocal density functional, Phys. Rev. Lett. 115, 036402 (2015).
  * Cohen  _et al._ [2012] A. J. Cohen, P. Mori-Sánchez, and W. Yang, Challenges for Density Functional Theory, Chemical Reviews 112, 289 (2012), pMID: 22191548, https://doi.org/10.1021/cr200107z .
  * Bohm and Pines [1953] D. Bohm and D. Pines, A Collective Description of Electron Interactins: III. Coulomb Interaction in a Degenerate Electron Gas, Phys. Rev. 92, 609 (1953).
  * Eshuis  _et al._ [2012] H. Eshuis, J. E. Bates, and F. Furche, Electron Correlation Methods Based on the Random Phase Approximation, Theor. Chem. Acc. 131, 1084 (2012).
  * Ren  _et al._ [2012a] X. Ren, P. Rinke, C. Joas, and M. Scheffler, Random-phase approximation and its applications in computational chemistry and materials science, J. Mater. Sci. 47, 7447 (2012a).
  * Ren  _et al._ [2012b] X. Ren, P. Rinke, V. Blum, J. Wieferink, A. Tkatchenko, A. Sanfilippo, K. Reuter, and M. Scheffler, New J. Phys. 14, 053020 (2012b).
  * Ren  _et al._ [2013] X. Ren, P. Rinke, G. E. Scuseria, and M. Scheffler, Renormalized second-order perturbation theory for the electron correlation energy: Concept, implementation, and benchmarks, Phys. Rev. B 88, 035120 (2013).
  * Langreth and Perdew [1977] D. C. Langreth and J. P. Perdew, Exchange-correlation energy of a metal surface: Wave-vector analysis, Phys. Rev. B 15, 2884 (1977).
  * Gunnarsson and Lundqvist [1976] O. Gunnarsson and B. I. Lundqvist, Exchange and correlation in atoms, molecules, and solids by the spin-density-functional formalism, Phys. Rev. B 13, 4274 (1976).
  * Zhao  _et al._ [2004] Y. Zhao, B. J. Lynch, and D. G. Truhlar, Development and Assessment of a New Hybrid Density Functional Model for Thermochemical Kinetics, J. Phys. Chem. A 108, 2715 (2004).
  * Grimme [2006] S. Grimme, Semiempirical hybrid density functional with perturbative second-order correlation, J. Chem. Phys. 124, 034108 (2006).
  * Zhang  _et al._ [2009] Y. Zhang, X. Xu, and W. A. Goddard III, Doubly hybrid density functional for accurate descriptions of nonbond interactions, thermochemistry, and thermochemical kinetics, Proc. Nat. Acad. Sci. U. S. A. 106, 4963 (2009).
  * Zhang and Xu [2014] I. Y. Zhang and X. Xu,  _A New-Generation Density Functional_ (Springer, New York, 2014).
  * Møller and Plesset [1934] C. Møller and M. S. Plesset, Phys. Rev. 46, 618 (1934).
  * Lebègue  _et al._ [2010] S. Lebègue, J. Harl, T. Gould, J. G. Ángyán, G. Kresse, and J. F. Dobson, Cohesive Properties and Asymptotics of the Dispersion Interaction in Graphite by the Random Phase Approximation, Phys. Rev. Lett. 105, 196401 (2010).
  * Harl  _et al._ [2010] J. Harl, L. Schimka, and G. Kresse, Assessing the quality of the random phase approximation for lattice constants and atomization energies of solids, Phys. Rev. B 81, 115126 (2010).
  * Eshuis and Furche [2012] H. Eshuis and F. Furche, Basis set convergence of molecular correlation energy differences with the random phase approximation, J. Chem. Phys. 136, 084105 (2012).
  * Janesko  _et al._ [2009] B. G. Janesko, T. M. Henderson, and G. E. Scuseria, Long-range-corrected hybrids including random phase approximation correlation, J. Chem. Phys. 130, 081105 (2009).
  * Tew  _et al._ [2007] D. P. Tew, W. Klopper, and T. Helgaker, Electron correlation: The many-body problem at the heart of chemistry, Journal of computational chemistry 28, 1307 (2007).
  * Noga  _et al._ [1992] J. Noga, W. Kutzelnigg, and W. Klopper, CC-R12, a correlation cusp corrected coupled-cluster method with a pilot application to the Be2 potential curve, Chemical physics letters 199, 497 (1992).
  * Grüneis  _et al._ [2009] A. Grüneis, M. Marsman, J. Harl, L. Schimka, and G. Kresse, Making the random phase approximation to electronic correlation accurate, J. Chem. Phys. 131, 154115 (2009).
  * Furche [2001] F. Furche, Molecular tests of the random phase approximation to the exchange-correlation energy functional, Phys. Rev. B 64, 195120 (2001).
  * Fabiano and Dalla Sala [2012] E. Fabiano and F. Dalla Sala, Accuracy of basis-set extrapolation schemes for DFT-RPA correlation energies in molecular calculations, Theor. Chem. Acc. 131, 1278 (2012).
  * Grüneis  _et al._ [2010] A. Grüneis, M. Marsman, and G. Kresse, Second-order Møller–Plesset perturbation theory applied to extended systems. II. Structural and energetic properties, J. Chem. Phys. 133, 074107 (2010).
  * Hattig  _et al._ [2012] C. Hattig, W. Klopper, A. Kohn, and D. P. Tew, Explicitly correlated electrons in molecules, Chemical reviews 112, 4 (2012).
  * Kong  _et al._ [2012] L. Kong, F. A. Bischoff, and E. F. Valeev, Explicitly correlated R12/F12 methods for electronic structure, Chemical reviews 112, 75 (2012).
  * Kutzelnigg and Klopper [1991] W. Kutzelnigg and W. Klopper, Wave functions with terms linear in the interelectronic coordinates to take care of the correlation cusp. I. General theory, J. Chem. Phys. 94, 1985 (1991).
  * T. H. Dunning [1989] J. T. H. Dunning, Gaussian basis sets for use in correlated molecular calculations. I. The atoms boron through neon and hydrogen, J. Chem. Phys. 90, 1007 (1989).
  * Wang  _et al._ [2021] Y. Wang, Y. Li, J. Chen, I. Y. Zhang, and X. Xu, Doubly Hybrid Functionals Close to Chemical Accuracy for Both Finite and Extended Systems: Implementation and Test of XYG3 and XYGJ-OS, JACS Au 1, 543 (2021), publisher: American Chemical Society.
  * Zhang and Grüneis [2019] I. Y. Zhang and A. Grüneis, Coupled Cluster Theory in Materials Science, Front. Mater. 6, Article 123 (2019).
  * Ye and Berkelbach [2022] H.-Z. Ye and T. C. Berkelbach, Correlation-consistent Gaussian basis sets for solids made simple, Journal of Chemical Theory and Computation 18, 1595 (2022).
  * Goedecker  _et al._ [1996] S. Goedecker, M. Teter, and J. Hutter, Separable dual-space Gaussian pseudopotentials, Phys. Rev. B 54, 1703 (1996).
  * Hartwigsen  _et al._ [1998] C. Hartwigsen, S. Goedecker, and J. Hutter, Relativistic separable dual-space Gaussian pseudopotentials from H to Rn, Phys. Rev. B 58, 3641 (1998).
  * Marsman  _et al._ [2009] M. Marsman, A. Grüneis, J. Paier, and G. Kresse, Second-order Møller–Plesset perturbation theory applied to extended systems. I. Within the projector-augmented-wave formalism using a plane wave basis set, The Journal of Chemical Physics 130, 184103 (2009), https://pubs.aip.org/aip/jcp/article-pdf/doi/10.1063/1.3126249/9452255/184103_1_online.pdf .
  * Shepherd  _et al._ [2012] J. J. Shepherd, A. Grüneis, G. H. Booth, G. Kresse, and A. Alavi, Convergence of many-body wave-function expansions using a plane-wave basis: From homogeneous electron gas to solid state systems, Phys. Rev. B 86, 035111 (2012).
  * Klimeš  _et al._ [2015] J. Klimeš, M. Kaltak, E. Maggio, and G. Kresse, Singles correlation energy contributions in solids, The Journal of Chemical Physics 143, 102816 (2015).
  * Kohn [1996] W. Kohn, Phys. Rev. Lett. 76, 3168 (1996).
  * Rojas  _et al._ [1995] H. N. Rojas, R. W. Godby, and R. J. Needs, Space-Time Method for ab initio Calculations of Self-Energies and Dielectric Response Functions of Solids, Phys. Rev. Lett. 74, 1827 (1995).
  * Klimeš  _et al._ [2014] J. K. Klimeš, M. Kaltak, and G. Kresse, Predictive GW calculations using plane waves and pseudopotentials, Phys. Rev. B 90, 075125 (2014).
  * Zhang  _et al._ [2013a] I. Y. Zhang, X. Ren, P. Rinke, V. Blum, and M. Scheffler, Numeric atom-centered-orbital basis sets with valence-correlation consistency from H to Ar, New J. Phys. 15, 123033 (2013a).
  * Zhang  _et al._ [2019] I. Y. Zhang, A. J. Logsdail, X. Ren, S. V. Levchenko, L. Ghiringhelli, and M. Scheffler, Test set for materials science and engineering with user-friendly graphic tools for error analysis: Systematic benchmark of the numerical and intrinsic errors in state-of-the-art electronic-structure approximations, New J. Phys. 21, 013025 (2019).
  * Blum  _et al._ [2009] V. Blum, F. Hanke, R. Gehrke, P. Havu, V. Havu, X. Ren, K. Reuter, and M. Scheffler, Ab-initio molecular simulations with numeric atom-centered orbitals, Comp. Phys. Comm. 180, 2175 (2009).
  * Zhang  _et al._ [2013b] I. Y. Zhang, X. Ren, P. Rinke, V. Blum, and M. Scheffler, Numeric atom-centered-orbital basis sets with valence-correlation consistency from H to Ar, New Journal of Physics 15, 123033 (2013b).
  * Van Lenthe and Baerends [2003] E. Van Lenthe and E. J. Baerends, Optimized slater-type basis sets for the elements 1–118, Journal of Computational Chemistry 24, 1142 (2003).
  * Te Velde  _et al._ [2001] G. t. Te Velde, F. M. Bickelhaupt, E. J. Baerends, C. Fonseca Guerra, S. J. van Gisbergen, J. G. Snijders, and T. Ziegler, Chemistry with ADF, Journal of Computational Chemistry 22, 931 (2001).
  * Raffenetti [1973] R. C. Raffenetti, Even-tempered atomic orbitals. II. Atomic SCF wavefunctions in terms of even-tempered exponential bases, J. Chem. Phys. 59, 5936 (1973).
  * Nelder and Mead [1965] J. A. Nelder and R. Mead, A simplex method for function minimization, The computer journal 7, 308 (1965).
  * Ren  _et al._ [2011] X. Ren, A. Tkatchenko, P. Rinke, and M. Scheffler, Beyond the Random Phase Approximation for the Electron Correlation Energy: The Importance of Single Excitations, Phys. Rev. Lett. 106, 153003 (2011).
  * Yang and Ren [2022] S. Yang and X. Ren, Phase stability of the argon crystal: first-principles study based on random phase approximation plus renormalized single excitation corrections, New J. Phys. 24, 033049 (2022).
  * Tahir  _et al._ [2022] M. N. Tahir, T. Zhu, H. Shang, J. Li, V. Blum, and X. Ren, Localized resolution of identity approach to the analytical gradients of random-phase approximation ground-state energy: algorithm and benchmarks, Journal of Chemical Theory and Computation 18, 5297 (2022).
  * Halkier  _et al._ [1998] A. Halkier, T. Helgaker, P. Jørgensen, W. Klopper, H. Koch, J. Olsen, and A. K. Wilson, Chem. Phys. Lett. 286, 243 (1998).
  * Helgaker  _et al._ [1997] T. Helgaker, W. K. H. Koch, and J. Noga, Basis-set convergence of correlated calculations on water, J. Chem. Phys. 106, 9639 (1997).
  * van Lenthe  _et al._ [1994] E. van Lenthe, E. Baerends, and J. Snijders, Relativistic total energy using regular approximations, J. Chem. Phys. 101, 9783 (1994).
  * Togo and Tanaka [2015] A. Togo and I. Tanaka, First principles phonon calculations in materials science, Scr. Mater. 108, 1 (2015).
  * Schwalbe  _et al._ [1977] L. A. Schwalbe, R. K. Crawford, H. H. Chen, and R. A. Aziz, Thermodynamic consistency of vapor pressure and calorimetric data for argon, krypton, and xenon, J. Chem. Phys. 66, 4493 (1977).
  * Losee and Simmons [1968] D. Losee and R. Simmons, Thermal-expansion measurements and thermodynamics of solid krypton, Phys. Rev. 172, 944 (1968).
  * Rościszewski  _et al._ [2000] K. Rościszewski, B. Paulus, P. Fulde, and H. Stoll, Ab initio coupled-cluster calculations for the fcc and hcp structures of rare-gas solids, Phys. Rev. B 62, 5482 (2000).
  * Perdew  _et al._ [1996b] J. P. Perdew, M. Ernzerhof, and K. Burke, Rationale of mixing exact exchange with density functional approximations, J. Chem. Phys. 105, 9982 (1996b).
  * Heyd  _et al._ [2003] J. Heyd, G. E. Scuseria, and M. Ernzerhof, Hybrid functionals based on a screened coulomb potential, J. Chem. Phys. 118, 8207 (2003).
  * Tkatchenko and Scheffler [2009] A. Tkatchenko and M. Scheffler, Phys. Rev. Lett. 102, 073005 (2009).
  * Tkatchenko  _et al._ [2012] A. Tkatchenko, R. A. DiStasio Jr., R. Car, and M. Scheffler, Phys. Rev. Lett. 108, 236402 (2012).
  * Harl and Kresse [2008] J. Harl and G. Kresse, Cohesive energy curves for noble gas solids calculated by adiabatic connection fluctuation-dissipartion theory, Phys. Rev. B 77, 045136 (2008).
  * Ke  _et al._ [2020] F. Ke, L. Zhang, Y. Chen, K. Yin, C. Wang, Y.-K. Tzeng, Y. Lin, H. Dong, Z. Liu, J. S. Tse, _et al._ , Synthesis of atomically thin hexagonal diamond with compression, Nano letters 20, 5916 (2020).
  * Nayak  _et al._ [2014] A. P. Nayak, S. Bhattacharyya, J. Zhu, J. Liu, X. Wu, T. Pandey, C. Jin, A. K. Singh, D. Akinwande, and J.-F. Lin, Pressure-induced semiconducting to metallic transition in multilayered molybdenum disulphide, Nature communications 5, 3731 (2014).
  * Chi  _et al._ [2018] Z. Chi, X. Chen, F. Yen, F. Peng, Y. Zhou, J. Zhu, Y. Zhang, X. Liu, C. Lin, S. Chu, _et al._ , Superconductivity in pristine 2Ha\- MoS2 at ultrahigh pressure, Phys. Rev. Lett. 120, 037002 (2018).
  * Jephcoat  _et al._ [1987] A. Jephcoat, H. Mao, L. Finger, D. Cox, and R. Hemley, Pressure-induced structural phase transitions in solid xenon, Phys. Rev. Lett. 59, 2670 (1987).
  * Cynn  _et al._ [2001] H. Cynn, C. S. Yoo, B. Baer, V. Iota-Herbei, A. K. McMahan, M. Nicol, and S. Carlson, Martensitic fcc-to-hcp Transformation Observed in Xenon at High Pressure, Phys. Rev. Lett. 86, 4552 (2001).
  * Goettel  _et al._ [1989] K. A. Goettel, J. H. Eggert, I. F. Silvera, and W. C. Moss, Optical evidence for the metallization of xenon at 132 (5) gpa, Phys. Rev. Lett. 62, 665 (1989).
  * Errandonea  _et al._ [2002] D. Errandonea, B. Schwager, R. Boehler, and M. Ross, Phase behavior of krypton and xenon to 50 GPa, Phys. Rev. B 65, 214110 (2002).
  * Errandonea  _et al._ [2006] D. Errandonea, R. Boehler, S. Japel, M. Mezouar, and L. R. Benedetti, Structural transformation of compressed solid Ar: An X-ray diffraction study to 114 GPa, Phys. Rev. B 73, 092106 (2006).
  * Reichlin  _et al._ [1989] R. Reichlin, K. E. Brister, A. K. McMahan, M. Ross, S. Martin, Y. K. Vohra, and A. L. Ruoff, Evidence for the Insulator-Metal Transition in Xenon from Optical, X-Ray, and Band-Structure Studies to 170 GPa, Phys. Rev. Lett. 62, 669 (1989).
  * Ross  _et al._ [2005] M. Ross, R. Boehler, and P. Söderlind, Xenon melting curve to 80 GPa and $5p\mathrm{\text{$-$}}d$ hybridization, Phys. Rev. Lett. 95, 257801 (2005).
  * Rosa  _et al._ [2018] A. Rosa, G. Garbarino, R. Briggs, V. Svitlyk, G. Morard, M. A. Bouhifd, J. Jacobs, T. Irifune, O. Mathon, and S. Pascarelli, Effect of the fcc-hcp martensitic transition on the equation of state of solid krypton up to 140 GPa, Phys. Rev. B 97, 094115 (2018).
  * Schwerdtfeger  _et al._ [2017] P. Schwerdtfeger, K. G. Steenbergen, and E. Pahl, Relativistic coupled-cluster and density-functional studies of argon at high pressure, Phys. Rev. B 95, 214116 (2017).
  * Kim  _et al._ [2006] E. Kim, M. Nicol, H. Cynn, and C.-S. Yoo, Martensitic fcc-to-hcp Transformations in Solid Xenon under Pressure: A First-Principles Study, Phys. Rev. Lett. 96, 035504 (2006).
  * Kwon  _et al._ [1995] I. Kwon, L. A. Collins, J. D. Kress, and N. Troullier, First-principles study of solid Ar and Kr under high compression, Phys. Rev. B 52, 15165 (1995).
  * Dewhurst  _et al._ [2002] J. K. Dewhurst, R. Ahuja, S. Li, and B. Johansson, Lattice dynamics of solid xenon under pressure, Phys. Rev. Lett. 88, 075504 (2002).
  * Polian  _et al._ [1989] A. Polian, J. Itie, E. Dartyge, A. Fontaine, and G. Tourillon, X-ray absorption spectroscopy on solid krypton up to 20 GPa, Phys. Rev. B 39, 3369 (1989).
  * Jephcoat [1998] A. P. Jephcoat, Rare-gas solids in the earth’s deep interior, Nature 393, 355 (1998).


# Surface Recognition for e-Scooter
Using Smartphone IMU Sensor

Areej Eweida, Nimord Segol, Maxim Freydin, , Niv Sfaradi, and Barak Or
Submitted: Jan 2023.All authors are with ALMA Technologies Ltd, Haifa, 340000,
Israel (e-mail: barak@almatechnologies.com).

###### Abstract

In recent years, as the use of micromobility gained popularity, technological
challenges connected to e-scooters became increasingly important. This paper
focuses on road surface recognition, an important task in this area. A
reliable and accurate method for road surface recognition can help improve the
safety and stability of the vehicle. Here a data-driven method is proposed to
recognize if an e-scooter is on a road or a sidewalk. The proposed method uses
only the widely available inertial measurement unit (IMU) sensors on a
smartphone device. deep neural networks (DNNs) are used to infer whether an
e-scooter is driving on a road or on a sidewalk by solving a binary
classification problem. A data set is collected and several different deep
models as well as classical machine learning approaches for the binary
classification problem are applied and compared. Experiment results on a route
containing the two surfaces are presented demonstrating the DNNs ability to
distinguish between them.

###### Index Terms:

Micromobility, Surface recognition, Deep Neural Network, Inertial Measurement
Unit, Machine Learning

##  I Introduction

Road surface recognition is the process of identifying and classifying
different types of road surfaces that a vehicle may encounter during operation
[1, 2, 3, 4]. This work focuses on road surface recognition for e-scooters. A
mode of transportation that has become increasingly popular in urban areas in
recent years. In particular, this work presents a method for distinguishing
between paved asphalt roads and sidewalks for e-scooters using only a mobile
device.

Accurate road surface recognition is important for several reasons. First,
different road surfaces can have different impacts on the stability and
handling of an e-scooter [5]. For example, a smooth, paved road may provide a
stable and comfortable ride, while a textured, rough sidewalk may be more
difficult to navigate. By identifying the road surface, e-scooters can adjust
their performance accordingly, ensuring a safe and enjoyable ride for the
user. Additionally, road surface recognition can be used to improve the
efficiency of e-scooters. Different road surfaces can have different levels of
resistance, which can impact the amount of energy required to travel a given
distance. By accurately identifying the road surface, e-scooters can optimize
their power output and energy consumption, extending the range of the vehicle
and improving the overall user experience. Finally, driving an e-scooter on a
sidewalk is illegal in many places. An accurate method that determines if a
rider is on a sidewalk or on a road can provide fleet managers and
micromobility companies with a valuable tool.

In recent years, several approaches were developed to address the problem of
road surface recognition for different transportation modes. In[6] the authors
present a road surface recognition system for bicycles using several sensors
attached to the bicycle including a camera an inertial measurement unit (IMU)
and a speedometer. In [7], the authors described a method to estimate the
roughness of a surface on which a car is traversing using LiDAR as well as a
global positioning system (GPS) receiver and IMU. In [8] the authors describe
a method of surface recognition for pedestrians using an IMU and a microphone
attached to the users’ shoes. A method to estimate the condition of the road
surface using lasers mounted on a vehicle was proposed in [9]. In [10] a deep
classifier that estimates the road surface was trained using data collected by
cameras. In [1] an algorithm is presented that extracts features from an IMU
sensor on a car and uses K nearest neighbors algorithm to detect abnormal road
surface conditions.

Other related works where machine learning algorithms, and in particular deep
neural networks, were applied to IMU sensor data include [11] where DNN-based
multi-models together with an EKF were proposed to deal with the global
navigation satellite system (GNSS) outage. In [12] it was demonstrated that
integration of DNN into classical signal processing methods can boost
navigation performance. In [13] the authors demonstrate how a DNN model can
estimate a mounting yaw angle of a smartphone strapped to a vehicle. Another
related work is [14] where DNNs were used to learn linear velocities from IMU
sensors in cars. To the best of the author’s knowledge, this work is the first
to present a road surface recognition algorithm, including a DNN, for
e-scooters.

In this paper, a novel approach to road surface recognition for e-scooters
using IMU sensors on smartphones is developed. The method is data-driven and
relies on a DNN binary classifier. IMU sensors are widely available in
smartphones. By mounting a smartphone on an e-scooter and collecting data from
the IMU sensors, and applying a DNN-based algorithm, it is possible to infer
the type of road surface the e-scooter is traveling on. The DNN is able to
detect the signature associated with each type of road surface and discern
between roads and sidewalks. The proposed approach has several advantages: It
is cost-effective, as it utilizes sensors that are already present in many
smartphones; It is also lightweight and compact, making it well-suited for use
on e-scooters; Additionally, smartphones are increasingly being used as
onboard computers for e-scooters, providing a convenient platform for
integrating the suggested road surface recognition system.

The rest of the paper is organized as follows: In Section II the data
collection and dataset generation processes are presented. In Section III the
proposed ML algorithms are detailed. In Section IV the experimental results
are described. In particular, the results of several algorithms on a
validation set are presented. An experiment on a route containing both types
of surfaces is also presented. Finally, in Section V conclusions are given.

##  II Data set creation

Any data-driven method requires an offline collection of data and model
training. This section contains a description of the data collection and
processing. A schematic view of the data process can be seen in Fig. (2).

###  II-A Raw data collection

IMU acceleration and gyroscope data was collected with an iPhone 13 strapped
to an e-scooter’s handle. Data was recorded at the maximum rate of 100 [Hz]. A
total of 100 minutes were collected for this task in 12 distinct driving
sessions conducted by the same rider. Different driving sessions had distinct
routes to increase diversity in the data set. See Fig. (1) for an image of the
route and the collecting device (smartphone). Of the routes, 7 were entirely
on a sidewalk, and the rest were on the road. To create a validation set, 2 of
the drives were set aside. All in all, 10 drives totaling 80 minutes were used
for training and the two drives, one on a sidewalk and one on the road, for
validation.

Figure 1: Left: an example of a street with both sidewalk and an asphalt road.
Right: a smartphone device mounted on a scooter. Note the device frame of
reference.

###  II-B Pre-processing

The collected data consists of acceleration $[m/s^{2}]$ and angular velocity
$[rad/s]$ in three axes. Before passing the acceleration and gyroscope data to
the classifier, the following pre-processing tasks are applied:

  1. 1.

Apply a low-pass filter with a cut-off frequency of 10 [Hz].

  2. 2.

Down-sample the data to 20 [Hz].

  3. 3.

Rotate (only roll and pitch) the measured acceleration and angular velocity
vectors to the navigation frame where the first component points to the
e-scooter forward direction and the third component is parallel to the gravity
vector.

  4. 4.

Subtract the gravity component from the third acceleration component.

These steps are applied both offline when training the model and online as a
part of the classification algorithm. The goal of the described processing is
to simplify the signal and avoid unnecessary noise and other information that
can affect the performance of the classifier. It is noted that the above
processing procedure has appeared in various prior works e.g. [15, 16, 13]. It
is provided in Figure 2 for this work to be self-contained.

Figure 2: Schematic description of the dataset generation process.

###  II-C Data generation

To create the dataset the recorded drives were partitioned to windows of $S$
seconds. In this work the value of $S$ was considered a hyper-parameter in the
set $\\{1,2,3\\}$. The overall data set consists of pairs
$X\in\mathrm{R}^{20\cdot S\times 6}$ of IMU measurements of $S[s]$ sampled in
20 [Hz] and $Y\in\\{0,1\\}$ with $1$ indicating a sidewalk and $0$ indicating
a paved road.

##  III Classifiers

Using the dataset collected as described above we train a total of 5 ML
models. Each of the models outputs one or zero signifying sidewalk or road.

###  III-A Models

Five ML models where tested. The most common DNN-based models for this type of
task are convolutional neural network (CNN) [17] and long short-term memory
(LSTM) [18]. Each one of them were used as a basis for a DNN architecture and
both of them were combined in a fused architecture. In addition, two classical
ML models were tested: support vector machine (SVM) [19] and random forest
[20]. A description of each model is provided:

  1. 1.

CNN \- A convolutional neural network with two convolution layers with 8
channels and 3 kernels followed by dropout and mean pooling with a rectified
linear unit (ReLU). One convolution layer followed by dropout and max pooling
with ReLU. Two linear layers with ReLU with 20 and 1 channels, respectively.

  2. 2.

LSTM \- A simplified LSTM layer with hidden size of 100 followed by a dense
layer.

  3. 3.

LSTM-CNN \- An LSTM layer with hidden size 100 followed by two convolution
layers with 8 channels and 3 kernels each, followed by dropout and mean
pooling with ReLU. The model ends with a single linear layer, as presented in
Fig.3.

  4. 4.

SVM \- Each sample $X\in\mathbb{R}^{20\cdot S\times 6}$ was flattened to a
vector $\tilde{X}\in\mathbb{R}^{120\cdot S}$ and a classifier was trained on
it.

  5. 5.

Random Forest \- Each sample $X\in\mathbb{R}^{20\cdot S\times 6}$ was
flattened to a vector $\tilde{X}\in\mathbb{R}^{120\cdot S}$ and an ensemble
random forest classifier was trained. The maximum depth allowed per tree was
set to 100.

###  III-B DNNs loss function and training hyperparameters

The loss function and training criteria for the DNN-based models are described
in this section.

####  III-B1 Loss function for DNN-based models

All of the DNNs considered in this work were trained using the binary cross
entropy loss, the most common loss function used for binary classification
problems. It is calculated using the following formula:

$\ell(p_{i},Y_{i})=-\left[{Y_{i}\cdot\log{(p_{i})}+(1-Y_{i})\cdot\log{(1-p_{i})}}\right],$
(1)

where the output of the DNN is $p_{i}=p(X_{i},{\bf{\Theta})}\in[0,1]$, which
depends on the sample $X_{i}$ and the weights $\bf{\Theta}$ and
$Y_{i}\in\\{0,1\\}$ is the true label of example $i$. The entire loss is

$L=\frac{1}{N}\sum_{i=1}^{N}\ell(p_{i},Y_{i}),$ (2)

where $N$ is the number of samples in the train set.

####  III-B2 Training DNN-based models

The DNN-based models were all trained for 30 epochs using the ADAM optimizer
[21] with a learning rate of $0.001$, $\beta$ coefficients of $\beta_{1}=0.9$,
$\beta_{2}=0.999$, and batch size $6$. We trained the model using a Tesla K80
GPU. We trained for approximately 3.5 minutes on average in each combination
of architecture and window size, 8 seconds per epoch.

##  IV Results

This section contains the results of each of the considered classifiers on the
validation set. The models were evaluated for windows of $S=\\{1,2,3\\}$
seconds. The model which obtained the best results on the validation set was
then tested on a separate drive containing both road and sidewalk.

###  IV-A Evaluation Metrics

Three evaluation metrics for the performance of the models were considered.

  1. 1.

Accuracy: accuracy is the proportion of true results among the total number of
cases examined.

  2. 2.

F1 Score: the F1 score is a number between 0 and 1 and represents the harmonic
mean of precision and recall.

  3. 3.

False Positive Rate (FPR): FPR is the proportion of samples classified as 1
when the true label was 0.

Each of these metrics was computed for the validation set for all the ML
models.

Model |  Window Size [s] |  Accuracy |  F1 |  FPR
---|---|---|---|---
CNN |  1 |  0.78 |  0.72 |  0.22
2 |  0.83 |  0.83 |  0.15
3 |  0.78 |  0.78 |  0.16
LSTM |  1 |  0.86 |  0.88 |  0.15
2 |  0.85 |  0.86 |  0.17
3 |  0.52 |  0.47 |  0.23
LSTM-CNN |  1 |  0.88 |  0.88 |  0.14
2 |  0.90 |  0.90 |  0.13
3 |  0.93 |  0.93 |  0.13
SVM |  1 |  0.73 |  0.70 |  0.16
2 |  0.72 |  0.68 |  0.16
3 |  0.71 |  0.66 |  0.18
Random Forest |  1 |  0.79 |  0.79 |  0.23
2 |  0.78 |  0.79 |  0.22
3 |  0.80 |  0.82 |  0.20
TABLE I: Performance metrics for each model on the validation set for
different window sizes.

###  IV-B Best performing model

The optimal solution for the classification problem on the validation set is
the LSTM-CNN model. The benefit of this model is its both spatially (IMU
signals information) and temporally (time) properties. The architecture layers
are presented in Fig. 3. The best results on the validation set were obtained
for a window of $3[s]$. Table I presents performance metrics for each of the
ML models considered in this work. Note that the model outperforms all other
models where it achieves maximal accuracy of $0.93$, maximal F1 of $0.93$, and
minimal FPR of $0.13$.

Figure 3: DNN architecture: LSTM, two convolutional layers are followed by one
dense layer with a connecting average-pooling layer. Finally, a linear layer
eventually outputs the probability of getting each one of the classes.

###  IV-C The performance of classical models

The two classical models considered were the SVM and Random Forest. In
contrast to the LSTM and CNN-based models, the classical approaches do not
treat the input signal as time series but as a flattened vector of features.
Table I shows a surprising result that both classical models performed well
and closely approached the metrics of the DNNs. This may suggest that the data
is separable as is even without significant effort in handcrafting features.
However, it is important to note that the IMU data was filtered, downsampled,
rotated, and the gravity component was subtracted.

###  IV-D Field experiment

An experiment was conducted on a single drive containing both road and
sidewalk segments. The route is shown in Fig. (4) and it was not used in the
training or validation sets. The light green segment is a sidewalk, the light
gray is an asphalt road, and the red points show false classification. The
route begins on a sidewalk and at a certain point, switches to an asphalt
road.

The trained LSTM-CNN model with a $3[s]$ window was chosen for this task as it
performed best on the validation set. A recording of raw IMU measurement was
passed through the pre-processing described in Section II-B and fed to the
classifier. The model was activated every $5[s]$. The classifier has
successfully classified $90\%$ of the route. In addition, a strong wind was
present during the recording of this experiment, which was not the case for
the training data. The wind has led to elevated levels of noise and vibration
which were not present in the training or validation datasets. However, the
results of this experiment agree well with the performance observed on the
validation set.

Figure 4: Results of the selected classifier on a route containing both
sidewalk and road.

##  V Conclusions

A data-driven approach to the road surface recognition problem for e-scooters
was suggested. For e-scooters it is specifically important to distinguish
between roads and sidewalks. A binary classifier was trained to perform this
task. Several models both based on deep learning and classical approaches were
tested. Experimental results showed that a DNN-based model produces good
results. Future research directions include developing a multi-class
classifier that distinguishes between several types of surfaces and developing
a real-time system that incorporates data from several sensors using a mobile
device to deduce the road surface.

## References

  * [1] R. Du, G. Qiu, K. Gao, L. Hu, and L. Liu, “Abnormal road surface recognition based on smartphone acceleration sensor,” _Sensors_ , vol. 20, no. 2, p. 451, 2020.
  * [2] A. Bystrov, M. Abbas, E. Hoare, T.-Y. Tran, N. Clarke, M. Gashinova, and M. Cherniakov, “Analysis of classification algorithms applied to road surface recognition,” in _2015 IEEE Radar Conference (RadarCon)_. IEEE, 2015, pp. 0907–0911.
  * [3] M. Aki, T. Rojanaarpa, K. Nakano, Y. Suda, N. Takasuka, T. Isogai, and T. Kawai, “Road surface recognition using laser radar for automatic platooning,” _IEEE Transactions on Intelligent Transportation Systems_ , vol. 17, no. 10, pp. 2800–2810, 2016.
  * [4] S. Cafiso, A. Di Graziano, V. Marchetta, and G. Pappalardo, “Urban road pavements monitoring and assessment using bike and e-scooter as probe vehicles,” _Case Studies in Construction Materials_ , vol. 16, p. e00889, 2022.
  * [5] J. Leoni, M. Tanelli, S. C. Strada, and S. M. Savaresi, “Assessing e-scooters safety and drivability: a quantitative analysis,” _IFAC-PapersOnLine_ , vol. 55, no. 24, pp. 260–265, 2022.
  * [6] M. Springer and C. Ament, “A mobile and modular low-cost sensor system for road surface recognition using a bicycle,” in _2020 IEEE International Conference on Multisensor Fusion and Integration for Intelligent Systems (MFI)_ , 2020, pp. 360–366.
  * [7] K. Chen, “Road roughness recognition based on lidar,” _Journal of Physics: Conference Series_ , vol. 2278, no. 1, p. 012008, may 2022.
  * [8] H. Mitake, H. Watanabe, and M. Sugimoto, “Footsteps and inertial data-based road surface condition recognition method,” in _Proceedings of the 18th International Conference on Mobile and Ubiquitous Multimedia_ , ser. MUM ’19. New York, NY, USA: Association for Computing Machinery, 2019.
  * [9] K. Urano, K. Hiroi, S. Kato, N. Komagata, and N. Kawaguchi, “Road surface condition inspection using a laser scanner mounted on an autonomous driving car,” in _2019 IEEE International Conference on Pervasive Computing and Communications Workshops (PerCom Workshops)_ , 2019, pp. 826–831.
  * [10] E. Šabanovič, V. Žuraulis, O. Prentkovskis, and V. Skrickij, “Identification of road-surface type using deep neural networks for friction coefficient estimation,” _Sensors_ , vol. 20, no. 3, 2020.
  * [11] J. Liu and G. Guo, “Vehicle localization during gps outages with extended kalman filter and deep learning,” _IEEE Transactions on Instrumentation and Measurement_ , vol. 70, pp. 1–10, 2021.
  * [12] B. Or and I. Klein, “A Hybrid Model and Learning-Based Adaptive Navigation Filter,” _IEEE Transactions on Instrumentation and Measurement_ , pp. 1–1, 2022.
  * [13] M. Freydin, N. Sfaradi, N. Segol, A. Eweida, and B. Or, “Mountnet: Learning an inertial sensor mounting angle with deep neural networks,” _arXiv preprint arXiv:2212.11120_ , 2022.
  * [14] H. Yan, Q. Shan, and Y. Furukawa, “RIDI: Robust IMU double integration,” in _Proceedings of the European Conference on Computer Vision (ECCV)_ , 2018, pp. 621–636.
  * [15] J. Farrell, _Aided navigation: GPS with high rate sensors_. McGraw-Hill, Inc., 2008.
  * [16] M. Freydin and B. Or, “Learning car speed using inertial sensors for dead reckoning navigation,” _IEEE Sensors Letters_ , vol. 6, no. 9, pp. 1–4, 2022\.
  * [17] K. O’Shea and R. Nash, “An introduction to convolutional neural networks,” _arXiv preprint arXiv:1511.08458_ , 2015.
  * [18] S. Hochreiter and J. Schmidhuber, “Long short-term memory,” _Neural computation_ , vol. 9, no. 8, pp. 1735–1780, 1997.
  * [19] C. Cortes and V. Vapnik, “Support-vector networks,” _Machine learning_ , vol. 20, no. 3, pp. 273–297, 1995.
  * [20] L. Breiman, “Random forests,” _Machine learning_ , vol. 45, no. 1, pp. 5–32, 2001.
  * [21] D. P. Kingma and J. Ba, “Adam: A method for stochastic optimization,” _arXiv preprint arXiv:1412.6980_ , 2014.