# Approaching the Basis Set Limit in Gaussian-Orbital-Based Periodic
Calculations with Transferability: Performance of Pure Density Functionals for
Simple Semiconductors

Joonho Lee<EMAIL_ADDRESS>Department of Chemistry, Columbia University,
New York, NY, USA  Xintian Feng  Q-Chem Inc., Pleasanton, CA, USA  Leonardo A.
Cunha  Department of Chemistry, University of California, Berkeley, CA, USA
Jérôme F. Gonthier  Department of Chemistry, University of California,
Berkeley, CA, USA  Evgeny Epifanovsky  Q-Chem Inc., Pleasanton, CA, USA
Martin Head-Gordon  Department of Chemistry, University of California,
Berkeley, CA, USA

###### Abstract

Simulating solids with quantum chemistry methods and Gaussian-type orbitals
(GTOs) has been gaining popularity. Nonetheless, there are few systematic
studies that assess the basis set incompleteness error (BSIE) in these GTO-
based simulations over a variety of solids. In this work, we report a GTO-
based implementation for solids, and apply it to address the basis set
convergence issue. We employ a simple strategy to generate large uncontracted
(unc) GTO basis sets, that we call the unc-def2-GTH sets. These basis sets
exhibit systematic improvement towards the basis set limit as well as good
transferability based on application to a total of 43 simple semiconductors.
Most notably, we found the BSIE of unc-def2-QZVP-GTH to be smaller than 0.7
m$E_{h}$ per atom in total energies and 20 meV in band gaps for all systems
considered here. Using unc-def2-QZVP-GTH, we report band gap benchmarks of a
combinatorially designed meta generalized gradient functional (mGGA), B97M-rV,
and show that B97M-rV performs similarly (a root-mean-square-deviation (RMSD)
of 1.18 eV) to other modern mGGA functionals, M06-L (1.26 eV), MN15-L (1.29
eV), and SCAN (1.20 eV). This represents a clear improvement over older pure
functionals such as LDA (1.71 eV) and PBE (1.49 eV) though all these mGGAs are
still far from being quantitatively accurate. We also provide several
cautionary notes on the use of our uncontracted bases and on future research
on GTO basis set development for solids.

##  I Introduction

Condensed phase simulations using quantum chemistry tools originally developed
for molecules have gained popularity over many years,Ayala  _et al._ (2001);
Katagiri (2005); Izmaylov and Scuseria (2008); Pisani  _et al._ (2008, 2012);
Del Ben  _et al._ (2012); Usvyat  _et al._ (2015); McClain  _et al._ (2017);
Shang and Yang (2020) with the hope of enabling development of new
systematically improvable tools that can go beyond standard density functional
approaches,Maurer  _et al._ (2019) as well as existing Green’s function
methodsHybertsen and Louie (1986); Aryasetiawan and Gunnarsson (1998) in the
field. These simulations can be broadly categorized into two classes: (1)
large $\Gamma$-point calculations to describe spatial inhomegeneity as found
in gas, liquid, and surface simulations and (2) calculations with a relatively
small unit cell and a large number of $\mathbf{k}$-points as relevant for
simulations of solids. The former category resembles large cluster
calculations that are routinely performed in the molecular community and the
use of Gaussian-type orbitals (GTOs) as a computational basis is not uncommon
and numerically well-behaved. The use of GTOs to reach the thermodynamic limit
(TDL) of (dense) solids often faces numerical difficulties associated with
overcompleteness of GTOs that leads to a severe linear dependency among basis
functions towards the TDL.Löwdin (1970); Klahn and Bingel (1977); VandeVondele
and Hutter (2007); Peintinger  _et al._ (2013) Nonetheless, many studies have
employed Gaussian basis sets either using those developed for molecular
calculations, those developed for periodic mean-field
calculations,VandeVondele and Hutter (2007); Peintinger  _et al._ (2013); Laun
_et al._ (2018); Oliveira  _et al._ (2019) or those optimized system-
specifically without much in the way of transferability guarantees.Morales and
Malone (2020); Daga  _et al._ (2020); Zhou  _et al._ (2021) The use of GTOs
for solid-state calculations has been growing as well exemplified by many
existing GTO-based quantum chemistry programs with the periodic boundary
condition capability.Li  _et al._ (1995); Soler  _et al._ (2002); Artacho  _et
al._ (2008); Verzijl and Thijssen (2012); Kudin and Scuseria (2000);
Balasubramani  _et al._ (2020); Dovesi  _et al._ (2020); Kühne  _et al._
(2020); Prentice  _et al._ (2020)

The development of compact GTO basis setsJensen (2013); Nagy and Jensen (2017)
has a long history in molecular quantum chemistry.Shavitt (1993); Pritchard
_et al._ (2019) Since McWeeny’s first proposalMcweeny (1950) and Boys’ early
attemptBoys S. (1950) to use GTOs for molecular systems, many developments on
contracted Gaussian basis sets such as atomic natural orbitalAlmlöf and Taylor
(1987), correlation-consistentDunning (1989) and polarization-consistentJensen
(2001, 2002) basis sets have made high-accuracy quantum chemistry calculations
practical. However, these highly optimized contracted basis sets are usually
not considered applicable to solids due to emerging linear
dependencies.Peintinger  _et al._ (2013) In the early days of basis set
development, even-temperedBardo and Ruedenberg (1973); Feller and Ruedenberg
(1979) and well-temperedHuzinaga  _et al._ (1985) bases were explored as a
means to obtain high-quality results using only primitive GTOs reducing the
complications in sophisticated optimization procedures for exponents and
contraction coefficients. In the even-tempered bases, one employs three
parameters for each angular momentum shell $l$ to define a set of “even-
tempered” primitive GTOs by

$\phi_{lmk}(\mathbf{r})\propto\exp(-\zeta_{lk}r^{2})r^{l}S_{lm}(\Omega)$ (1)

where $\phi_{lmk}$ is an atomic orbital, $l$ and $m$ are angular momentum
quantum number, $S_{lm}(\Omega)$ are the real spherical harmonics at a solid
angle $\Omega$, $k$ sets the total number of primitive GTOs for $l,m$, and
$\zeta_{lk}$ is parameterized by a geometric series,

$\zeta_{lk}=\alpha_{l}\beta_{l}^{k-1},\>\>\alpha_{l},\beta_{l}>0,\>\>\beta_{l}\neq
1$ (2)

In the well-tempered variants, a more sophisticated form is used for
$\zeta_{lk}$. In the even-tempered basis, one needs to pick a total of three
parameters $k$, $\alpha_{lk}$, and $\beta_{lk}$. The appropriate values may be
obtained by looking at atoms and small molecules though finding these values
can generally be tedious.Feller and Ruedenberg (1979) Even-tempered basis sets
are generally much larger than contracted GTOs and thus they are rarely used
in modern quantum chemistry calculations. Nonetheless, these bases have not
yet been explored in the context of solid-state applications.

In this work, we propose an even simpler basis set generation protocol than
that of even-tempered bases which does not involve any optimizations. Our
procedure is to generate large uncontracted GTO bases that yield density
functional theory (DFT) total energies per cell within 0.7 m$E_{h}$ per atom
in the unit cell from the complete basis set limit obtained by planewave (PW)
basis. The idea is to take two existing GTO bases (one from the
def2-seriesWeigend and Ahlrichs (2005) and SZV-MOLOPT-SR-GTHVandeVondele and
Hutter (2007)), uncontract these bases, and take the union of the resulting
primitive GTOs while removing core orbitals that are treated by the underlying
GTH pseudopotential. Like the even-tempered bases, our sets are much larger
than typical contracted GTOs available in the literature, but they are not
optimized for specific systems and/or mean-field methods so they should
naturally bear transferability.

As an application of these bases, we focus on the basic goal of quantifying
the basis set error of Gaussian-based DFT calculations. This goal is even more
important to reach when considering correlated wavefunction calculations.
However, the basis set incompleteness error (BSIE) in correlation energies can
be quantified and characterized only after the underlying mean-field energy is
converged to the basis set limit. The BSIE was directly quantified by
employing the same pseudopotential proposed by Hutter and co-workers (called
the GTH pseudopotential)Goedecker  _et al._ (1996); Hartwigsen  _et al._
(1998) in both the new Gaussian-based program developed in this work and a PW-
based code, Quantum Espresso (QE).Giannozzi  _et al._ (2020)

Furthermore, we also apply our basis set to validating the performance of ten
selected pure exchange-correlation (XC) functionals. These ten XC functionals
consist of one local density approximation (LDA) functional,Slater (1951);
Perdew and Zunger (1981) five generalized gradient approximation (GGA)
functionals (PBE,Perdew  _et al._ (1996) PBEsol,Perdew  _et al._ (2008)
revPBE,Zhang and Yang (1998) BLYP,Becke (1988); Lee  _et al._ (1988)
B97-DGrimme (2006)), and four meta GGA (mGGA) functionals (SCAN,Sun  _et al._
(2015) M06-L,Zhao and Truhlar (2006) MN15-L,Yu  _et al._ (2016)
B97M-rVMardirossian and Head-Gordon (2015); Mardirossian  _et al._ (2017)).
Our benchmark set has a total of 43 semiconductors where 40 of them were taken
from the SC40 setHeyd  _et al._ (2005) and the remaining 3 (LiH,Vidal and
Vidal-Valat (1986); Nolan  _et al._ (2009); Hoang and Van de Walle (2013)
LiF,Matsushita  _et al._ (2011); Garza and Scuseria (2016) and LiClMatsushita
_et al._ (2011); Garza and Scuseria (2016)) were taken from other places. The
performance of LDA and PBE on the majority of these systems using GTOs was
already documented in ref. 61 though the underlying BSIE of the associated GTO
basis sets is unclear. Many PW-based codes including QE have LDA, GGA, and
SCAN functionals available so it is not very difficult to assess their
performance using PW-based codes.Giannozzi  _et al._ (2020) In fact, the
performance of LDA and GGA functionals, as well as the SCAN mGGA, is
relatively well understood for band gap problems.Perdew  _et al._ (2017); Yang
_et al._ (2016) However, the recently developed functionals that were
combinatorially optimized for main group molecular chemistry,
$\omega$B97X-V,Mardirossian and Head-Gordon (2014) $\omega$B97M-V,Mardirossian
and Head-Gordon (2016) and B97M-V,Mardirossian and Head-Gordon (2015);
Mardirossian  _et al._ (2017) have rarely appeared in condensed phase
studiesPestana  _et al._ (2017); Ruiz Pestana  _et al._ (2018); Pestana  _et
al._ (2019); Lininger  _et al._ (2021); Li  _et al._ (2021) and are relatively
less common and less used in PW-based codes. The same is true for the
Minnesota functionals (M06-L and MN15-L). Replacing the -V tail with the -rV
tail (the rVV10 van der Waals (vdW) correctionSabatini  _et al._ (2013)
instead of the VV10 vdW correctionVydrov and Van Voorhis (2010)), an efficient
implementation of the -rV tail is now available in some planewave-based
codes.Giannozzi  _et al._ (2020) Aside from the computational cost associated
with the long-range exact exchange, an efficient implementation of these
functionals should be readily possible. These combinatorially optimized
functionals were found to be statistically the best XC functionals at each
rung of Jacob’s ladder for main group chemistry problems,Mardirossian and
Head-Gordon (2017), and they have performed very well in other molecular
benchmarks also.Goerigk  _et al._ (2017); Najibi and Goerigk (2018) In the
condensed phase, the mGGA, B97M-rV appears to describe properties of liquid
water as accurately as far more computationally demanding hybrid
functions.Ruiz Pestana  _et al._ (2018) However, the performance of B97M-rV
for band gap problems is largely unknown at present. Motivated by this, we
report the performance of B97M-rV for band gaps here.

This paper is organized as follows: (1) we first review basic formalisms of
periodic mean-field calculations, the gaussian planewave (GPW) density fitting
scheme, and an efficient implementation of rVV10, (2) we then describe our
strategies for generating transferable GTO bases for simulating solids towards
the TDL, (3) we discuss computational details, (4) we present results for
basis set convergence of DFT total energies and band gaps using the proposed
bases, (5) we assess the performance of pure XC functionals comparing against
experimental band gaps, (6) we deliver cautionary notes on using our bases and
on the future basis set development for solids featuring striking failures of
existing GTH bases, and (7) we then conclude.

##  II Theory

Periodic mean-field calculations using a linear combination of atomic orbitals
have been well-documented in many places.Evarestov (2007, 2012) Nonetheless,
we aim to give a pedagogical review of the relevant theories on periodic DFT
calculations within the GPW implementation and the implementation of rVV10
since these are the key compute kernels in our new implementation. Experienced
readers may skip some of the subsequent sections and start from Section II.4.

###  II.1 Periodic Mean-Field Calculations

As a consequence of real-space translational symmetry, crystal momentum
($\mathbf{k}$) is a good quantum number. Periodic mean-field (PMF)
calculations with GTOs are hence done with crystalline molecular orbitals
(CMOs), $\\{\psi_{i}^{\mathbf{k}}\\}$,Kudin and Scuseria (1998)

$\psi_{i}^{\mathbf{k}}(\mathbf{r})=\sum_{\mu}C_{\mu
i}^{\mathbf{k}}\phi_{\mu}^{\mathbf{k}}(\mathbf{r})$ (3)

where crystalline atomic orbitals (CAOs) are defined with a lattice summation,

$\phi_{\mu}^{\mathbf{k}}(\mathbf{r})=\sum_{\mathbf{R}}\phi_{\mu}^{\mathbf{R}}(\mathbf{r})e^{i\mathbf{k}\cdot\mathbf{R}}$
(4)

where $\mathbf{R}$ denotes a lattice vector represented by a sum of integer
multiples of primitive vectors of the direct lattice. In PMF calculations,
analogously to their molecular counterparts, the PMF energy is minimized when
the CMO coefficient matrix obeys a self-consistent Roothaan-Hall equation,

$\mathbf{F}^{\mathbf{k}}\mathbf{C}^{\mathbf{k}}=\mathbf{S}^{\mathbf{k}}\mathbf{C}^{\mathbf{k}}\mathbf{\epsilon}^{\mathbf{k}}$
(5)

where $\mathbf{F}^{\mathbf{k}}$ is the Fock matrix at $\mathbf{k}$,
$\mathbf{S}^{\mathbf{k}}$ is the overlap matrix of CAOs at $\mathbf{k}$
defined as

$S_{\mu\nu}^{\mathbf{k}}=\sum_{\mathbf{R}}\langle{\phi}_{\mu}^{\mathbf{0}}|{\phi}_{\nu}^{\mathbf{R}}\rangle
e^{i\mathbf{k}\cdot\mathbf{R}}=\sum_{\mathbf{R}}S_{\mu\nu}^{\mathbf{0}\mathbf{R}}e^{i\mathbf{k}\cdot\mathbf{R}},$
(6)

and $\mathbf{\epsilon}^{\mathbf{k}}$ is the band energy at $\mathbf{k}$.

In periodic calculations with GTOs, it is very common to observe linear
dependencies of the CAOs which makes the metric (overlap) matrix
$\mathbf{S}^{\mathbf{k}}$ poorly conditioned. Within finite precision computer
arithmetic, the resulting truncation error in the inverse metric can lead to
numerical instability, convergence issues, and non-trivial errors in the PMF
energies. Therefore, handling exact and near linear dependencies is crucial in
GTO-based periodic calculations especially when one attempts to get to the
basis set limit where linear dependencies become progressively severe. In this
work, we adopted the canonical orthogonalization procedure.Szabo and Ostlund
(1996) The canonical orthogonalization procedure is defined as follows:

  1. 1.

The diagonalization of $\mathbf{S}^{\mathbf{k}}$ is performed for each
$\mathbf{k}$:

$\mathbf{S}^{\mathbf{k}}=\mathbf{U}^{\mathbf{k}}\mathbf{s}^{\mathbf{k}}(\mathbf{U}^{\mathbf{k}})^{\dagger}$
(7)

  2. 2.

For a given threshold $\epsilon_{\text{lindep}}$, one retains the
$N_{\text{CMO}}^{\mathbf{k}}$ eigenvalues in $\mathbf{s}^{\mathbf{k}}$ above
$\epsilon_{\text{lindep}}$ along with their corresponding eigenvectors. We
refer these subsets of eigenvalues and eigenvectors to as
$\tilde{\mathbf{s}}^{\mathbf{k}}$ and $\tilde{\mathbf{U}}^{\mathbf{k}}$,
respectively.

  3. 3.

We then define the orthogonalization matrix $\mathbf{X}^{\mathbf{k}}$,

$\mathbf{X}^{\mathbf{k}}=\tilde{\mathbf{U}}^{\mathbf{k}}(\tilde{\mathbf{s}}^{\mathbf{k}})^{-1/2}$
(8)

The dimension of $\mathbf{X}^{\mathbf{k}}$ is
$N_{\text{CAO}}$-by-$N_{\text{CMO}}^{\mathbf{k}}$ and
$N_{\text{CMO}}^{\mathbf{k}}$ is the dimension of the effective variational
space after removing numerical linear dependencies. We note that we then have

$(\mathbf{X}^{\mathbf{k}})^{\dagger}\mathbf{S}^{\mathbf{k}}\mathbf{X}^{\mathbf{k}}=\mathbf{I}_{N_{\text{CMO}}^{\mathbf{k}}}$
(9)

The choice of $\epsilon_{\text{lindep}}$ should be made so as to balance
between numerical stability (i.e., removing enough basis functions to avoid
excessive roundoff error and precision loss) and quality of the resulting
basis set (i.e., keeping as many basis functions as possible). We picked
$\epsilon_{\text{lindep}}$ to be $10^{-6}$ which is the default value of our
molecular computations in Q-Chem.Shao  _et al._ (2015); Epifanovsky  _et al._
(2021) We note that this linear dependency threshold is chosen to be
reasonable for double precision, and could be tightened up if one could afford
quadruple or higher precision arithmetic.

###  II.2 Review of the GPW algorithm

The GPW density fitting algorithm was first proposed by Hutter and co-
workersLippert  _et al._ (1997) and has been popularized via the
implementation in CP2K.VandeVondele  _et al._ (2005); Kühne  _et al._ (2020)
The central idea of the algorithm is that one employs planewaves as the
auxiliary basis set for density-fitting while using GTOs as the primary
computational basis set. This strategy is particularly well-suited for solid-
state calculations since periodic boundary conditions are naturally imposed
and planewave density fitting can be done efficiently. While applying GPW to
three-dimensional (3D) systems is the most straightforward, lower-dimensional
systems (0D, 1D, and 2D) need special attention to remove spurious image-image
interactions. The application of GPW was successfully carried out by Füsti-
Molnár and PulayFüsti-Molnar and Pulay (2002); Füsti-Molnár and Pulay (2002)
for molecules (i.e., 0D) where spurious image-image interactions were removed
exactly by using a truncated Coulomb potential. A similar idea can be
generalized to 1D and 2D.Rozzi  _et al._ (2006)

Among various terms in $\mathbf{F}^{\mathbf{k}}$, in this work, we focus on
the Coulomb matrix, $\mathbf{J}^{\mathbf{k}}$, because this contribution is
typically the computational bottleneck in pure DFT calculations. We want to
compute the Coulomb matrix element between a basis function $\phi_{\mu}$
located in a unit cell $\mathbf{R}=\mathbf{0}$ (denoted as
$\phi_{\mu}^{\mathbf{0}}$) and a basis function $\phi_{\nu}$ located in a unit
cell $\mathbf{R}$ (denoted as $\phi_{\nu}^{\mathbf{R}}$),

$\displaystyle J_{\mu\nu}^{\mathbf{0}\mathbf{R}}$
$\displaystyle\equiv\int_{\mathbf{r}}\phi_{\mu}^{\mathbf{0}}(\mathbf{r})V_{J}(\mathbf{r})\phi_{\nu}^{\mathbf{R}}(\mathbf{r})$
$\displaystyle=\sum_{\mathbf{R^{\prime}}}\int_{\mathbf{r}\in\mathbf{R^{\prime}}}\phi_{\mu}^{\mathbf{0}}(\mathbf{r})V_{J}(\mathbf{r}-\mathbf{R^{\prime}})\phi_{\nu}^{\mathbf{R}}(\mathbf{r})$
(10)

where $V_{J}(\mathbf{r})$ is the Coulomb potential defined as

$V_{J}(\mathbf{r})=\int_{\mathbf{r^{\prime}}}\frac{\rho(\mathbf{r^{\prime}})}{\mathbf{r}-\mathbf{r^{\prime}}}$
(11)

and we used the fact that $V_{J}(\mathbf{r})$ is periodic in the unit cell
displacements. We note that $\mathbf{r}\in\mathbf{R}^{\prime}$ implies that
the domain of the integration is restricted to the unit cell centered at
$\mathbf{R}^{\prime}$.

The evaluation of $V_{J}(\mathbf{r})$ can be done with $\mathcal{O}(N_{g}\log
N_{g})$ complexity via the fast Fourier transform (FFT) algorithm for discrete
Fourier transform where $N_{g}$ is the number of grid points within the
simulation cell. In reciprocal space,

$V_{J}(\mathbf{G})=\frac{4\pi}{|\mathbf{G}|^{2}}\tilde{\rho}(\mathbf{G})$ (12)

where

$\tilde{\rho}(\mathbf{G})=\frac{1}{\Omega}\int_{\mathbf{r}}\rho(\mathbf{r})e^{i\mathbf{G}\cdot\mathbf{r}}$
(13)

with $\Omega$ being the volume of the computational unit cell. Using these,
the GPW algorithm computes $\mathbf{J}^{\mathbf{k}}$ as follows:

  1. 1.

Compute $\rho(\mathbf{r})$ within a unit cell via

$\rho(\mathbf{r})=\frac{1}{N_{k}}\sum_{\mathbf{k}}\sum_{i}\sum_{\mu\nu}C_{\mu
i}^{\mathbf{k}}(C_{\nu
i}^{\mathbf{k}})^{*}\phi_{\mu}^{\mathbf{k}}(\mathbf{r})(\phi_{\nu}^{\mathbf{k}}(\mathbf{r}))^{*}$
(14)

where $N_{k}$ is the number of k-points.

  2. 2.

Fourier transform $\rho(\mathbf{r})$ to obtain $\tilde{\rho}(\mathbf{G})$.
This is the “density-fitting” step using a planewave auxiliary basis set.

  3. 3.

Compute the Coulomb potential in reciprocal space via Eq. 12 and inverse
Fourier transform to obtain $V_{J}(\mathbf{r})$. Note that we ignore the
$|\mathbf{G}|=0$ component.

  4. 4.

Compute $\mathbf{J}^{\mathbf{k}}$ via

$J_{\mu\nu}^{\mathbf{k}}=\int_{\mathbf{r}\in\text{U.C.}}(\phi_{\mu}^{\mathbf{k}}(\mathbf{r}))^{*}V_{J}(\mathbf{r})\phi_{\nu}^{\mathbf{k}}(\mathbf{r})$
(15)

where the quadrature is performed only within the unit cell (U.C.).

Our implementation computes $\phi_{\mu}^{\mathbf{k}}(\mathbf{r})$ once in the
beginning and stores these in memory. Therefore, our GPW implementation for
the J-build has $\mathcal{O}(N_{k}N_{g})$ storage cost (due to storing
$\phi_{\mu}^{\mathbf{k}}(\mathbf{r})$) and $\mathcal{O}(N_{k}N_{g}+N_{g}\log
N_{g})$ compute cost assuming sparsity of CAOs. Since $N_{g}$ scales with the
unit cell volume while $N_{k}$ does not, this algorithm approaches
$\mathcal{O}(N)$ scaling. Diagonalization is performed by dense linear algebra
with $\mathcal{O}(N^{3})$ scaling.

###  II.3 Summary of implementation of rVV10

Some of the more modern density functionals use the VV10 vdW correction, but
the cost of evaluating VV10 scales quadratically with system size. Using ideas
from the work of Román-Pérez and Soler,Román-Pérez and Soler (2009) Sabatini
and others proposed an alternative functional form called rVV10Sabatini  _et
al._ (2013) which can be implemented efficiently with linear complexity for
planewave codes while retaining similar accuracy as VV10. Subsequently, the
use of rVV10 was verified for combinatorially optimized density functionals
(B97M-V, $\omega$B97X-V, and $\omega$B97M-V) leading to B97M-rV,
$\omega$B97X-rV, and $\omega$B97M-rV.Mardirossian  _et al._ (2017) We are
interested in investigating the performance of these combinatorially optimized
functionals for band gaps so an efficient implementation of rVV10 is highly
desirable.

The rVV10 energy functional readsVydrov and Van Voorhis (2010); Sabatini  _et
al._ (2013)

$E_{\text{rVV10}}=E^{\text{local}}_{\text{rVV10}}+E^{\text{non-
local}}_{\text{rVV10}}$ (16)

where the local part can be absorbed into the local density approximation
terms and the non-local part poses implementational challenges with a naíve
quadratic scaling cost. The non-local contribution is defined as

$E^{\text{non-
local}}_{\text{rVV10}}=\frac{1}{2}\int_{\mathbf{r}}\int_{\mathbf{r^{\prime}}}\rho(\mathbf{r})\kappa(\mathbf{r})^{-3/2}\rho(\mathbf{r^{\prime}})\kappa(\mathbf{r^{\prime}})^{-3/2}\Phi(\mathbf{r,r^{\prime}})$
(17)

where $\rho(\mathbf{r})$ is the electron density, and the kernel
$\Phi(\mathbf{r},\mathbf{r}^{\prime})$ is

$\Phi(\mathbf{r,r^{\prime}})=\frac{-1.5}{(q(\mathbf{r})R^{2}+1)(q(\mathbf{r^{\prime}})R^{2}+1)(q(\mathbf{r})R^{2}+q(\mathbf{r^{\prime}})R^{2}+2)}$
(18)

with $R=|\mathbf{r}-\mathbf{r}^{\prime}|$. The remaining terms are

$q(\mathbf{r})=\kappa(\mathbf{r})^{-1}\sqrt{C|\frac{\nabla\rho(\mathbf{r})}{\rho(\mathbf{r})}|^{4}+\frac{4}{3}\pi\rho(\mathbf{r})},$
(19)

and

$\kappa(\mathbf{r})=1.5b\pi(\frac{\rho(\mathbf{r})}{9\pi})^{1/6}.$ (20)

The fixed parameters $b$ and $C$ are a part of the definition of each XC
functional that includes the rVV10 contribution. The evaluation of this leads
to an overall quadratic scaling in $N_{g}$ due to its six-dimensional double
integral in Eq. 17.

As discussed in ref. 92, we first use cubic splines to interpolate $\Phi$ such
that

$\Phi(\mathbf{r,r^{\prime}})\approx\sum_{\alpha,\beta}\Phi(q_{\alpha},q_{\beta},R)p_{\alpha}(q(\mathbf{r}))p_{\beta}(q(\mathbf{r^{\prime}}))$
(21)

where $q_{\alpha}$ and $q_{\beta}$ are interpolation points and $p_{\alpha}$
and $p_{\beta}$ are interpolating polynomials. This makes the evaluation of
$\Phi$ computationally convenient because $\Phi$ becomes a function of only
$|\mathbf{r-r^{\prime}}|$. Its dependence on $q_{\alpha}$ and $q_{\beta}$ is
easy to handle as $q_{\alpha}$ and $q_{\beta}$ are fixed interpolation points.
The number of the interpolation points is also very manageable as it is
typically set to 20.Sabatini  _et al._ (2013) We now define an intermediate,

$\theta_{\alpha}(\mathbf{r})=\rho(\mathbf{r})\kappa(\mathbf{r})^{-3/2}p_{\alpha}(q(\mathbf{r}))$
(22)

and use it to recast the non-local energy contribution into a convolution
form:

$\displaystyle E^{\text{non-local}}_{\text{rVV10}}$
$\displaystyle=\frac{1}{2}\sum_{\alpha,\beta}\int_{\mathbf{r}}\int_{\mathbf{r^{\prime}}}\theta_{\alpha}(\mathbf{r})\theta_{\beta}(\mathbf{r^{\prime}})\Phi_{\alpha\beta}(|\mathbf{r-r^{\prime}}|)$
$\displaystyle=\frac{1}{2}\sum_{\alpha,\beta}\int_{\mathbf{G}}\tilde{\theta}_{\alpha}(\mathbf{G})\tilde{\theta}_{\beta}(\mathbf{G})\tilde{\Phi}_{\alpha\beta}(|\mathbf{G}|)$
(23)

Since $\Phi_{\alpha\beta}(|\mathbf{r-r^{\prime}}|)=\Phi_{\alpha\beta}(R)$ is
spherically symmetric, its Fourier transform can be computed by one-
dimensional Fourier-sine transformation. The values of
$\tilde{\Phi}_{\alpha\beta}(|\mathbf{G}|)$ on a pre-defined set of
$|\mathbf{G}|$ points can be precomputed and tabulated. These tabulated values
are then used for interpolation to perform the convolution in Eq. 23 for a
specific set of $|\mathbf{G}|$. We note that the convolution in Eq. 23 can be
performed in $\mathcal{O}(N_{g}\log N_{g})$ time as opposed to the quadratic-
scaling runtime of the naíve algorithm. A similar approach can be used to
compute the Fock matrix contribution associated with rVV10.

###  II.4 Strategies for assessing the basis set error for simple solids and
generating transferable GTOs

Our goal in this work is to access the near basis set limit of pure density
functionals for solids using GTOs. For this purpose, it is critical to have
well-defined basis set limit reference values. A popular planewave code,
Quantum Espresso (QE), also implements the GTH pseudopotentialGoedecker  _et
al._ (1996) developed by Hutter and co-workers, which was originally developed
for use in the CP2K program.Kühne  _et al._ (2020) We have adopted the same
GTH pseudopotential for use in our code as well. This allows for a direct
comparison between QE and our code, which is particularly useful because QE
can converge the total energy to the basis set limit almostly completely by
increasing the planewave cutoff.

We considered the 40 semiconductors benchmark set (SC40) first proposed by
Scuseria and co-workersHeyd  _et al._ (2005) along with three rocksalt solids
(LiH,Vidal and Vidal-Valat (1986); Nolan  _et al._ (2009); Hoang and Van de
Walle (2013) LiF,Matsushita  _et al._ (2011); Garza and Scuseria (2016)
LiClMatsushita  _et al._ (2011); Garza and Scuseria (2016)). For these
compounds, all-electron GTO basis sets have been proposed but their accuracy
remains largely unknown.Heyd  _et al._ (2005) Moreover, to be used with the
GTH pseudopotential, we need a basis set without core electrons.
Unfortunately, the standard GTH basis set series does not have a broad
coverage of the periodic table beyond its minimal basis set (SZV-
GTH).VandeVondele and Hutter (2007) To access the basis set limit for a
variety of solids considered in this work, we propose a simple way to generate
a large basis set which yields the total energy per cell close to the basis
set limit (errors smaller than 0.7 $mE_{h}$ per atom for DFT calculations
performed here, as will be shown later). We also note that this same strategy
of uncontracting existing GTO bases can be applied to the generation of all-
electron bases as well.

To generate the basis set, we follow a straightforward procedure:

  1. 1.

We take the existing def2-bases and uncontract the contracted GTOs therein. We
then remove GTOs with an exponent greater than 20 since they correspond to
core electrons that are already covered by the GTH pseudopotential. This
cutoff of 20 was empirically determined and we expect that the results are not
sensitive to the precise value of the cutoff given the large size of our final
basis set (see below for more discussion).

  2. 2.

We take the union of these uncontracted def2 bases and the uncontracted SZV-
MOLOPT-SR-GTH basis set to enhance the resolution within the minimal basis set
space defined by the GTH pseudopotential. This final basis set will be
referred to as unc-def2-X-GTH where X can be SVP, TZVP, QZVP, etc.

One may think that having a fixed cutoff of 20 for all elements could be
unphysical because increasing the atomic number tends to increase all of GTO
exponents. In our case, however, the GTOs from def2 bases with an exponent
larger than the largest exponent in SZV-MOLOPT-SR-GTH belong to the core
region that is already treated by the GTH pseudopotential. Inspecting the
range of exponents in SZV-MOLOPT-SR-GTH basis, one finds that the largest ones
are smaller than 20 with the exception of Na (23.5) and Mg (30.7) up to atomic
number 86. Based on our results on solids involving Mg, the cutoff of 20 works
well for this element as well. Overall, the contraction of electron density
due to the increase in the nuclear charge is reflected appropriately and there
is no sensitivity stemming from this cutoff. We also note that when taking the
union of two bases some of the exponents can be very close in value, but for
simplicity we do not remove those obvious near-linear-dependencies. Instead we
let the canonical orthogonalization procedure take care of them. We report
these unc-def-GTH bases (unc-def2-SVP-GTH, unc-def2-SVPD-GTH, unc-def2-TZVP-
GTH, unc-def2-TZVPD-GTH, unc-def2-TZVPP-GTH, unc-def2-TZVPPD-GTH, unc-
def2-QZVP-GTH, unc-def2-QZVPD-GTH, unc-def2-QZVPP-GTH, unc-def2-QZVPPD-GTH)
through the Zenodo repository,zen (2021) as well as in the text files included
in the final publication.

With regard to the existing GTH-based contracted GTO basis sets, at present
neither the range of Gaussian exponents nor the contraction coefficients have
been specifically optimized to approach the basis set limit: rather they have
been designed to offer a good trade-off between compute cost and accuracy for
solid-state applications. The use of uncontracted basis functions in this work
is an attempt to probe the suitability of using a broad range of Gaussian
exponents and angular momenta while obtaining the contraction coefficients via
variational energy minimization (i.e., the MO coefficients are the contraction
coefficients in our case). As a consequence of decontraction, our proposed
basis sets range from quite large to very large and are heavily linearly
dependent. Nonetheless, this brute force approach will permit us to assess
systematic convergence of our total energies towards the basis set limit
energies obtained through QE. We emphasize that potential numerical
instability issues are quite well handled by the simple canonical
orthogonalization procedure.

Last but not least, we note that our Gaussian basis set generation procedure
does not utilize any system-dependent parameters or optimization protocols. As
evidenced by even-tempered bases,Bardo and Ruedenberg (1973); Feller and
Ruedenberg (1979) this is particularly important for ensuring transferability.
Our exponents retain both tight exponents that are effective for condensed
phase and relatively diffuse exponents that are effective for atomic (or
molecular) limits. Therefore, we expect that the BSIE is relatively
insensitive to the underlying geometry. Nonetheless, when atoms come too close
together, GTOs are expected to perform more poorly due to the higher degree of
linear dependence as will be shown later in Section V.3.

##  III Computational details

We implemented a GPW-based periodic DFT code in a development version of
Q-Chem.Shao  _et al._ (2015); Epifanovsky  _et al._ (2021) Our implementation
assumes overall $\mathcal{O}(N^{2})$ memory storage that amounts to storing
Fock, density, and CMO coefficient matrices. For the systems examined in this
work, our memory strorage is dominated by keeping the GTO basis function
values evaluated on the FFT grid despite its formal linear scaling based on
sparsity. In the future, this practical memory bottleneck can be removed by
computing these on the fly. We also note that our GPW algorithms scale
linearly with the system size to produce the Fock matrix and our SCF program
scales cubically with system size due to linear algebra functions such as
matrix diagonalization. We control the resolution of the PW density fitting
basis with a single parameter: the kinetic energy cutoff ($E_{\text{cut}}$).
Each auxiliary basis PW can be indexed by 3 integers, $(n_{1},n_{2},n_{3})$
which reside on a
$(2n_{1}^{\text{max}}-1)\times(2n_{2}^{\text{max}}-1)\times(2n_{3}^{\text{max}}-1)$
grid where each integer
$n_{i}\in\\{-n_{i}^{\text{max}},\cdots,n_{i}^{\text{max}}\\}$ with

$n_{i}^{\text{max}}=\frac{\sqrt{8E_{\text{cut}}}}{||\mathbf{b}_{i}||}$ (24)

where $\mathbf{b}_{i}$ denotes one of the reciprocal vectors. For our GPW
calculations, we used $E_{\text{cut}}$ of 1500 eV for everything except those
that contain Ba (2000 eV) and Mg (4500 eV). The resulting density fitting
error was found to be smaller than 100 $\mu E_{h}$ per atom in the unit cell,
which is negligible for the purpose of this paper.

The reference planewave basis calculations were all performed with QE where we
used $E_{\text{cut}}$ (for the wavefunction itself) of 1200 Ry for total
energy calculations. For the band structure calculations, we used
$E_{\text{cut}}$ of 1200 Ry for systems containing Mg and $E_{\text{cut}}$ of
750 Ry for everything else.

The lattice constants were fixed at experimental values Heyd  _et al._ (2005)
and experimental band gaps for the SC40 set were taken from 61. The
experimental band gaps and lattice constants of LiH, LiF, and LiCl were taken
from refs. 62; 63; 64; 65; 66. We used the GTH-LDA pseudopotential in all
calculations for both GPW and PW (through QE) calculations to enable direct
comparison of total energies. We used the Monkhorst-PackMonkhorst and Pack
(1976) $\mathbf{k}$-mesh to sample the first Brillouin zone and ensured the
convergence of the total energy per cell to the TDL for all solids examined
here. We found that a $6\times 6\times 6$ $\mathbf{k}$-mesh is enough to
converge the total energy per cell to an error of smaller than 0.1 mH for all
solids considered. Therefore, for band structure calculations and cold curve
calculations, we used a $6\times 6\times 6$ Monkhorst-Pack $\mathbf{k}$-mesh.
Since the GPW implementation is also available in other open-source packages
such as CP2KKühne  _et al._ (2020) and PySCF,Sun  _et al._ (2020) we also used
these two packages to validate our implementation in the initial stage of this
work.

We examined a total of ten XC functionals, LDA (Slater exchangeSlater (1951)
and PZ81 correlationPerdew and Zunger (1981)), PBE,Perdew  _et al._ (1996)
PBEsol,Perdew  _et al._ (2008) revPBE,Zhang and Yang (1998) BLYP,Becke (1988);
Lee  _et al._ (1988) B97-D,Grimme (2006) SCAN,Sun  _et al._ (2015) M06-L,Zhao
and Truhlar (2006) MN15-L,Yu  _et al._ (2016) and B97M-rV.Mardirossian and
Head-Gordon (2015); Mardirossian  _et al._ (2017) For 11 solids in our
benchmark set (C, Si, SiC, BN, BP, AlN, MgO, MgS, LiH, LiF, LiCl), widely used
GTH basis setsVandeVondele and Hutter (2007) are available: DZVP-GTH, TZVP-
GTH, TZV2P-GTH, QZV2P-GTH, and QZV3P-GTH. We therefore assessed the accuracy
of those existing bases only over a smaller subset of our benchmark set, but
our proposed basis sets were examined for all 43 solids considered in this
work. The basis set convergence study against PW was carried out only for LDA
and PBE while the overall band gap accuracy was examined for all ten
functionals.

##  IV Results and Discussion

###  IV.1 Basis set convergence of total DFT energies

Figure 1: Root mean square deviation (RMSD) of DFT total energies (m$E_{h}$)
per cell with respect to that of QE over 11 solids as a function of the number
of $\mathbf{k}$-points for (a) LDA and (b) PBE functionals using GTH and unc-
def2-GTH bases.

We first examine the subset of 11 solids for LDA and PBE functionals as
presented in Fig. 1. In particular, Fig. 1 shows the root-mean-square-
deviation (RMSD) of total energies compared to QE total energies (namely total
energies in the basis set limit) as function of the size of the
$\mathbf{k}$-mesh. $N_{k}=216$ ($6\times 6\times 6$) is enough to reach the
TDL.  While smaller $\mathbf{k}$-mesh calculations such as the $\Gamma$-point
only calculations are unphysical, we are interested in how the basis set error
changes as a function of the $\mathbf{k}$-mesh size. Since both PW and our GTO
results share the same Hamiltonian for each size of $\mathbf{k}$-mesh, we can
assess the basis set error of GTO calculations as long as the reference PW
energies are fully converged to the basis set limit. For all $\mathbf{k}$-mesh
sizes, the GTH basis set series shows systematically more accurate results
relative to the basis set limit as cardinality (and the size of the basis set)
increases. We note that $N_{k}=1$ (just including the $\Gamma$-point) shows
the largest basis set error in all examples considered here. This is because
the local expressive power of GTOs also increases as one increases the size of
$\mathbf{k}$-mesh due to the non-orthogonality of GTOs. While the systematic
improvement of GTH bases is very appealing, we note that the residual basis
set error with QZV3P-GTH is still about 5 m$E_{h}$ which is quite large
considering how small the simulation cells are (2 or 4 atoms total).

The unc-def2-GTH basis series also shows a systematic improvement with
cardinal number, with much smaller errors than the corresponding contracted
GTH basis results. As an example, the performance of unc-def2-SVP-GTH is
nearly on par with TZV2P-GTH except at the $\Gamma$-point. The larger bases,
unc-def2-TZVP-GTH and unc-def2-QZVP-GTH, both perform excellently on this set,
including the $\Gamma$-point result. In particular, unc-def2-QZVP-GTH is able
to deliver total energies in the TDL that are all within 1 m$E_{h}$ of the
basis set limit. This shows the completeness of our proposed bases though of
course these are much bigger in size than standard GTH bases, and therefore
far more computationally demanding. We provide more detailed information on
selected elements in Section V.1. Finally, we note that Fig. 1 (a) for LDA and
(b) for PBE show virtually no difference, which suggests that our conclusions
do not depend on functional (of course functionals that depend particularly
strongly on fine details of the density may be far harder to converge to the
basis set limit using GTOsMardirossian and Head-Gordon (2013)).

Figure 2: Root mean square deviation (RMSD) of DFT total energies (m$E_{h}$)
per cell with respect to that of QE over 43 solids as a function of the number
of $\mathbf{k}$-points for (a) LDA and (b) PBE functionals using unc-def2-GTH
bases.

In Fig. 2, we repeat the same analysis but over the entire benchmark set of 43
solids. As before, unc-def2-GTH bases struggle for $N_{k}=1$ but work well for
larger $\mathbf{k}$-meshes. RMSD systematically decreases as we increase the
size of the basis set. With the largest basis set, unc-def2-QZVP-GTH, we
achieve better than 1 m$E_{h}$ accuracy in the TDL for the LDA and PBE
functionals, as measured by the RMSD values. Systems with the largest error in
the TDL are SrSe (1.2 m$E_{h}$) in the case of LDA and GaP (1.4 m$E_{h}$) in
the case of PBE. As observed in the case of even-tempered bases, we expect
that the result can be systematically made better by adding more exponents and
increasing the maximum angular momentum.Bardo and Ruedenberg (1973); Feller
and Ruedenberg (1979) For instance, in the case of SrSe/LDA, employing unc-
def2-QZVPP-GTH (adding two additional f functions to both Sr and Se), we
observe an error of 0.4 m$E_{h}$ which is three times smaller than that of
unc-def2-QZVP-GTH. While we can obtain overall better results by using unc-
def2-QZVPP-GTH, we will mainly focus on the use of the unc-def2-QZVP-GTH basis
set for the rest of the paper for simplicity. In summary, these benchmark
calculations suggest that unc-def2-GTH basis sets can achieve near basis set
limit DFT total energies reliably towards the TDL. This result implies that
the range of exponents and angular momenta in our bases is quite appropriate
for solids.

###  IV.2 Basis set convergence of DFT band gaps

In many materials applications, DFT calculations are used not just to compute
the ground state energy but to obtain spectral information through Kohn-Sham
orbital energies.Perdew  _et al._ (2017) In doing so, one uses information
from virtual orbitals in addition to that from occupied orbitals. In the case
of total energies presented in Section IV.1, only occupied orbitals affect the
results. Here, we are assessing the quality of the difference between the
lowest energy virtual orbital (i.e., the conduction band minimum) and the
higher energy occupied orbital. It is possible that some BSIEs may cancel when
taking energy differences.

Figure 3: Root mean square deviation (RMSD) of DFT band gaps (meV) with
respect to those of QE over 11 solids (a) LDA and (b) PBE functionals using
GTH and unc-def2-GTH bases.

In Fig. 3, we present the RMSD of band gaps using GTO bases compared to the
basis set limit results for LDA and PBE functionals. To compare unc-def2-GTH
bases with GTH bases, we limit ourselves to the subset of 11 solids for the
time being. Somewhat surprisingly, we observe almost no improvement in the
band gap when going from DZVP-GTH to TZVP-GTH. By contrast, Fig. 1 shows a
reduction in the RMSD of total energies of about 8 m$E_{h}$ when increasing
the basis set size from DZVP-GTH to TZVP-GTH. However, this total energy
improvement does not result in any band gap improvement. Nonetheless, past
TZVP-GTH, the GTH bases do show systematic improvement in the band gap
estimation. With the largest GTH basis set (QZV3P-GTH), RMSD in the band gap
is 18-20 meV depending on the XC functional. Consistent with the total energy
benchmark presented in Section IV.1, unc-def2-GTH bases also exhibit
systematic improvement. While the quality of unc-def2-SVP-GTH was on par with
TZV2P-GTH in Fig. 1, its band gap is clearly worse than that of TZV2P-GTH
highlighting favorable error cancellation in TZV2P-GTH. Nonetheless, unc-
def2-TZVP-GTH is similar to QZV3P-GTH and unc-def2-QZVP-GTH has RMSD of 5.8
meV and 4.2 meV, respectively for LDA and PBE, showing its ability to converge
band gaps to the basis set limit.

Figure 4: Root mean square deviation (RMSD) of DFT band gaps (meV) with
respect to those of QE over 43 solids (a) LDA and (b) PBE functionals using
unc-def2-GTH bases.

Encouraged by these results, we also analyzed the BSIEs in band gaps over all
43 solids using unc-def2-GTH bases as presented in Fig. 4. With unc-def2-SVP-
GTH, the RMSD value is about 40 meV and it becomes less than 20 meV when using
unc-def2-TZVP-GTH. Lastly, with unc-def2-QZVP-GTH, the RMSD value becomes 6.9
meV and 6.3 meV, respectively, for LDA and PBE. However, we note that the
largest deviation is about 20 meV in both functionals, which corresponds to
the band gap of SrSe. SrSe is the system with the largest total energy error
for LDA as noted in the discussion of Fig. 2 in Section IV.1. Again, this
remaining error can be further reduced by adding more GTOs to unc-def2-QZVP-
GTH (e.g., using unc-def2-QZVPP-GTH), but we do not pursue this here. The
central message of this section is that the BSIE in the band gap reported in
this paper using unc-def2-QZVP-GTH is smaller than 20 meV based on the
numerical data. This is about 50 times smaller than the intrinsic errors in
standard functionals for band gaps, so unc-def2-QZVP-GTH should be suitable
for benchmarking purposes.

###  IV.3 Performance of pure DFT functionals

| LDA | PBE | PBEsol | revPBE | BLYP | B97-D | SCAN | M06-L | MN15-L | B97M-rV | Exp.
---|---|---|---|---|---|---|---|---|---|---|---
C | 4.12 | 4.33 | 4.16 | 4.38 | 4.60 | 4.57 | 4.64 | 4.84 | 4.24 | 4.67 | 5.48
Si | 0.49 | 0.66 | 0.52 | 0.72 | 0.94 | 0.91 | 0.93 | 1.12 | 0.96 | 0.92 | 1.17
Ge | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.20 | 0.46 | 0.43 | 0.26 | 0.74
SiC | 1.33 | 1.50 | 1.35 | 1.54 | 1.85 | 1.95 | 1.81 | 1.83 | 1.92 | 2.07 | 2.42
BN | 4.36 | 4.64 | 4.42 | 4.71 | 5.03 | 5.07 | 5.03 | 4.94 | 4.98 | 5.28 | 6.22
BP | 1.20 | 1.38 | 1.23 | 1.42 | 1.66 | 1.62 | 1.67 | 1.95 | 1.61 | 1.70 | 2.4
BAs | 1.16 | 1.34 | 1.19 | 1.40 | 1.60 | 1.59 | 1.57 | 1.85 | 1.55 | 1.66 | 1.46
BSb | 0.76 | 0.91 | 0.78 | 0.96 | 1.15 | 1.17 | 1.06 | 1.21 | 1.08 | 1.15 | N/A
AlP | 1.47 | 1.67 | 1.50 | 1.75 | 1.98 | 2.04 | 1.99 | 2.20 | 2.12 | 2.16 | 2.51
AlAs | 1.36 | 1.58 | 1.40 | 1.66 | 1.89 | 1.95 | 1.86 | 2.00 | 1.99 | 2.03 | 2.23
AlSb | 1.17 | 1.36 | 1.20 | 1.45 | 1.57 | 1.58 | 1.56 | 1.78 | 1.63 | 1.63 | 1.68
bGaN | 1.61 | 1.79 | 1.68 | 1.85 | 1.86 | 1.96 | 1.86 | 1.88 | 1.43 | 1.86 | 3.3
GaP | 1.44 | 1.66 | 1.50 | 1.75 | 1.70 | 1.72 | 1.85 | 1.89 | 1.84 | 2.05 | 2.35
GaAs | 0.29 | 0.51 | 0.41 | 0.58 | 0.44 | 0.46 | 0.62 | 0.92 | 0.72 | 1.01 | 1.52
GaSb | 0.00 | 0.17 | 0.08 | 0.23 | 0.08 | 0.08 | 0.21 | 0.50 | 0.36 | 0.64 | 0.73
InP | 0.42 | 0.61 | 0.52 | 0.68 | 0.57 | 0.56 | 0.57 | 0.86 | 0.36 | 0.88 | 1.42
InAs | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.03 | 0.41
InSb | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.13 | 0.23
ZnS | 1.80 | 2.12 | 1.94 | 2.24 | 2.15 | 2.24 | 2.38 | 2.56 | 2.00 | 2.36 | 3.66
ZnSe | 1.01 | 1.33 | 1.15 | 1.45 | 1.35 | 1.43 | 1.60 | 1.80 | 1.31 | 1.67 | 2.7
ZnTe | 1.07 | 1.38 | 1.22 | 1.49 | 1.38 | 1.43 | 1.60 | 1.83 | 1.40 | 1.79 | 2.38
CdS | 0.85 | 1.15 | 0.98 | 1.27 | 1.16 | 1.22 | 1.22 | 1.34 | 0.77 | 1.23 | 2.55
CdSe | 0.34 | 0.64 | 0.48 | 0.76 | 0.66 | 0.70 | 0.72 | 0.88 | 0.34 | 0.82 | 1.9
CdTe | 0.52 | 0.81 | 0.66 | 0.93 | 0.81 | 0.82 | 0.84 | 1.07 | 0.52 | 1.06 | 1.92
MgS | 3.27 | 3.57 | 3.39 | 3.73 | 3.71 | 3.88 | 4.14 | 4.03 | 4.06 | 4.17 | 5.4
MgTe | 2.30 | 2.59 | 2.43 | 2.75 | 2.72 | 2.89 | 3.16 | 3.18 | 2.99 | 3.23 | 3.6
MgO | 4.68 | 4.93 | 4.79 | 5.06 | 5.17 | 5.49 | 5.59 | 5.01 | 5.69 | 5.55 | 7.22
MgSe | 1.70 | 2.01 | 1.89 | 2.16 | 2.10 | 2.47 | 2.58 | 2.66 | 3.20 | 2.70 | 2.47
CaS | 2.17 | 2.41 | 2.27 | 2.51 | 2.52 | 2.64 | 2.92 | 2.56 | 3.03 | 3.05 | N/A
CaSe | 1.90 | 2.14 | 2.00 | 2.24 | 2.27 | 2.38 | 2.63 | 2.28 | 2.73 | 2.78 | N/A
CaTe | 1.42 | 1.65 | 1.51 | 1.74 | 1.80 | 1.90 | 2.08 | 1.74 | 2.19 | 2.27 | N/A
SrS | 2.22 | 2.49 | 2.33 | 2.61 | 2.62 | 2.74 | 2.92 | 2.57 | 2.86 | 2.95 | N/A
SrSe | 2.01 | 2.28 | 2.12 | 2.40 | 2.43 | 2.54 | 2.69 | 2.35 | 2.63 | 2.74 | N/A
SrTe | 1.57 | 1.83 | 1.66 | 1.94 | 2.00 | 2.10 | 2.20 | 1.87 | 2.17 | 2.29 | N/A
BaS | 2.01 | 2.26 | 2.11 | 2.38 | 2.36 | 2.44 | 2.58 | 2.25 | 2.38 | 2.51 | 3.88
BaSe | 1.83 | 2.07 | 1.92 | 2.19 | 2.19 | 2.26 | 2.39 | 2.07 | 2.21 | 2.35 | 3.58
BaTe | 1.49 | 1.74 | 1.58 | 1.85 | 1.87 | 1.93 | 2.03 | 1.74 | 1.90 | 2.03 | 3.08
LiH | 2.64 | 3.01 | 2.78 | 3.15 | 3.44 | 3.69 | 3.61 | 3.87 | 4.52 | 4.39 | 4.9
LiF | 8.92 | 9.33 | 9.11 | 9.56 | 9.49 | 9.92 | 10.08 | 9.64 | 10.27 | 9.77 | 14.2
LiCl | 6.01 | 6.40 | 6.18 | 6.61 | 6.56 | 6.83 | 7.21 | 7.08 | 7.48 | 7.22 | 9.4
AlN | 4.25 | 4.38 | 4.25 | 4.43 | 4.66 | 4.78 | 4.87 | 4.75 | 4.94 | 5.24 | 6.13
GaN | 1.86 | 2.05 | 1.94 | 2.12 | 2.12 | 2.21 | 2.12 | 2.13 | 1.68 | 2.11 | 3.5
InN | 0.00 | 0.00 | 0.00 | 0.01 | 0.00 | 0.04 | 0.00 | 0.00 | 0.00 | 0.00 | 0.69
RMSD | 1.72 | 1.50 | 1.64 | 1.41 | 1.39 | 1.28 | 1.20 | 1.26 | 1.27 | 1.17 | N/A
MAD | -1.46 | -1.23 | -1.37 | -1.14 | -1.10 | -1.02 | -0.95 | -0.90 | -0.99 | -0.84 | N/A
MAX | 5.28 | 4.87 | 5.09 | 4.64 | 4.71 | 4.28 | 4.12 | 4.56 | 3.93 | 4.43 | N/A
Table 1:  Experimental and theoretical band gaps (or fundamental gaps) (eV)
from various functionals over 43 solids. N/A means “not available”. RMSD, MAD,
and MAX denote, respectively, root-mean-square-deviation, mean-average-
deviation, and maximum deviation in reference to experimental values. We took
experimental references for three rocksalt solids (LiH,LiF,LiCl) from refs 62;
63; 64, 65; 66, and 65; 66, respectively. The rest of experimental values were
taken from ref. 61. Figure 5: Band gap (eV) comparison over 36 solids between
DFT (ten different functionals) and experiments: Blue: root-mean-square-
deviation (RMSD) of DFT band gaps (eV) with respect to those of experiments
and Red: mean-average-deviation (MAD) of DFT band gaps (eV) with respect to
those of experiments

Having established the accuracy of unc-def2-GTH bases, we assess the
performance of pure DFT functionals over these simple solids. Unfortunately,
some of the 43 solids considered here do not have experimental band gaps.
These solids are BSb, CaS, CaSe, CaTe, SrS, SrSe, and SrTe. Leaving aside
these seven cases, we have a total of 36 experimental band gaps. unc-
def2-QZVP-GTH is used with all XC functionals considered in this section. The
DFT band gaps over 43 solids along with the available experimental gaps are
presented in Table 1. For an overall summary, it may be more instructuve to
look at statistics of the band gap results as shown in Fig. 5. Looking at the
mean-average-deviation (MAD), it is immediately evident that all pure
functionals examined here exhibit the infamous band gap underestimation
problem of pure functionals.Perdew  _et al._ (2017) Nonetheless, one can still
find systematic improvement for going from the simplest functional, LDA, to
more modern meta GGA functionals, SCAN, M06-L, MN15-L, and B97M-rV in terms of
the root-mean-square-deviation (RMSD) values. While the performance of B97M-rV
is not great for those band gaps, it still stays as one of the more accurate
pure functionals for these problems. This is encouraging because B97M-rV is
statistically the most accurate pure XC functional in main group chemistry
applications.Mardirossian and Head-Gordon (2017) Overall, all pure functionals
perform poorly in this benchmark study and the inclusion of exact exchange
seems necessary. This is not a new observation on its own and has been well-
documented even for nearly the same benchmark set that we study here.Civalleri
_et al._ (2012); Heyd  _et al._ (2005) In ref. 97, the authors consider many
different ways to analyze the statistical data of LDA, PBE, and PBEsol along
with other hybrid functionals and revealed that hybrid functionals always
perform the best in nine out of ten statistical analyses. It will be
interesting to revisit this benchmark set with modern hybrid functionals in
the future.

For simplicity and due to the unavailability of functional-specific GTH
pseudopotentials for most XC functionals considered here, we employed the GTH-
LDA pseudopotential for all functionals in this section. Since this is not
ideal, we checked the sensitivity of our conclusions with respect to the
choice of the pseudopotential by testing GTH-LDA, GTH-PBE, and GTH-BLYP
pseudopotentials with the BLYP functional. In all cases the RMSD and MAD are
affected by less than 0.1 eV, which is a smaller energy scale than that of the
band gap errors by roughly a factor of 10. For completeness, we provide the
relevant numerical data in the Supplementary Material (see Table S1).In the
future, all-electron calculations could be done with all-electron basis sets
generated via a similar protocol presented here. Alternatively, one could
generate functional-specific GTH pseudopotentials for the modern XC
functionals considered here.

##  V Outlook for future basis set design

In this section, we would like to deliver cautionary notes on using our
proposed bases and some discussion on future research in basis set design for
solids.

###  V.1 Our basis set is accurate but very large

| Si | C | O | Mg
---|---|---|---|---
SZV-GTH | 4 | 4 | 4 | 5
DZVP-GTH | 13 | 13 | 13 | 14
TZVP-GTH | 17 | 17 | 17 | 18
TZV2P-GTH | 22 | 22 | 22 | 23
QZV2P-GTH | 26 | 26 | 26 | 27
QZV3P-GTH | 31 | 31 | 31 | 32
unc-def2-SVP-GTH | 40 | 41 | 40 | 53
unc-def2-TZVP-GTH | 62 | 58 | 57 | 68
unc-def2-QZVP-GTH | 90 | 83 | 81 | 86
Table 2: Number of basis functions in the basis sets used in this work for
selected elements (Si, C, O, and Mg).

While our proposed unc-def2-GTH bases are of high quality, these bases are
very large due to the decontraction from the original contracted GTO bases.
This large size carries a significant computational cost. This is the major
drawback of even-tempered and well-tempered bases, and it is one that our unc-
def2-GTH bases also share. To be more concrete, we provide the number of basis
functions for selected elements (Si, C, O, Mg) in Table 2. unc-def2-SVP-GTH is
about three times bigger than DZVP-GTH while our unc-def2-TZVP-GTH is roughly
three times bigger than TZV2P-GTH. Similarly, our largest basis set unc-
def2-QZVP-GTH is about 2.5–3 times larger than QZV3P-GTH.

Because of compute cost and memory demand, there is a need to compress these
bases for practical calculations. Perhaps, the most difficult (but most
effective if done correctly) way to compress them is to obtain transferable
contraction coefficients. One could start by inspecting the molecular orbitals
(or Bloch orbitals) that our calculations produce for those simple solids.
Another strategy is to take these mean-field molecular orbitals and compress
the virtual space for subsequent correlation calculations, for instance using
the random phase approximation (RPA). The use of natural orbitals to compress
the virtual space was shown to be effective, and would be a good starting
point for making our basis more compact Ramberger  _et al._ (2019). We note
that it is also unclear whether our proposed bases exhibit any scaling
properties which will allow for higher accuracy by using basis set
extrapolation for correlation energy calculations, which could be further
investigated in the future.

###  V.2 Even low-lying virtual orbitals can be difficult to describe well

Figure 6: LDA band structure of MgO: (a) comparing QZV3P-GTH against QE and
(b) comparing unc-def2-QZVP-GTH against QE. The band energies are shifted such
that the highest valence band energy is located at zero. The red circle in (a)
highlights the qualitative failure of QZV3P-GTH virtual orbitals.

Basis sets that are optimized for mean-field calculations such as GTH bases
often behave erratically in correlated calculations.Morales and Malone (2020)
Since these bases tend to yield good occupied orbitals, the poor performance
of correlation calculations can be attributed to virtual orbitals.
Furthermore, low-lying virtual orbitals play important roles in describing
optical properties and related excited states. Therefore, high-quality basis
sets should produce qualitatively accurate virtuals. As an example, we present
the band structure of MgO using QZV3P-GTH and unc-def2-QZVP-GTH and compare
them against that of QE. MgO has a total of 8 occupied orbitals and we
computed up to the 16-th band in QE for comparison purposes. We note that the
challenge of MgO conduction bands for GTOs was noted before in ref. 99, but we
focus on a wider range of conduction bands here. The pertinent band structures
are presented in Fig. 6.

In both bases, the valence bands and the first few conduction bands are in an
excellent agreement with those of QE. However, the higher-lying virtuals of
QZV3P-GTH (in Fig. 6(a)) start to deviate significantly from those of QE. The
most striking failure is the lack of the 5-th virtual orbital highlighted
under a red circle in Fig. 6 (a). On the other hand, the virtuals from unc-
def2-QZVP-GTH have visually indistinguishable energies when compared to QE
highlighting its potential utility for correlated calculations as well. We
also tried a smaller unc-def2-GTH basis set, namely unc-def2-TZVP-GTH. It
turns out that even unc-def2-TZVP-GTH misses the same virtual that QZV3P-GTH
misses as well.  With further investigations, we found that the 5-th virtual
orbital is missed by basis sets without any f function on the Mg atom. To be
more concrete, we added one f-function to Mg in the QZV3P-GTH basis where the
exponent of 0.16 was taken from unc-def2-QZVP-GTH. With this basis set, we
recover the missing band at the $\Gamma$-point. This additional f-function
introduces only a 0.2 m$E_{h}$ energy lowering in the ground state, but it is
essential to capture one of the low-lying conduction bands.

This example emphasizes that more attention to the low-lying virtual orbitals
should be paid when designing GTO basis sets for applications such as
conduction band structure, time-dependent DFT and correlated methods such as
RPA. Existing GTO bases designed primarily to describe the occupied space may
likewise exhibit qualitative failures like this case.

###  V.3 Transferability across different lattice constants is challenging

Cold curves of solids are often of great interest for experimentalists. Cold
curves are analogous to potential energy curves (PECs) in molecular quantum
chemistry. Similar to PECs, as one shrinks the lattice constant and brings
atoms close to one another, a larger number of near linear dependencies occur,
and the quality of the underlying GTO basis degrades because of discarding
such functions by canonical orthogonalization. Furthermore, system-dependent
optimization strategies can struggle for cold curves because basis sets are
usually optimized for one specific geometry (usually equilibrium
geometries).Morales and Malone (2020); Daga  _et al._ (2020) As a result of
this, varying lattice constants can be challenging using GTO basis sets as the
system approaches its high-pressure configuration (shorter lattice constants)
or atomic limits (longer lattice constants).

Figure 7: Investigation of a SiC cold curve using PBE: (a) Basis set error
with respect to QE as a function of lattice constant for various basis sets.
(b) Average number of linearly independent AOs as a function of the lattice
constant for various basis sets.

As an example to illustrate this point, we computed a cold curve of SiC using
PBE and the GTH-LDA pseudopotential with TZV2P-GTH, QZV3P-GTH, unc-QZV3P-GTH,
unc-def2-TZVP, and unc-def2-QZVP. The Brillouin zone was sampled with $6\times
6\times 6$ $\mathbf{k}$-points via the Monkhorst-Pack algorithm. Here unc-
QZV3P-GTH is the basis obtained from decontracting QZV3P-GTH. Using unc-
QZV3P-GTH, we can quantify the error coming from the contraction coefficients
of QZV3P-GTH. As before, the QE results (the same functional and
pseudopotential) serve as the basis set limit reference values. The pertinent
results are presented in Fig. 7.

Fig. 7 (a) shows that TZV2P-GTH, QZV3P-GTH, and unc-QZV3P-GTH bases make a
large error especially when compressing the lattice. It is also instructive to
quantify the nonparallelity error (NPE) over the range of lattice constants
examined here as a means to measure error cancellation. The NPE is 25.7
m$E_{h}$, 10.6 m$E_{h}$, 8.3 m$E_{h}$, 3.0 m$E_{h}$, and 0.9 m$E_{h}$,
respectively, for TZV2P-GTH, QZV3P-GTH, unc-QZV3P-GTH, unc-def2-TZVP-GTH, and
unc-def2-QZVP-GTH. Interestingly, unc-QZV3P-GTH reduces the basis set error by
only a small amount, which implies that the contraction coefficients of
QZV3P-GTH for those elements are transferable over a wide range of lattice
constants. It also suggests that the range of exponents in QZV3P-GTH becomes
inappropriate for smaller lattice constants. Comparing the exponents of unc-
QZV3P-GTH and unc-def2-TZVP-GTH, we find that unc-def2-TZVP-GTH has more
compact exponents for spd shells and has an f shell for Si that is not present
in unc-QZV3P-GTH. These more compact GTOs likely become more important at
closer distances and hence explain the differences between two bases.

The main cause of these generally large NPEs is the fact that at closer
distances the quality of those GTO bases degrades as shown in Fig. 7 (b) which
quantifies the number of orthogonalized basis functions retained after
canonical orthogonalization. Since each $\mathbf{k}$-point has a slightly
different number of linearly independent AOs, we present the average values
over 216 $\mathbf{k}$-points as a function of lattice constant. Evidently, the
number of linearly independent AOs decreases as the lattice constant
decreases, which in turn increases the basis set incompleteness error.
Nonetheless, the largest basis set, unc-def2-QZVP-GTH, is able to achieve a
satisfying NPE in this case, which highlights the utility of this basis set
for accurate cold curve simulations.

##  VI Conclusions

In this manuscript, we discussed strategies for generating large and accurate
uncontracted Gaussian bases (unc-def2-GTH bases) which do not resort to
system- or method- specific optimizations. Using a new implementation of the
Gaussian atomic orbital plus planewave density fitting approach in Q-Chem, the
basis set incompleteness error in our proposed bases were then assessed over
43 prototypical semiconductors by comparing the pure density functional theory
total energies per cell and band gaps against those from fully converged
planewave results. We found that the basis set incompleteness error in total
energy and band gap with our largest GTO basis (unc-def2-QZVP-GTH) is smaller
than 0.7 m$E_{h}$ per atom in the unit cell and less than 20 meV,
respectively, verifying the validity of the range of exponents and angular
momenta in the proposed bases.

In the application of our bases, we focused on the assessment of ten pure
density functionals for predicting the band gaps of 36 semiconductors whose
experimental gaps are well documented. Not surprisingly, we found that all
examined pure functionals (LDA, PBE, PBEsol, revPBE, BLYP, B97-D, SCAN, M06-L,
MN15-L, B97M-rV) significantly underestimate the band gaps of these materials.
The combinatorially optimized mGGA functional, B97M-rV, performs as well as do
other modern mGGA functionals. Our work suggests that combinatorially
optimized range-separated hybrid functionals such as $\omega$B97X-rV and
$\omega$B97M-rV will be highly interesting to study since they may also
exhibit better accuracy compared to other relatively older range-separated
hybrid functionals or even short-range hybrid functionals.

We also made several cautionary remarks on our GTO bases as well as on the
future research in GTO basis design for solids:

  1. 1.

Our basis sets are accurate but large so there is a need for a way to compress
our basis sets further for both mean-field and correlation calculations.

  2. 2.

The widely used GTH bases may qualitatively fail for describing low-lying
virtual orbitals which will affect the subsequent correlation and optical
calculations. At much greater compute cost, our unc-def2-QZVP-GTH basis set
was shown to accurately capture all of the low-lying virtual orbitals of MgO
including the one missed by QZV3P-GTH.

  3. 3.

Reducing the non-parallelity error of the basis set incompleteness error is
challenging particularly due to the high pressure region of cold curves that
exhibits a higher number of near linear dependencies.

In the future, we will test several ways (e.g., finding universal contraction
coefficients and frozen natural orbitalsRamberger  _et al._ (2019)) to
compress our unc-def2-GTH bases and investigate the basis set convergence of
correlation and optical methods with these bases in the future. Furthermore,
simple solids like LiH have many numerical data of the total Hartree-Fock
energies towards the basis set limit,Gillan  _et al._ (2008); Marsman  _et
al._ (2009); Paier  _et al._ (2009); Civalleri  _et al._ (2010); Usvyat  _et
al._ (2011) which could be a good testbed for our basis sets in the future.

##  VII Supplementary Material

The supplementary material of this work is available, which contains the test
of the impact of different GTH pseudopotentials in the band gaps.

##  VIII Acknowledgement

This work was supported by the National Institutes of Health SBIR program
through Grant No. 2R44GM128480-02A1. We thank Eloy Ramos for initial
engagement with this project in 2017 and Kuan-Yu Liu for help with the
implementation of the $\mathbf{k}$-point parser used in this work. JL thanks
David Reichman for support.

##  IX Data Availability Statement

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

##  X Conflict of Interest

E.E. and M.H.-G. are part-owners of Q-Chem, Inc.

## References

  * Ayala  _et al._ (2001) Philippe Y. Ayala, Konstantin N. Kudin, and Gustavo E. Scuseria, “Atomic orbital Laplace-transformed second-order Møller–Plesset theory for periodic systems,” J. Chem. Phys. 115, 9698–9707 (2001).
  * Katagiri (2005) Hideki Katagiri, “Equation-of-motion coupled-cluster study on exciton states of polyethylene with periodic boundary condition,” J. Chem. Phys. 122, 224901 (2005).
  * Izmaylov and Scuseria (2008) Artur F. Izmaylov and Gustavo E. Scuseria, “Resolution of the identity atomic orbital Laplace transformed second order Møller–Plesset theory for nonconducting periodic systems,” Phys. Chem. Chem. Phys. 10, 3421–3429 (2008).
  * Pisani  _et al._ (2008) Cesare Pisani, Lorenzo Maschio, Silvia Casassa, Migen Halo, Martin Schütz, and Denis Usvyat, “Periodic local MP2 method for the study of electronic correlation in crystals: Theory and preliminary applications,” J. Comput. Chem. 29, 2113–2124 (2008).
  * Pisani  _et al._ (2012) Cesare Pisani, Martin Schütz, Silvia Casassa, Denis Usvyat, Lorenzo Maschio, Marco Lorenz, and Alessandro Erba, “C ryscor : a program for the post-Hartree–Fock treatment of periodic systems,” Phys. Chem. Chem. Phys. 14, 7615–7628 (2012).
  * Del Ben  _et al._ (2012) Mauro Del Ben, Jürg Hutter, and Joost VandeVondele, “Second-Order Møller–Plesset Perturbation Theory in the Condensed Phase: An Efficient and Massively Parallel Gaussian and Plane Waves Approach,” J. Chem. Theory Comput. 8, 4177–4188 (2012).
  * Usvyat  _et al._ (2015) Denis Usvyat, Lorenzo Maschio, and Martin Schütz, “Periodic local MP2 method employing orbital specific virtuals,” J. Chem. Phys. 143, 102805 (2015).
  * McClain  _et al._ (2017) James McClain, Qiming Sun, Garnet Kin-Lic Chan, and Timothy C. Berkelbach, “Gaussian-Based Coupled-Cluster Theory for the Ground-State and Band Structure of Solids,” J. Chem. Theory Comput. 13, 1209–1218 (2017).
  * Shang and Yang (2020) Honghui Shang and Jinlong Yang, “Implementation of Laplace Transformed MP2 for Periodic Systems With Numerical Atomic Orbitals,” Front. Chem. 8 (2020), 10.3389/fchem.2020.589992.
  * Maurer  _et al._ (2019) Reinhard J. Maurer, Christoph Freysoldt, Anthony M. Reilly, Jan Gerit Brandenburg, Oliver T. Hofmann, Torbjörn Björkman, Sébastien Lebègue, and Alexandre Tkatchenko, “Advances in Density-Functional Calculations for Materials Modeling,” Annu. Rev. Mater. Res. 49, 1–30 (2019).
  * Hybertsen and Louie (1986) Mark S. Hybertsen and Steven G. Louie, “Electron correlation in semiconductors and insulators: Band gaps and quasiparticle energies,” Phys. Rev. B 34, 5390–5413 (1986).
  * Aryasetiawan and Gunnarsson (1998) F. Aryasetiawan and O. Gunnarsson, “The GW method,” Rep. Prog. Phys. 61, 237–312 (1998).
  * Löwdin (1970) Per-Olov Löwdin, “On the Nonorthogonality Problem,” in  _Advances in Quantum Chemistry_, Vol. 5 (Academic Press, Cambridge, MA, USA, 1970) pp. 185–199.
  * Klahn and Bingel (1977) Bruno Klahn and Werner A. Bingel, “Completeness and linear independence of basis sets used in quantum chemistry,” Int. J. Quantum Chem. 11, 943–957 (1977).
  * VandeVondele and Hutter (2007) Joost VandeVondele and Jürg Hutter, “Gaussian basis sets for accurate calculations on molecular systems in gas and condensed phases,” J. Chem. Phys. 127, 114105 (2007).
  * Peintinger  _et al._ (2013) Michael F. Peintinger, Daniel Vilela Oliveira, and Thomas Bredow, “Consistent Gaussian basis sets of triple-zeta valence with polarization quality for solid-state calculations,” J. Comput. Chem. 34, 451–459 (2013).
  * Laun  _et al._ (2018) Joachim Laun, Daniel Vilela Oliveira, and Thomas Bredow, “Consistent gaussian basis sets of double- and triple-zeta valence with polarization quality of the fifth period for solid-state calculations,” J. Comput. Chem. 39, 1285–1290 (2018).
  * Oliveira  _et al._ (2019) Daniel Vilela Oliveira, Joachim Laun, Michael F. Peintinger, and Thomas Bredow, “BSSE-correction scheme for consistent gaussian basis sets of double- and triple-zeta valence with polarization quality for solid-state calculations,” J. Comput. Chem. 40, 2364–2376 (2019).
  * Morales and Malone (2020) Miguel A. Morales and Fionn D. Malone, “Accelerating the convergence of auxiliary-field quantum Monte Carlo in solids with optimized Gaussian basis sets,” J. Chem. Phys. 153, 194111 (2020).
  * Daga  _et al._ (2020) Loredana Edith Daga, Bartolomeo Civalleri, and Lorenzo Maschio, “Gaussian Basis Sets for Crystalline Solids: All-Purpose Basis Set Libraries vs System-Specific Optimizations,” J. Chem. Theory Comput. 16, 2192–2201 (2020).
  * Zhou  _et al._ (2021) Yanbing Zhou, Emanuel Gull, and Dominika Zgid, “Material-Specific Optimization of Gaussian Basis Sets against Plane Wave Data,” J. Chem. Theory Comput. 17, 5611–5622 (2021).
  * Li  _et al._ (1995) Y. S. Li, M. C. Wrinn, J. M. Newsam, and M. P. Sears, “Parallel implementation of a mesh-based density functional electronic structure code,” J. Comput. Chem. 16, 226–234 (1995).
  * Soler  _et al._ (2002) José M. Soler, Emilio Artacho, Julian D. Gale, Alberto García, Javier Junquera, Pablo Ordejón, and Daniel Sánchez-Portal, “The SIESTA method for ab initio order-N,” J. Phys.: Condens. Matter 14, 2745–2779 (2002).
  * Artacho  _et al._ (2008) Emilio Artacho, E. Anglada, O. Diéguez, J. D. Gale, A. García, J. Junquera, R. M. Martin, P. Ordejón, J. M. Pruneda, D. Sánchez-Portal, and J. M. Soler, “The SIESTA method; developments and applicability,” J. Phys.: Condens. Matter 20, 064208 (2008).
  * Verzijl and Thijssen (2012) C. J. O. Verzijl and J. M. Thijssen, “DFT-Based Molecular Transport Implementation in ADF/BAND,” J. Phys. Chem. C 116, 24393–24412 (2012).
  * Kudin and Scuseria (2000) Konstantin N. Kudin and Gustavo E. Scuseria, “Linear-scaling density-functional theory with Gaussian orbitals and periodic boundary conditions: Efficient evaluation of energy and forces via the fast multipole method,” Phys. Rev. B 61, 16440–16453 (2000).
  * Balasubramani  _et al._ (2020) Sree Ganesh Balasubramani, Guo P. Chen, Sonia Coriani, Michael Diedenhofen, Marius S. Frank, Yannick J. Franzke, Filipp Furche, Robin Grotjahn, Michael E. Harding, Christof Hättig, Arnim Hellweg, Benjamin Helmich-Paris, Christof Holzer, Uwe Huniar, Martin Kaupp, Alireza Marefat Khah, Sarah Karbalaei Khani, Thomas Müller, Fabian Mack, Brian D. Nguyen, Shane M. Parker, Eva Perlt, Dmitrij Rappoport, Kevin Reiter, Saswata Roy, Matthias Rückert, Gunnar Schmitz, Marek Sierka, Enrico Tapavicza, David P. Tew, Christoph van Wüllen, Vamsee K. Voora, Florian Weigend, Artur Wodyński, and Jason M. Yu, “TURBOMOLE: Modular program suite for ab initio quantum-chemical and condensed-matter simulations,” J. Chem. Phys. 152, 184107 (2020).
  * Dovesi  _et al._ (2020) Roberto Dovesi, Fabien Pascale, Bartolomeo Civalleri, Klaus Doll, Nicholas M. Harrison, Ian Bush, Philippe D’Arco, Yves Noël, Michel Rérat, Philippe Carbonnière, Mauro Causà, Simone Salustro, Valentina Lacivita, Bernard Kirtman, Anna Maria Ferrari, Francesco Silvio Gentile, Jacopo Baima, Mauro Ferrero, Raffaella Demichelis, and Marco De La Pierre, “The CRYSTAL code, 1976–2020 and beyond, a long story,” J. Chem. Phys. 152, 204111 (2020).
  * Kühne  _et al._ (2020) Thomas D. Kühne, Marcella Iannuzzi, Mauro Del Ben, Vladimir V. Rybkin, Patrick Seewald, Frederick Stein, Teodoro Laino, Rustam Z. Khaliullin, Ole Schütt, Florian Schiffmann, Dorothea Golze, Jan Wilhelm, Sergey Chulkov, Mohammad Hossein Bani-Hashemian, Valéry Weber, Urban Borštnik, Mathieu Taillefumier, Alice Shoshana Jakobovits, Alfio Lazzaro, Hans Pabst, Tiziano Müller, Robert Schade, Manuel Guidon, Samuel Andermatt, Nico Holmberg, Gregory K. Schenter, Anna Hehn, Augustin Bussy, Fabian Belleflamme, Gloria Tabacchi, Andreas Glöß, Michael Lass, Iain Bethune, Christopher J. Mundy, Christian Plessl, Matt Watkins, Joost VandeVondele, Matthias Krack, and Jürg Hutter, “CP2K: An electronic structure and molecular dynamics software package - Quickstep: Efficient and accurate electronic structure calculations,” J. Chem. Phys. 152, 194103 (2020).
  * Prentice  _et al._ (2020) Joseph C. A. Prentice, Jolyon Aarons, James C. Womack, Alice E. A. Allen, Lampros Andrinopoulos, Lucian Anton, Robert A. Bell, Arihant Bhandari, Gabriel A. Bramley, Robert J. Charlton, Rebecca J. Clements, Daniel J. Cole, Gabriel Constantinescu, Fabiano Corsetti, Simon M.-M. Dubois, Kevin K. B. Duff, José María Escartín, Andrea Greco, Quintin Hill, Louis P. Lee, Edward Linscott, David D. O’Regan, Maximillian J. S. Phipps, Laura E. Ratcliff, Álvaro Ruiz Serrano, Edward W. Tait, Gilberto Teobaldi, Valerio Vitale, Nelson Yeung, Tim J. Zuehlsdorff, Jacek Dziedzic, Peter D. Haynes, Nicholas D. M. Hine, Arash A. Mostofi, Mike C. Payne, and Chris-Kriton Skylaris, “The ONETEP linear-scaling density functional theory program,” J. Chem. Phys. 152, 174111 (2020).
  * Jensen (2013) Frank Jensen, “Atomic orbital basis sets,” WIRES: Comput. Mol. Sci. 3, 273–295 (2013).
  * Nagy and Jensen (2017) Balazs Nagy and Frank Jensen, “Basis sets in quantum chemistry,” Rev. Computat. Chem. 30, 93–149 (2017).
  * Shavitt (1993) Isaiah Shavitt, “The History and Evolution of Gaussian Basis Sets,” Isr. J. Chem. 33, 357–367 (1993).
  * Pritchard  _et al._ (2019) Benjamin P. Pritchard, Doaa Altarawy, Brett Didier, Tara D. Gibson, and Theresa L. Windus, “New Basis Set Exchange: An Open, Up-to-Date Resource for the Molecular Sciences Community,” J. Chem. Inf. Model. 59, 4814–4820 (2019).
  * Mcweeny (1950) R. Mcweeny, “Gaussian Approximations, to Wave Functions - Nature,” Nature 166, 21–22 (1950).
  * Boys S. (1950) F. Boys S., “Electronic wave functions - I. A general method of calculation for the stationary states of any molecular system,” Proc. R. Soc. London A - Math. Phys. Sci. 200, 542–554 (1950).
  * Almlöf and Taylor (1987) Jan Almlöf and Peter R. Taylor, “General contraction of Gaussian basis sets. I. Atomic natural orbitals for first- and second-row atoms,” J. Chem. Phys. 86, 4070–4077 (1987).
  * Dunning (1989) Thom H. Dunning, “Gaussian basis sets for use in correlated molecular calculations. I. The atoms boron through neon and hydrogen,” J. Chem. Phys. 90, 1007–1023 (1989).
  * Jensen (2001) Frank Jensen, “Polarization consistent basis sets: Principles,” J. Chem. Phys. 115, 9113–9125 (2001).
  * Jensen (2002) Frank Jensen, “Polarization consistent basis sets. ii. estimating the kohn–sham basis set limit,” J. Chem. Phys. 116, 7372–7379 (2002).
  * Bardo and Ruedenberg (1973) Richard D. Bardo and Klaus Ruedenberg, “Even-tempered atomic orbitals. III. Economic deployment of Gaussian primitives in expanding atomic SCF orbitals,” J. Chem. Phys. 59, 5956–5965 (1973).
  * Feller and Ruedenberg (1979) David F. Feller and Klaus Ruedenberg, “Systematic approach to extended even-tempered orbital bases for atomic and molecular calculations,” Theor. Chim. Acta 52, 231–251 (1979).
  * Huzinaga  _et al._ (1985) S. Huzinaga, M. Klobukowski, and H. Tatewaki, “The well-tempered GTF basis sets and their applications in the SCF calculations on N2, CO, Na2, and P2,” Can. J. Chem. 63, 1812 (1985).
  * Weigend and Ahlrichs (2005) Florian Weigend and Reinhart Ahlrichs, “Balanced basis sets of split valence, triple zeta valence and quadruple zeta valence quality for H to Rn: Design and assessment of accuracy,” Phys. Chem. Chem. Phys. 7, 3297–3305 (2005).
  * Goedecker  _et al._ (1996) S. Goedecker, M. Teter, and J. Hutter, “Separable dual-space Gaussian pseudopotentials,” Phys. Rev. B 54, 1703–1710 (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–3662 (1998).
  * Giannozzi  _et al._ (2020) Paolo Giannozzi, Oscar Baseggio, Pietro Bonfà, Davide Brunato, Roberto Car, Ivan Carnimeo, Carlo Cavazzoni, Stefano De Gironcoli, Pietro Delugas, Fabrizio Ferrari Ruffino, _et al._ , “Quantum espresso toward the exascale,” J. Chem. Phys. 152, 154105 (2020).
  * Slater (1951) John C Slater, “A simplification of the hartree-fock method,” Phys. Rev. 81, 385 (1951).
  * Perdew and Zunger (1981) John P Perdew and Alex Zunger, “Self-interaction correction to density-functional approximations for many-electron systems,” Phys. Rev. B 23, 5048 (1981).
  * Perdew  _et al._ (1996) John P Perdew, Kieron Burke, and Matthias Ernzerhof, “Generalized gradient approximation made simple,” Phys. Rev. Lett. 77, 3865 (1996).
  * Perdew  _et al._ (2008) John P. Perdew, Adrienn Ruzsinszky, Gábor I. Csonka, Oleg A. Vydrov, Gustavo E. Scuseria, Lucian A. Constantin, Xiaolan Zhou, and Kieron Burke, “Restoring the Density-Gradient Expansion for Exchange in Solids and Surfaces,” Phys. Rev. Lett. 100, 136406 (2008).
  * Zhang and Yang (1998) Yingkai Zhang and Weitao Yang, “Comment on “Generalized Gradient Approximation Made Simple”,” Phys. Rev. Lett. 80, 890 (1998).
  * Becke (1988) Axel D Becke, “Density-functional exchange-energy approximation with correct asymptotic behavior,” Phys. Rev. A 38, 3098 (1988).
  * Lee  _et al._ (1988) Chengteh Lee, Weitao Yang, and Robert G Parr, “Development of the colle-salvetti correlation-energy formula into a functional of the electron density,” Phys. Rev. B 37, 785 (1988).
  * Grimme (2006) Stefan Grimme, “Semiempirical GGA-type density functional constructed with a long-range dispersion correction,” J. Comput. Chem. 27, 1787–1799 (2006).
  * Sun  _et al._ (2015) Jianwei Sun, Adrienn Ruzsinszky, and John P Perdew, “Strongly constrained and appropriately normed semilocal density functional,” Phys. Rev. Lett. 115, 036402 (2015).
  * Zhao and Truhlar (2006) Yan Zhao and Donald G. Truhlar, “A new local density functional for main-group thermochemistry, transition metal bonding, thermochemical kinetics, and noncovalent interactions,” J. Chem. Phys. 125, 194101 (2006).
  * Yu  _et al._ (2016) Haoyu S. Yu, Xiao He, and Donald G. Truhlar, “MN15-L: A New Local Exchange-Correlation Functional for Kohn–Sham Density Functional Theory with Broad Accuracy for Atoms, Molecules, and Solids,” J. Chem. Theory Comput. 12, 1280–1293 (2016).
  * Mardirossian and Head-Gordon (2015) Narbe Mardirossian and Martin Head-Gordon, “Mapping the genome of meta-generalized gradient approximation density functionals: The search for b97m-v,” J. Chem. Phys. 142, 074111 (2015).
  * Mardirossian  _et al._ (2017) Narbe Mardirossian, Luis Ruiz Pestana, James C Womack, Chris-Kriton Skylaris, Teresa Head-Gordon, and Martin Head-Gordon, “Use of the rvv10 nonlocal correlation functional in the b97m-v density functional: Defining b97m-rv and related functionals,” J. Phys. Chem. Lett. 8, 35–40 (2017).
  * Heyd  _et al._ (2005) Jochen Heyd, Juan E Peralta, Gustavo E Scuseria, and Richard L Martin, “Energy band gaps and lattice parameters evaluated with the heyd-scuseria-ernzerhof screened hybrid functional,” J. Chem. Phys. 123, 174101 (2005).
  * Vidal and Vidal-Valat (1986) JEAN PIERRE Vidal and G Vidal-Valat, “Accurate debye–waller factors of 7lih and 7lid by neutron diffraction at three temperatures,” Acta Crystallogr., Sect. B: Struct. Sci. 42, 131–137 (1986).
  * Nolan  _et al._ (2009) SJ Nolan, MJ Gillan, D Alfè, NL Allan, and FR Manby, “Calculation of properties of crystalline lithium hydride using correlated wave function theory,” Phys. Rev. B 80, 165109 (2009).
  * Hoang and Van de Walle (2013) Khang Hoang and Chris G Van de Walle, “Lih as a li+ and h- ion provider,” Solid State Ion. 253, 53–56 (2013).
  * Matsushita  _et al._ (2011) Yu-ichiro Matsushita, Kazuma Nakamura, and Atsushi Oshiyama, “Comparative study of hybrid functionals applied to structural and electronic properties of semiconductors and insulators,” Phys. Rev. B 84, 075205 (2011).
  * Garza and Scuseria (2016) Alejandro J Garza and Gustavo E Scuseria, “Predicting band gaps with hybrid density functionals,” J. Phys. Chem. Lett. 7, 4165–4170 (2016).
  * Perdew  _et al._ (2017) John P Perdew, Weitao Yang, Kieron Burke, Zenghui Yang, Eberhard KU Gross, Matthias Scheffler, Gustavo E Scuseria, Thomas M Henderson, Igor Ying Zhang, Adrienn Ruzsinszky, _et al._ , “Understanding band gaps of solids in generalized kohn–sham theory,” Proc. Natl. Acad. Sci. U.S.A 114, 2801–2806 (2017).
  * Yang  _et al._ (2016) Zeng-hui Yang, Haowei Peng, Jianwei Sun, and John P Perdew, “More realistic band gaps from meta-generalized gradient approximations: Only in a generalized kohn-sham scheme,” Phys. Rev. B 93, 205205 (2016).
  * Mardirossian and Head-Gordon (2014) Narbe Mardirossian and Martin Head-Gordon, “$\omega$b97x-v: A 10-parameter, range-separated hybrid, generalized gradient approximation density functional with nonlocal correlation, designed by a survival-of-the-fittest strategy,” Phys. Chem. Chem. Phys. 16, 9904–9924 (2014).
  * Mardirossian and Head-Gordon (2016) Narbe Mardirossian and Martin Head-Gordon, “$\omega$ b97m-v: A combinatorially optimized, range-separated hybrid, meta-gga density functional with vv10 nonlocal correlation,” J. Chem. Phys. 144, 214110 (2016).
  * Pestana  _et al._ (2017) Luis Ruiz Pestana, Narbe Mardirossian, Martin Head-Gordon, and Teresa Head-Gordon, “Ab initio molecular dynamics simulations of liquid water using high quality meta-gga functionals,” Chem. Sci. 8, 3554–3565 (2017).
  * Ruiz Pestana  _et al._ (2018) Luis Ruiz Pestana, Ondrej Marsalek, Thomas E Markland, and Teresa Head-Gordon, “The quest for accurate liquid water properties from first principles,” J. Phys. Chem. Lett. 9, 5009–5016 (2018).
  * Pestana  _et al._ (2019) Luis Ruiz Pestana, Hongxia Hao, and Teresa Head-Gordon, “Diels–alder reactions in water are determined by microsolvation,” Nano Lett. 20, 606–611 (2019).
  * Lininger  _et al._ (2021) Christianna N Lininger, Joseph A Gauthier, Wan-Lu Li, Elliot Rossomme, Valerie Vaissier Welborn, Zhou Lin, Teresa Head-Gordon, Martin Head-Gordon, and Alexis T Bell, “Challenges for density functional theory: calculation of co adsorption on electrocatalytically relevant metals,” Phys. Chem. Chem. Phys. 23, 9394–9406 (2021).
  * Li  _et al._ (2021) Wan-Lu Li, Christianna N Lininger, Kaixuan Chen, Valerie Vaissier Welborn, Elliot Rossomme, Alexis T Bell, Martin Head-Gordon, and Teresa Head-Gordon, “Critical role of thermal fluctuations for co binding on electrocatalytic metal surfaces,” JACS Au , (in press) (2021).
  * Sabatini  _et al._ (2013) Riccardo Sabatini, Tommaso Gorni, and Stefano De Gironcoli, “Nonlocal van der waals density functional made simple and efficient,” Phys. Rev. B 87, 041108 (2013).
  * Vydrov and Van Voorhis (2010) Oleg A Vydrov and Troy Van Voorhis, “Nonlocal van der waals density functional: The simpler the better,” J. Chem. Phys. 133, 244103 (2010).
  * Mardirossian and Head-Gordon (2017) Narbe Mardirossian and Martin Head-Gordon, “Thirty years of density functional theory in computational chemistry: an overview and extensive assessment of 200 density functionals,” Mol. Phys. 115, 2315–2372 (2017).
  * Goerigk  _et al._ (2017) L. Goerigk, A. Hansen, C. Bauer, S. Ehrlich, A. Najibi, and S. Grimme, “A look at the density functional theory zoo with the advanced GMTKN55 database for general main group thermochemistry, kinetics and noncovalent interactions,” Phys. Chem. Chem. Phys. 19, 32184–32215 (2017).
  * Najibi and Goerigk (2018) A. Najibi and L. Goerigk, “The nonlocal kernel in van der Waals density functionals as an additive correction: An extensive analysis with special emphasis on the B97M-V and $\omega$B97M-V approches,” J. Chem. Theory Comput. 14, 5725–5738 (2018).
  * Evarestov (2007) Robert Evarestov,  _Quantum Chemistry of Solids_ (Springer-Verlag, Berlin, Germany, 2007).
  * Evarestov (2012) R. A. Evarestov,  _Quantum Chemistry of Solids_ (Springer-Verlag, Berlin, Germany, 2012).
  * Kudin and Scuseria (1998) Konstantin N Kudin and Gustavo E Scuseria, “A fast multipole algorithm for the efficient treatment of the coulomb problem in electronic structure calculations of periodic systems with gaussian orbitals,” Chem. Phys. Lett. 289, 611–616 (1998).
  * Szabo and Ostlund (1996) Attila Szabo and Neil S. Ostlund,  _Modern Quantum Chemistry: Introduction to Advanced Electronic Structure Theory_ (Courier Corporation, 1996).
  * Shao  _et al._ (2015) Yihan Shao, Zhengting Gan, Evgeny Epifanovsky, Andrew T. B. Gilbert, Michael Wormit, Joerg Kussmann, Adrian W. Lange, Andrew Behn, Jia Deng, Xintian Feng, Debashree Ghosh, Matthew Goldey, Paul R. Horn, Leif D. Jacobson, Ilya Kaliman, Rustam Z. Khaliullin, Tomasz Kuś, Arie Landau, Jie Liu, Emil I. Proynov, Young Min Rhee, Ryan M. Richard, Mary A. Rohrdanz, Ryan P. Steele, Eric J. Sundstrom, H. Lee Woodcock, Paul M. Zimmerman, Dmitry Zuev, Ben Albrecht, Ethan Alguire, Brian Austin, Gregory J. O. Beran, Yves A. Bernard, Eric Berquist, Kai Brandhorst, Ksenia B. Bravaya, Shawn T. Brown, David Casanova, Chun-Min Chang, Yunqing Chen, Siu Hung Chien, Kristina D. Closser, Deborah L. Crittenden, Michael Diedenhofen, Robert A. DiStasio, Hainam Do, Anthony D. Dutoi, Richard G. Edgar, Shervin Fatehi, Laszlo Fusti-Molnar, An Ghysels, Anna Golubeva-Zadorozhnaya, Joseph Gomes, Magnus W. D. Hanson-Heine, Philipp H. P. Harbach, Andreas W. Hauser, Edward G. Hohenstein, Zachary C. Holden, Thomas-C. Jagau, Hyunjun Ji, Benjamin Kaduk, Kirill Khistyaev, Jaehoon Kim, Jihan Kim, Rollin A. King, Phil Klunzinger, Dmytro Kosenkov, Tim Kowalczyk, Caroline M. Krauter, Ka Un Lao, Adèle D. Laurent, Keith V. Lawler, Sergey V. Levchenko, Ching Yeh Lin, Fenglai Liu, Ester Livshits, Rohini C. Lochan, Arne Luenser, Prashant Manohar, Samuel F. Manzer, Shan-Ping Mao, Narbe Mardirossian, Aleksandr V. Marenich, Simon A. Maurer, Nicholas J. Mayhall, Eric Neuscamman, C. Melania Oana, Roberto Olivares-Amaya, Darragh P. O’Neill, John A. Parkhill, Trilisa M. Perrine, Roberto Peverati, Alexander Prociuk, Dirk R. Rehn, Edina Rosta, Nicholas J. Russ, Shaama M. Sharada, Sandeep Sharma, David W. Small, Alexander Sodt, Tamar Stein, David Stück, Yu-Chuan Su, Alex J. W. Thom, Takashi Tsuchimochi, Vitalii Vanovschi, Leslie Vogt, Oleg Vydrov, Tao Wang, Mark A. Watson, Jan Wenzel, Alec White, Christopher F. Williams, Jun Yang, Sina Yeganeh, Shane R. Yost, Zhi-Qiang You, Igor Ying Zhang, Xing Zhang, Yan Zhao, Bernard R. Brooks, Garnet K. L. Chan, Daniel M. Chipman, Christopher J. Cramer, William A. Goddard, Mark S. Gordon, Warren J. Hehre, Andreas Klamt, Henry F. Schaefer, Michael W. Schmidt, C. David Sherrill, Donald G. Truhlar, Arieh Warshel, Xin Xu, Alán Aspuru-Guzik, Roi Baer, Alexis T. Bell, Nicholas A. Besley, Jeng-Da Chai, Andreas Dreuw, Barry D. Dunietz, Thomas R. Furlani, Steven R. Gwaltney, Chao-Ping Hsu, Yousung Jung, Jing Kong, Daniel S. Lambrecht, WanZhen Liang, Christian Ochsenfeld, Vitaly A. Rassolov, Lyudmila V. Slipchenko, Joseph E. Subotnik, Troy Van Voorhis, John M. Herbert, Anna I. Krylov, Peter M. W. Gill, and Martin Head-Gordon, “Advances in molecular quantum chemistry contained in the Q-Chem 4 program package,” Mol. Phys. 113, 184–215 (2015).
  * Epifanovsky  _et al._ (2021) Evgeny Epifanovsky, Andrew T. B. Gilbert, Xintian Feng, Joonho Lee, Yuezhi Mao, Narbe Mardirossian, Pavel Pokhilko, Alec F. White, Marc P. Coons, Adrian L. Dempwolff, Zhengting Gan, Diptarka Hait, Paul R. Horn, Leif D. Jacobson, Ilya Kaliman, Jörg Kussmann, Adrian W. Lange, Ka Un Lao, Daniel S. Levine, Jie Liu, Simon C. McKenzie, Adrian F. Morrison, Kaushik D. Nanda, Felix Plasser, Dirk R. Rehn, Marta L. Vidal, Zhi-Qiang You, Ying Zhu, Bushra Alam, Benjamin J. Albrecht, Abdulrahman Aldossary, Ethan Alguire, Josefine H. Andersen, Vishikh Athavale, Dennis Barton, Khadiza Begam, Andrew Behn, Nicole Bellonzi, Yves A. Bernard, Eric J. Berquist, Hugh G. A. Burton, Abel Carreras, Kevin Carter-Fenk, Romit Chakraborty, Alan D. Chien, Kristina D. Closser, Vale Cofer-Shabica, Saswata Dasgupta, Marc de Wergifosse, Jia Deng, Michael Diedenhofen, Hainam Do, Sebastian Ehlert, Po-Tung Fang, Shervin Fatehi, Qingguo Feng, Triet Friedhoff, James Gayvert, Qinghui Ge, Gergely Gidofalvi, Matthew Goldey, Joe Gomes, Cristina E. González-Espinoza, Sahil Gulania, Anastasia O. Gunina, Magnus W. D. Hanson-Heine, Phillip H. P. Harbach, Andreas Hauser, Michael F. Herbst, Mario Hernández Vera, Manuel Hodecker, Zachary C. Holden, Shannon Houck, Xunkun Huang, Kerwin Hui, Bang C. Huynh, Maxim Ivanov, Ádám Jász, Hyunjun Ji, Hanjie Jiang, Benjamin Kaduk, Sven Kähler, Kirill Khistyaev, Jaehoon Kim, Gergely Kis, Phil Klunzinger, Zsuzsanna Koczor-Benda, Joong Hoon Koh, Dimitri Kosenkov, Laura Koulias, Tim Kowalczyk, Caroline M. Krauter, Karl Kue, Alexander Kunitsa, Thomas Kus, István Ladjánszki, Arie Landau, Keith V. Lawler, Daniel Lefrancois, Susi Lehtola, Run R. Li, Yi-Pei Li, Jiashu Liang, Marcus Liebenthal, Hung-Hsuan Lin, You-Sheng Lin, Fenglai Liu, Kuan-Yu Liu, Matthias Loipersberger, Arne Luenser, Aaditya Manjanath, Prashant Manohar, Erum Mansoor, Sam F. Manzer, Shan-Ping Mao, Aleksandr V. Marenich, Thomas Markovich, Stephen Mason, Simon A. Maurer, Peter F. McLaughlin, Maximilian F. S. J. Menger, Jan-Michael Mewes, Stefanie A. Mewes, Pierpaolo Morgante, J. Wayne Mullinax, Katherine J. Oosterbaan, Garrette Paran, Alexander C. Paul, Suranjan K. Paul, Fabijan Pavošević, Zheng Pei, Stefan Prager, Emil I. Proynov, Ádám Rák, Eloy Ramos-Cordoba, Bhaskar Rana, Alan E. Rask, Adam Rettig, Ryan M. Richard, Fazle Rob, Elliot Rossomme, Tarek Scheele, Maximilian Scheurer, Matthias Schneider, Nickolai Sergueev, Shaama M. Sharada, Wojciech Skomorowski, David W. Small, Christopher J. Stein, Yu-Chuan Su, Eric J. Sundstrom, Zhen Tao, Jonathan Thirman, Gábor J. Tornai, Takashi Tsuchimochi, Norm M. Tubman, Srimukh Prasad Veccham, Oleg Vydrov, Jan Wenzel, Jon Witte, Atsushi Yamada, Kun Yao, Sina Yeganeh, Shane R. Yost, Alexander Zech, Igor Ying Zhang, Xing Zhang, Yu Zhang, Dmitry Zuev, Alán Aspuru-Guzik, Alexis T. Bell, Nicholas A. Besley, Ksenia B. Bravaya, Bernard R. Brooks, David Casanova, Jeng-Da Chai, Sonia Coriani, Christopher J. Cramer, György Cserey, A. Eugene DePrince, Robert A. DiStasio, Andreas Dreuw, Barry D. Dunietz, Thomas R. Furlani, William A. Goddard, Sharon Hammes-Schiffer, Teresa Head-Gordon, Warren J. Hehre, Chao-Ping Hsu, Thomas-C. Jagau, Yousung Jung, Andreas Klamt, Jing Kong, Daniel S. Lambrecht, WanZhen Liang, Nicholas J. Mayhall, C. William McCurdy, Jeffrey B. Neaton, Christian Ochsenfeld, John A. Parkhill, Roberto Peverati, Vitaly A. Rassolov, Yihan Shao, Lyudmila V. Slipchenko, Tim Stauch, Ryan P. Steele, Joseph E. Subotnik, Alex J. W. Thom, Alexandre Tkatchenko, Donald G. Truhlar, Troy Van Voorhis, Tomasz A. Wesolowski, K. Birgitta Whaley, H. Lee Woodcock, Paul M. Zimmerman, Shirin Faraji, Peter M. W. Gill, Martin Head-Gordon, John M. Herbert, and Anna I. Krylov, “Software for the frontiers of quantum chemistry: An overview of developments in the Q-Chem 5 package,” J. Chem. Phys. 155, 084801 (2021).
  * Lippert  _et al._ (1997) Gerald Lippert, Jürg Hutter, and Michele Parrinello, “A hybrid gaussian and plane wave density functional scheme,” Mol. Phys. 92, 477–488 (1997).
  * VandeVondele  _et al._ (2005) Joost VandeVondele, Matthias Krack, Fawzi Mohamed, Michele Parrinello, Thomas Chassaing, and Jürg Hutter, “Quickstep: Fast and accurate density functional calculations using a mixed gaussian and plane waves approach,” Comput. Phys. Commun. 167, 103–128 (2005).
  * Füsti-Molnar and Pulay (2002) László Füsti-Molnar and Peter Pulay, “Accurate molecular integrals and energies using combined plane wave and Gaussian basis sets in molecular electronic structure theory,” J. Chem. Phys. 116, 7795–7805 (2002).
  * Füsti-Molnár and Pulay (2002) László Füsti-Molnár and Peter Pulay, “The Fourier transform Coulomb method: Efficient and accurate calculation of the Coulomb operator in a Gaussian basis,” J. Chem. Phys. 117, 7827–7835 (2002).
  * Rozzi  _et al._ (2006) Carlo A. Rozzi, Daniele Varsano, Andrea Marini, Eberhard K. U. Gross, and Angel Rubio, “Exact Coulomb cutoff technique for supercell calculations,” Phys. Rev. B 73, 205119 (2006).
  * Román-Pérez and Soler (2009) Guillermo Román-Pérez and José M Soler, “Efficient implementation of a van der waals density functional: application to double-wall carbon nanotubes,” Phys. Rev. Lett. 103, 096102 (2009).
  * zen (2021) “Data repository for ‘approaching the basis set limit in gaussian-orbital-based periodic calculations with transferability: Performance of pure density functionals for simple semiconductors’,”  (2021).
  * Monkhorst and Pack (1976) Hendrik J Monkhorst and James D Pack, “Special points for brillouin-zone integrations,” Phys. Rev. B 13, 5188 (1976).
  * Sun  _et al._ (2020) Qiming Sun, Xing Zhang, Samragni Banerjee, Peng Bao, Marc Barbry, Nick S Blunt, Nikolay A Bogdanov, George H Booth, Jia Chen, Zhi-Hao Cui, _et al._ , “Recent developments in the pyscf program package,” J. Chem. Phys. 153, 024109 (2020).
  * Mardirossian and Head-Gordon (2013) Narbe Mardirossian and Martin Head-Gordon, “Characterizing and understanding the remarkably slow basis set convergence of several minnesota density functionals for intermolecular interaction energies,” J. Chem. Theory Comput. 9, 4453–4461 (2013).
  * Civalleri  _et al._ (2012) Bartolomeo Civalleri, Davide Presti, Roberto Dovesi, and Andreas Savin, “On choosing the best density functional approximation,” in  _Chemical Modelling: Applications and Theory Volume 9_, Vol. 9 (The Royal Society of Chemistry, 2012) pp. 168–185.
  * Ramberger  _et al._ (2019) Benjamin Ramberger, Zoran Sukurma, Tobias Schäfer, and Georg Kresse, “Rpa natural orbitals and their application to post-hartree-fock electronic structure methods,” J. Chem. Phys. 151, 214106 (2019).
  * Irmler  _et al._ (2018) Andreas Irmler, Asbjörn Burow, and Fabian Pauly, “Robust Periodic Fock Exchange with Atom-Centered Gaussian Basis Sets,” J. Chem. Theory Comput. 14, 4567–4580 (2018).
  * Gillan  _et al._ (2008) M. J. Gillan, D. Alfè, S. de Gironcoli, and F. R. Manby, “High-precision calculation of Hartree-Fock energy of crystals,” J. Comput. Chem. 29, 2098–2106 (2008).
  * 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,” J. Chem. Phys. 130, 184103 (2009).
  * Paier  _et al._ (2009) Joachim Paier, Cristian V. Diaconu, Gustavo E. Scuseria, Manuel Guidon, Joost VandeVondele, and Jürg Hutter, “Accurate Hartree-Fock energy of extended systems using large Gaussian basis sets,” Phys. Rev. B 80, 174114 (2009).
  * Civalleri  _et al._ (2010) Bartolomeo Civalleri, Roberto Orlando, Claudio M. Zicovich-Wilson, Carla Roetti, Victor R. Saunders, Cesare Pisani, and Roberto Dovesi, “Comment on “Accurate Hartree-Fock energy of extended systems using large Gaussian basis sets”,” Phys. Rev. B 81, 106101 (2010).
  * Usvyat  _et al._ (2011) Denis Usvyat, Bartolomeo Civalleri, Lorenzo Maschio, Roberto Dovesi, Cesare Pisani, and Martin Schütz, “Approaching the theoretical limit in periodic local MP2 calculations with atomic-orbital basis sets: The case of LiH,” J. Chem. Phys. 134, 214105 (2011).