Basis set construction for molecular electronic structure theory: Natural orbital and Gauss-Slater basis for smooth pseudpotentials

0.49 MB·English

by  F. R. Petruzielo

Basis set construction for molecular electronic structure theory:
Natural orbital and Gauss-Slater basis for smooth pseudopotentials
F. R. Petruzielo1,∗ Julien Toulouse2,† and C. J. Umrigar1‡
1Laboratory of Atomic and Solid State Physics,
Cornell University, Ithaca, New York 14853, USA
2Laboratoire de Chimie Th´eorique,
Universit´e Pierre et Marie Curie and CNRS,
75005 Paris, France
1
(Dated: January 28, 2011)
1
0 A simple yet general method for constructing basis sets for molecular electronic structure calcu-
2 lationsispresented. Thesebasissetsconsistofatomicnaturalorbitalsfrom amulti-configurational

n self-consistentfieldcalculation supplementedwithprimitivefunctions,chosensuchthattheasymp-
a toticsareappropriateforthepotentialofthesystem. Primitivesareoptimizedforthehomonuclear
J diatomic molecule to produce a balanced basis set. Two general features that facilitate this basis
7 constructionaredemonstrated. First,weakcouplingexistsbetweentheoptimalexponentsofprim-
2 itiveswithdifferentangularmomenta. Second,theoptimalprimitiveexponentsforachosensystem
depend weakly on the particular level of theory employed for optimization. The explicit case con-

] sidered here is a basis set appropriate for the Burkatzki-Filippi-Dolg pseudopotentials. Since these
ci pseudopotentials are finite at nuclei and have a Coulomb tail, the recently proposed Gauss-Slater
s functions are the appropriate primitives. Double- and triple-zeta bases are developed for elements
- hydrogenthrough argon. These newbases offer significant gains overthecorresponding Burkatzki-
l
r Filippi-Dolg bases at various levels of theory. Using a Gaussian expansion of the basis functions,
t
m these bases can be employed in any electronic structure method. Quantum Monte Carlo provides
an added benefit: expansions are unnecessary since theintegrals are evaluated numerically.
.
t
a
m I. INTRODUCTION exactly reproduce the correct electron-nucleus cusp and
- long-range asymptotic behavior of the orbitals. For cal-
d
n In quantum chemistry (QC) calculations, molecular culations on systems with a potential that is finite at
o orbitals are traditionally expanded in a combination of the nucleus and has a Coulomb tail, Gauss-Slater (GS)
c primitive Gaussian basis functions and linear combina- primitives[5]arethe appropriatechoicesincetheyintro-
[ tionsofGaussianprimitivescalledcontractedbasisfunc- ducenocuspattheoriginandreproducetheexponential

2 tions [1]. These basis sets cannot express the correct long-range asymptotic behavior of the orbitals.
v molecular orbital asymptotic behavior but are used in Despite shortcomings, traditional QC basis sets have
8 QCcalculationstopermitanalyticevaluationofthetwo- yielded good results. The natural orbitals (NOs) from a
1 electron integrals [2]. post Hartree-Fock (HF) method are a particularly suc-
3 Analyticintegralevaluationsignificantlylimitsflexibil- cessful form of contracted function [6–9]. The simplest
3 ity in basis set choice but is essential for computational NO construction involves diagonalizing the one-particle
.
5 efficiencyinQCcalculations. However,inpractice,other density matrix from a ground state atomic calculation
0 basisfunctionformscanbeconsideredsinceanarbitrary [6]. This construction is unbalanced due to obvious bias
0 function can be expanded in Gaussians. Of course, the favoring the atom. More complicated constructions in-
1
fidelity ofthis representationis limited. An expansionin volve diagonalizing the average one-particle density ma-
:
v afinite number ofGaussianscannotreproducethe expo- trixofseveralsystems: atomicgroundandexcitedstates,
Xi nential decay of the wavefunction at large distances or ions, diatomic molecules, and atoms in an external elec-
the Kato cusp conditions [3] at nuclei, but it can mimic tric field [7–9]. These constructions produce excellent
r
a these features over a finite range. results, but they are complex.
Quantum Monte Carlo (QMC) calculations [4] offer Asimplebutgeneralmethodforconstructingbasissets
greater freedom in choice of basis functions because ma- formolecularelectronicstructurecalculationsisproposed
trix elements are evaluated using Monte Carlo integra- and tested here. The bases are combinations of the NOs
tion. Consequently, the correct short- and long-distance obtainedfromdiagonalizingthe one-particledensity ma-
asymptotics can be satisfied exactly. For systems with trix from an atomic multiconfigurational self-consistent
a divergent nuclear potential, Slater basis functions can field (MCSCF) calculation and primitive functions ap-
propriate for the potential in the system. The primi-
tives are optimized for the homonuclear dimer in cou-
pled cluster calculations with single and double excita-
∗ [email protected] tions(CCSD),withtheintentionofproducingabalanced
† [email protected] basis set. Importantly, optimal exponents for the prim-
‡ [email protected] itive functions are shown to depend weakly on the level
2
of theory used in the optimization. Additionally, results for the BFD pseudopotentials, a single S orbital is the
show that coupling is weak between primitive functions ground state wavefunction, and this can be obtained ex-
of different angular momenta. This enables efficient de- actly in HF. Thus, the basis for each element in Group
termination of optimal exponents. 1AincludesasingleScontraction,noPcontractions,and
The utility of the above construction is demonstrated an appropriate number of GS primitives.
for the elements hydrogen through argon with the non-
divergent pseudopotentials of Burkatzki, Filippi, and
Dolg (BFD) [10]. Since these pseudopotentials are fi- A. Contracted Functions
nite at the nuclei and have a Coulomb tail, the GS
functions are the appropriate primitives. These pseu-
A contracted basis function is a linear combination of
dopotentials are chosenfor demonstrated accuracyin all
Gaussian primitives:
cases tested and because they are accompanied by a ba-
sis set. The BFD basis [10] serves as a metric for test-
ing the new basis. The benefits of our bases extend to ϕ (r,θ,φ)= c 2(2αi)n+12 rn−1e−αir2Zm(θ,φ),
allelectronicstructuremethods tested,including CCSD, nlm i s Γ(n+ 1) l
i 2
HF,theBeckethree-parameterhybriddensityfunctional X
(1)
(B3LYP) [11], and QMC.
where r,θ,φ are the standard spherical coordinates, n
ThemainareaofinterestfortheauthorsisQMC.Since is the principal quantum number, l is the azimuthal
QMCresultsdependlessonbasissetthantraditionalQC quantum number, m is the magnetic quantum number,
methods [5], only double-zeta (2z) and triple-zeta (3z) Zm(θ,φ) isarealsphericalharmonic,c istheith expan-
l i
bases are presented. sion coefficient, and α is the ith Gaussian exponent. In
i
Thispaperisorganizedasfollows. Basisfunctionform practice, the restriction n=l+1 applies.
and properties are demonstrated in Sec. II. Results for
The exponents of the primitive functions that form
calculationswiththe newbasesarediscussedinSec. III.
the contractedbasis functions are determined as follows.
Concluding remarks are provided in Sec. IV. Supple-
For each angular momentum for which a contraction is
mentary material is provided on EPAPS [12].
desired, an uncontracted basis consisting of nine even-
tempered primitive Gaussians is generated. For each
set of uncontracted Gaussians, the minimum exponent
II. BASIS SET
and even-tempering coefficient are varied to minimize
the CCSD energy of the atom using a Python wrapper
The number of basis functions for each angular mo- around GAMESS [13].
mentumfollowsthecorrelationconsistentpolarizedbasis An assumptionof weakcoupling between the different
set prescription of Dunning [1]. 2z and 3z bases appro- angular momenta underlies the optimization procedure.
priate for the BFD pseudopotentials are generated for Consequently, the uncontracted basis for each angular
the elements hydrogen through argon. Since the BFD momentum is optimized separately. This optimization
pseudopotential removes no core for hydrogen and he- is performed by calculating the CCSD energy on an ini-
lium, the 2z basis for these elements consists of two S tially coarse grid composed of different minimum expo-
functions and one P function, while the 3z basis consists nents and even-tempering coefficients. Once regions of
of three S functions, two P functions, and one D func- lowCCSDenergyareidentified,afinergridisusedtoob-
tion. Since the BFD pseudopotential removes a helium tainthefinalminimumexponentandeven-temperingco-
core for the first row atoms and a neon core for the sec- efficient. Inadditiontotheassumptionofweakcoupling,
ond row atoms, the remaining elements lithium through twoother propertiesofthe problemmakethis globalop-
argon have the same number of basis functions. In par- timizationpossiblewith modestcomputerresources;low
ticular, the 2z basis consists of two S functions, two P dimensionality of search space and efficiency of atomic
functions,andoneDfunction,whilethe3z basisconsists CCSD calculations.
of three S functions, three P functions, two D functions, Next, anatomic MCSCFcalculationin acomplete ac-
and one F function. tive space (CAS) with the optimized uncontracted basis
The bases consist of a combination of contracted and isperformedinGAMESS.Forthesecalculations,allelec-
primitive functions. Since the BFD pseudopotentialsare trons not removedby the pseudopotentialare allowedto
finiteattheoriginandhaveaCoulombtail,theGSfunc- excite. For helium, the active space consists of the or-
tions arethe appropriateprimitives. With the exception bitals from the n = 1 and n = 2 shells. For beryllium
oftheelementsinGroup1Aoftheperiodictable(i.e. H, throughneon,the activespaceincludesthe orbitalsfrom
Li, and Na), the basis for each element includes a single the n = 2 and n = 3 shells. For magnesium through
S contraction and a single P contraction combined with argon, the active space is composed of the orbitals from
an appropriate number of GS primitives. Only two con- the n=3 and n=4 shells, with the exception of the 4D
tractions are employedto reduce the computational cost and4Forbitals. Asubsetofthenaturalorbitalsfromthe
of using this basis in QC calculations. Since elements in MCSCFcalculationsareusedasthecontractedfunctions
Group 1A of the periodic table have only one electron of our basis.
3
All atomic calculations are performed in D symme- expansion is
2h
try since GAMESS does not permit imposition of full
rsaomtaetioantoamlsiycmsumbesthreyl.laHreenncoe,tdnieffceersesnatriclyomeqpuoinveanletnsto.fAthde- ϕζ (r,θ,φ)= cζ 2(2αζi)l+23 rle−αζir2Zm(θ,φ),
nlm i s Γ(l+ 3) l
ditionally, mixing may occur among orbitals of different i 2
X
angular momenta. For instance, there is mixing of S (6)
orbitals with both D and D orbitals. This where cζ is the ith expansion coefficient and αζ is the
3z2−r2 x2−y2 i i
anisotropy can be removed by averaging the different ith Gaussian exponent. Notice that the expansion per-
components of a particular subshell and zeroing out the mits the case for which n 6= l+1 for the GS function.
off-diagonalblocks ofthe one-particledensity matrix [7]. Additionally, the following scaling relations hold for the
Asimplerapproachtakeninthis workisfoundto pro- expansion coefficients and Gaussian exponents:
duce results of similar quality. For each angular mo-
mentum for which a contraction is desired, the NO with αζ =ζ2α1 (7)
i i
that angular momentum which has the largest occupa- cζ =c1. (8)
tion number is chosen. Additionally, NO elements which i i
do not correspond to the dominant character of the or-
Once the Gaussian expansions are found for unit expo-
bital are zeroed out. For instance, an NO with large co-
nents, expansions of arbitrary GSs follow immediately
efficients on the S basis functions and small coefficients
from the scaling relations. For QC calculations in this
on the D basis functions is considered to be dominated
paper, GSs are expanded in six Gaussians. However, if
by S character, so the D coefficients are zeroed out. Fi-
the purpose of the initial QC calculation is to generate
nally, the NOs are normalized. The NOs selected in this
crude starting orbitals for QMC calculations in which
proceduregeneratethecontractedfunctionsforthebasis
orbital optimization is performed, it is only necessary to
set. The expansions of the contractions are given in the
expand GS primitives in a single Gaussian. In this case,
supplementary material [12].
the costof QCcalculations is the samefor Gaussianand
GSprimitives. TheexpansionsofGSfunctionswithunit
exponent in both one and six Gaussians are given in the
B. Gauss-Slater Primitives supplementary material [12].
As mentioned above, the restriction n ≥ l+1 is im-
GS functions [5] are defined as posed for GS functions, instead of the more familiar
n = l + 1 restriction imposed for Gaussian primitives.
ϕζ (r,θ,φ)=Nζ rn−1e−(1ζ+rζ)2r Zm(θ,φ), (2) Thismotivatesconstructionoftwotypesofbases. Inthe
nlm n l first,ANO-GS,therestrictionn=l+1isenforced. Inthe
second,ANO-GSn, for eachl there canbe atmosta sin-
where ζ is the GS exponent and Nζ is the normaliza-
n gleGSprimitivewithaparticularn. Foreachadditional
tion factor. The restriction n ≥ l+1 is imposed for GS
primitive with a particular l, n must be incremented.
functions. For r ≪1, the GS behaves like a Gaussian:
For example, consider lithium. The 2z ANO-GS basis
has one S contraction, one GS-1S function, two GS-2P
ϕζnlm(r,θ,φ)∼=Nnζ rn−1e−(ζr)2 Zlm(θ,φ), (3) functions, and one GS-3D function. On the other hand,
the2z ANO-GSnbasishasoneScontraction,oneGS-1S
and for r≫1, the GS behaves like a Slater: function, one GS-2P function, one GS-3P function, and
one GS-3D function.
ϕζnlm(r,θ,φ)∼=Nnζ rn−1e−ζr Zlm(θ,φ). (4) A caveat to the above definition of the ANO-GSn ba-
sis is that GS-2S functions are not permitted since a
Consequently,GSfunctionsintroducenocuspatthe ori- single GS-2S function will introduce an undesired cusp
ginandcanreproducecorrectlong-rangeasymptoticbe- in the wavefunction. Additionally, the 2z ANO-GS and
havior of the orbitals. ANO-GSn basis sets are identicalfor allelements except
UnlikeGaussiansandSlaters,normalizationofGSshas lithium and sodium. When the 2z ANO-GS and ANO-
no closed form expression. Nevertheless, normalizing an GSn basis sets are identical, the basis sets are referred
arbitraryGS is trivial with the following scaling relation to as a 2z ANO-GS/GSn basis. For both lithium and
between Nnζ and Nn1: sodium, the basis sets differ because these systems have
no P contractions and instead have a second P primi-
Nζ =ζn+1/2 N1. (5) tive for the 2z basis. This primitive is a GS-2P for the
n n
ANO-GSbasisandaGS-3PfortheANO-GSnbasis. Ad-
Values for N1 are given in the supplementary material ditionally, weak coupling between functions of different
n
[12]. angular momentum causes the GS-1S and GS-3D func-
Since GSs are not analytically integrable, the radial tions in the ANO-GS bases for lithium and sodium to
part must be expanded in Gaussians for use in QC pro- differ from their counterparts in the ANO-GSn bases.
grams that evaluate matrix elements analytically. The However,the optimal exponents differ by less than 0.01.
4
OptimalexponentselectionfortheGSprimitivesisdis- theory using a Python wrapper around GAMESS. For
cussednow. Insteadofoptimizingexponentsfortheatom each angular momentum, an energy landscape is defined
as was done to generate the contractions, optimization by a grid of primitive exponents ranging from 0.1 to 6.0
of the GS exponents is performed for the homonuclear with 0.1 spacing. Thorough investigation has revealed
diatomic molecule at experimental bond length [14–22]. that exponents larger than 6.0 are not optimal for the
This advantageously produces a balanced basis set. systems considered. Low lying minima of this energy
WeakcouplingbetweenGSfunctionsofdifferentangu- landscape are then handled with increasingly finer grids
lar momenta is assumed, so the initial optimization for until energy changes are less than 0.01 mH. During this
each angular momentum is performed separately. This investigation of local minima, all angular momenta are
assumptionisvalidatedinFigure1,whichcontainsplots handled simultaneously to account for any coupling ef-
of the CCSD energy for Si while varying individual GS fects. Results of this optimization are shown in Figure
2
exponentsinthe2z ANO-GS/GSnbasis. Boththecurve 2. Optimalexponents for ANO-GS andANO-GSn bases
shape and exponent value which minimizes the energy exhibit a linear trend across each row of the periodic ta-
varylittlewithfixedexponentvalue,signifyingweakcou- ble. For nearly degenerate minima, the exponent follow-
pling between GS functions of different angular momen- ing the trendin the figure is chosenasoptimal, resulting
tum. in energy increase no greater than several 0.1 mH. The
optimal GS exponents are given in the supplementary
-7.54 material [12].
GS-2P=0.8
In some cases, the optimal exponents for primitives
GS-2P=1.3
with the same n and l are very close. This can lead to
-7.55 GS-2P=2.8
large equal and opposite coefficients on these basis func-
tions when constructing molecular orbitals. Numerical
)
H -7.56 problems could result, providing further motivation for
(
y the ANO-GSn basis, in which each pair of n and l is
g
r unique. However,allofourtestswiththeANO-GSbasis
e
En -7.57 have had no numerical problems.
Finally, the optimal primitive exponents are found to
-7.58 depend weakly on the electronic structure method em-
ployed in the optimization, as demonstrated in Figure 3
for Si with the 2z ANO-GS/GSn basis. The globally
2
-7.59
minimizingexponentsarenearlyequalindifferentmeth-
0 0.5 1 1.5 2 2.5 3 3.5 4
ods. This exponent transferability to different levels of
GS-1S Exponent
theory is extremely attractive for a basis set.
-7.40
GS-1S=0.5
GS-1S=1.5
-7.44 GS-1S=3.5 III. RESULTS
)
H -7.48 Section II demonstrates that the ANO-GS and ANO-
(
y GSn bases exhibit desirable properties. However, it re-
g
er mains to be shown that these basis sets produce accu-
n -7.52
E rateresults. Fortunately,thebasissetaccompanyingthe
BFDpseudopotentialservesasametricfortestingANO-
-7.56 GS and ANO-GSn basis quality. The BFD basis for el-
ements in Groups 1A and 2A of the periodic table has
recently been updated [23], but the number of functions
-7.60
in the new basis is inconsistent with the correlationcon-
0 0.5 1 1.5 2 2.5 3 3.5 4
sistentpolarizedbasisprescription[1]. Since comparison
GS-3D Exponent
would be difficult, their published functions are consid-
FIG. 1. Change in Si2 CCSD energy for 2z ANO-GS/GSn ered in this work.
basis shows weak coupling between GS functions of different Figure 4 shows the CCSD total energy gain per elec-
angular momenta. TOP: Energy versus GS-1S exponent for tron of the ANO-GS and ANO-GSn bases over the BFD
threevaluesoftheGS-2PexponentwiththeGS-3Dexponent bases [10] for atoms and homonuclear dimers of hydro-
fixed at its optimal value. Bottom: Energy versus GS-3D genthroughargon. Energygainsperelectrontendtoin-
exponent for three values of the GS-1S exponent with the creaseacrosseachrow of the periodic table. BothANO-
GS-2P exponentfixed at its optimal value.
GS and ANO-GSn bases yield energy gains for most
molecules and atoms. The energy gains per electron are
generally larger for molecules than for atoms, and larger
The optimization is performed at the CCSD level of fortheANO-GSnbasisthanfortheANO-GSbasis. The
5
6 0
CCSD
3z ANO-GSn GS-1S -2 HF
5 GGSS--32SP H) -4 B3LYP
m
GS-3P
(
4 GS-3D y -6
nt GS-4D rg
e GS-4F e -8
on 3 En
Exp 2 e in --1120
g
n
a -14
h
1 C
-16
0 -18
HHeLiBeBCNOFNeNaMgAlSiPSClAr 0 0.5 1 1.5 2 2.5 3 3.5 4
GS-1S Exponent
6 0 0
3z ANO-GS GS-1S
5 GS-1S ) -1
GS-2P H -5
m
GS-2P
(
4 GS-3D y -2
nt GS-3D rg -10
e GS-4F e
on 3 En -3
Exp 2 e in -15 -4
g
n
a
h -20
1 C -5
0 -25 -6
HHeLiBeBCNOFNeNaMgAlSiPSClAr 0 0.5 1 1.5 2 2.5 3 3.5 4
GS-2P Exponent
5 0 0
2z ANO-GS/GSn GS-1S
GS-2P -20
)
4 GS-3D H -10
m
( -40
y
ent 3 erg -60 -20
n n
o E
Exp 2 e in -80 -30
g
n -100
a
1 h -40
C -120
0 -140 -50
HHeLiBeBCNOFNeNaMgAlSiPSClAr 0 0.5 1 1.5 2 2.5 3 3.5 4
GS-3D Exponent
FIG.2. OptimalexponentsforANO-GSandANO-GSnbases FIG. 3. Change in Si2 energy for 2z ANO-GS/GSn basis
exhibit a linear trend across each row of the periodic table. shows optimal exponents depend weakly on electronic struc-
The 2z ANO-GS and ANO-GSn bases are identical for all turemethod(CCSD,HF,andB3LYP).Top: GS-1Sexponent
elements except lithium and sodium. The GS-1S and GS- isvariedwithGS-2PandGS-3Dexponentsfixedattheirop-
3D exponentsfor theseelementseach differbyless than 0.01 timal values. Middle: GS-2P exponent is varied with GS-1S
between 2z ANO-GS and ANO-GSn bases, so 2z ANO-GS andGS-3Dexponentsfixedattheiroptimalvalues. Thelarge
and ANO-GSnare shown together as 2z ANO-GS/GSn. Ex- increaseinenergyaroundanexponentof1.0occurssincethe
ponents for GS functions of P angular momentum are not Pprimitive and Pcontraction become nearly linearly depen-
includedforlithiumandsodiumsincetheseelementshavean dent. Bottom: GS-3DexponentisvariedwithGS-1SandGS-
extra primitiveof P angular momentum. 2P exponents fixed at their optimal values. For Middle and
Bottom, HF and B3LYP energy scale is on the right y-axis.
This difference in energy scale occurs since higher angular
momentum functions are less important in these effectively
single-determinant theories.
6
energygainsforthe2z basesaregenerallylargerthanfor ANO-GSnbasesisdramatic: the2zANO-GS/GSnresult
the3z bases,asexpected,sincetheenergylefttorecover isnearlyhalfwaybetweenthe2zand3zBFDresults,and
becomes smaller as the basis size increases. the3zANO-GS/GSnresultisnearlyhalfwaybetweenthe
3z and5z BFDresults. ANO-GSandANO-GSnbenefits
aremoreprominentinHFandB3LYP:formostsystems,
14
H) 13 Atom 2z ANO-GS/GSn the2zANO-GS/GSnresultisclosertothe3zBFDresult
m 12 Dimer 2z ANO-GS/GSn thanthe2z BFDresult,andthe3z ANO-GS/GSnresult
n ( 11 Atom 3z ANO-GS is closer to the 5z BFD result than the 3z BFD result.
o Atom 3z ANO-GSn
ctr 10 Dimer 3z ANO-GS DifferencesbetweenresultswiththeANO-GSandANO-
e 9 Dimer 3z ANO-GSn GSn bases are small.
El
r 8 Figure 7 shows the fraction of experimental atomiza-
e
p tionenergyrecoveredusingdiffusionMonteCarlo(DMC)
n 4 with the BFD, ANO-GS, andANO-GSnbases. For each
Gai 3 system, the DMC calculations are performed with both
y 2 a single-configuration state function (single-CSF) refer-
g
er 1 ence(DMC-1CSF)andfull-valencecompleteactivespace
n 0
E reference (DMC-FVCAS). However, for each of the con-
-1
HHeLiBeBCNOFNeNaMgAlSiPSClAr sCtSitFuernetfearteonmcessinarteheeqseuimvaolleenctu.leAs,llthDeMFCVCcaAlcSualantdiosninsgalree-
performedwitha0.01H−1 time stepandtrialwavefunc-
FIG. 4. CCSD total energy gains per electron of ANO-GS tionobtainedbyoptimizing Jastrow,orbital,andconfig-
and ANO-GSn relative to the corresponding BFD basis [10] urationstatefunction(CSF) parameters(whereapplica-
for atoms and homonuclear dimers of hydrogen through ar- ble) via the linear method [28–30] in variational Monte
gon. Energy gains per electron tend to increase across each Carlo. The DMC-1CSF and DMC-FVCAS calculations
row of the periodic table. The 2z ANO-GS and ANO-GSn exhibit similar trends to the HF and B3LYP calculation
basesareidenticalforallelementsexceptlithiumandsodium.
for most systems: the 2z ANO-GS/GSn result is closer
Differences between 2z ANO-GS and ANO-GSn results for
tothe3z BFDresultthanthe2z BFDresult,andthe3z
theseelementsis∼0.01mH,sotheyareshowntogetheras2z
ANO-GS/GSn resultis closer to the 5z BFD resultthan
ANO-GS/GSn.
the 3z BFD result. Again, differences between results
with the ANO-GS and ANO-GSn bases are small.
The ANO-GS and ANO-GSn bases also produce more There are several important points that can be made
accurateCCSDatomizationenergiesthantheBFDbasis by comparing the DMC calculations of Figure 7 to the
for the homonuclear dimers of hydrogen through argon. CCSD calculations of Figure 6. First, the DMC results
Figure 5 shows the fraction of experimental atomization for the atomization energies have a weaker dependence
energy recovered in CCSD for the homonuclear dimers on basis size than the CCSD results. Second, for a given
which are not weakly bound. The 2z ANO-GS/ANO- basisset,themostbasicDMCcalculations,DMC-1CSF,
GSn basis recoversmore atomizationenergythan the 2z yield superior results comparedto CCSD. In addition to
BFD basis for all dimers except those of Group 1A ele- yielding superior results, DMC-1CSF calculations have
ments. Similarly, the 3z ANO-GSn basis recovers more better computational cost scaling than CCSD calcula-
atomization energy than the 3z BFD basis for the same tions. Undercertainassumptions,thecostofDMC-1CSF
systems,butthedifferencesaresmall. The3z ANO-GSn calculationsscalesasO(N3)[31],whilethecostofCCSD
is on average slightly better than the 3z ANO-GS basis, calculations scales as O(N6) [32], where N is the num-
the largest gains being for F and Cl . ber of electrons. However, it is important to note that
2 2
For Group 1A elements, the BFD bases recover more theprefactorofthescalingissignificantlysmallerforthe
atomization energy in CCSD than do their ANO-GS or CCSD calculations.
ANO-GSn counterparts. This occurs due to inaccurate Finally,ourresultsarenotthefirsttoshowthatDMC
BFD energies for the atoms, as can be seen in Figure calculations can produce accurate atomization energies.
4. However, as described above, we used the published In particular, DMC-1CSF calculations of the entire G2
BFD bases for these elements rather than the updated set have been performed for both pseudopotential and
BFD bases [23] to maintain consistency. all-electron systems [33, 34] and produced excellent re-
Finally, improvements of the ANO-GS and ANO-GSn sults. Additionally, there is goodagreementbetweenthe
bases extend to other systems and methods. Figure 6 pseudopotential and all-electronresults with a mean ab-
shows the fraction of experimental atomization energy solutedeviationofabout2.0kcal/molovertheentireG2
recovered for five systems in the G2 set [26] with the set [34]. Although these previous results are very good,
BFD, ANO-GS, and ANO-GSn bases in three quantum there is room for improvement, particularly for the open
chemistry methods. For CCSD, the ANO-GS and ANO- shell systems. A systematic study with DMC-FVCAS
GSnbasesoutperformtheBFDbasisforallsystems. For calculations is currently underway in our group, which
sulfur dioxide the improvementdue to the ANO-GS and should produce results to (near) chemical accuracy for
7
y 1.2 2z BFD 3z ANO-GS
g
r 2z ANO-GS/GSn 3z ANO-GSn
e
n 1.1 3z BFD 5z BFD
E
g 1
n
di
n 0.9
Bi
p. 0.8
x
E
0.7
of
n 0.6
o
cti 0.5
a
r
F
0.4
H Li B C N O F Na Al Si P S Cl
2 2 2 2 2 2 2 2 2 2 2 2 2
FIG. 5. Fraction of experimental atomization energy recovered in CCSD with BFD, ANO-GS, and ANO-GSn bases for the
homonuclear dimers of hydrogen through argon which are not weakly bound. The 2z ANO-GS and ANO-GSn bases are
identicalfor all elementsexceptlithium andsodium. Differencesbetween 2z ANO-GSandANO-GSnatomization energies for
these elements is ∼ 0.01mH, so they are shown together as 2z ANO-GS/GSn. Calculated values are corrected for zero point
energy [17, 24] to compare with experiment [14, 16,17, 25].
all systems in the G2 set. the ANO-GSn basis is slightly better than the ANO-GS
basis, but either is a sound choice.
In the future, these basis sets will be extended to
IV. CONCLUSION include the transition metals, and, bases will be con-
structed for all-electron calculations, for which Slater
functions are the appropriate primitives.
Asimpleyetgeneralmethodforconstructingbasissets
formolecularelectronicstructuretheorycalculationshas
beenpresented. Thesebasissetsconsistofacombination V. ACKNOWLEDGMENTS
of atomic natural orbitals from an MCSCF calculation
withprimitivefunctions optimizedforthecorresponding
WethankClaudiaFilippiforveryvaluablediscussions.
homonucleardimer. Thefunctionalformoftheprimitive
ThisworkwassupportedbytheNSF(GrantNos. DMR-
functionsischosentohavethecorrectasymptoticsforthe
0908653 and CHE-1004603). Computations were per-
nuclear potential of the system.
formedinpartattheComputationCenterforNanotech-
Itwasshownthatoptimalexponentsofprimitiveswith
nology Innovation at Rensselaer Polytechnic Institute.
different angular momenta are weakly coupled. This en-
ables efficient determination of optimal exponents. Ad-
ditionally, it was demonstrated that the particular elec-
tronic structure method employed in optimization has
little effect on the optimal values of the primitive expo-
nents.
Twosets of2z and3z bases,ANO-GSandANO-GSn,
appropriate for the Burkatzki, Filippi, and Dolg non-
divergentpseudopotentialswereconstructedforelements
hydrogen through argon. Since these pseudopotentials
do not diverge at nuclei and have a Coulomb tail, GS
functions are the appropriate primitives.
It was demonstrated that both ANO-GS and ANO-
GSn basis sets offer significant gains over the Burkatzki,
Filippi and Dolg basis sets for CCSD, HF, B3LYP [11],
and QMC calculations. Improvements were observed in
both total energies and atomization energies. The latter
indicatesthatbasissetsprovidingabalanceddescription
ofatomsandmoleculeswereproducedbyusingboththe
atom and the dimer in the optimization. On average,
8
1.10 1.10
y y
g CCSD g DMC-FVCAS
er er
En 1.00 En 1.00
g g
n n
ndi 0.90 ndi 0.90
Bi Bi
p. p.
Ex 0.80 Ex 0.80
of of
n n
ctio 0.70 ctio 0.70
a a
Fr Fr
0.60 0.60
LiF O2 P2 S2 SO2 LiF O2 P2 S2 SO2
1.10 1.10
y y
g B3LYP g DMC-1CSF
er er
En 1.00 En 1.00
g g
n n
ndi 0.90 ndi 0.90
Bi Bi
p. p.
Ex 0.80 Ex 0.80
of of
n n
ctio 0.70 ctio 0.70
a a
Fr Fr
0.60 0.60
LiF O2 P2 S2 SO2 LiF O2 P2 S2 SO2
1.00
y 2z BFD FIG. 7. Fraction of experimental atomization energy recov-
rg HF 2z ANO-GS eredindiffusionMonteCarlo(DMC)forLiF,O2,P2,S2,and
e
n 0.80 2z ANO-GSn SO2 with the BFD, ANO-GS, and ANO-GSn bases. DMC
E 3z BFD calculations are performed with both a single-CSF reference
ng 3z ANO-GS (DMC-1CSF)andfull-valencecompleteactivespacereference
ndi 0.60 3z ANO-GSn (DMC-FVCAS).The2z ANO-GSandANO-GSnbasesyield
Bi 5z BFD* different results only for LiF. The 5z BFD* calculations do
p. not include the G or H functions from the 5z BFD basis.
Ex 0.40 Calculated atomization energies are corrected for zero point
of energy [17, 24] to compare with experiment [14, 16, 17, 27].
n The legend for this plot is identical to that of Figure 6.
o 0.20
cti
a
r
F
0.00
LiF O P S SO
2 2 2 2
FIG. 6. Fraction of experimental atomization energy recov-
ered inHF,B3LYP,and CCSD forLiF, O2,P2,S2,andSO2
with BFD, ANO-GS, and ANO-GSn bases. The 2z ANO-
GS and ANO-GSn bases yield different results only for LiF.
The 5z BFD* calculations do not include the G or H func-
tionsfromthe5zBFDbasis. Calculatedatomizationenergies
are corrected for zero point energy [17, 24] to compare with
experiment [14, 16, 17, 27].
9
[1] T. DunningJr, J. Chem. Phys.90, 1007 (1989). Number 101, edited by R. D. Johnson (NIST, Gaithers-
[2] S.Boys, Proc. R.Soc. A 200, 542554 (1950). burg, MD, 2010).
[3] T. Kato, Comm. Pure Appl.Math 10, 151177 (1957). [18] R. Grisenti, W. Schollkopf, J. Toennies, G. Hegerfeldt,
[4] W.Foulkes,L.Mitas, R.Needs,and G.Rajagopal, Rev. T.Kohler,andM.Stoll,Phys.Rev.Lett.85,2284(2000).
Mod. Phys. 73, 3383 (2001). [19] V. E. Bondybey and J. H. English, J. Chem. Phys. 80,
[5] F. R. Petruzielo, J. Toulouse, and C. J. Umrigar, 568 (1984).
J. Chem. Phys. 132, 094109 (2010). [20] R. Aziz and M. Slaman, Chem. Phys.130, 187 (1989).
[6] J.Alml¨ofandP.Taylor,J.Chem.Phys.86,4070(1987). [21] Z. Fu, G. W. Lemire, G. A. Bishea, and M. D. Morse,
[7] P.-O. Widmark, P.-A. Malmqvist, and B. Roos, Theor. J. Chem. Phys. 93, 8420 (1990).
Chim. Acta. 77, 291 (1990). [22] P. R. Herman, P. E. Larocque, and B. P. Stoicheff, J.
[8] P.-O. Widmark, B. Joakim, and B. Roos, Theor. Chim. Chem. Phys.89, 4535 (1988).
Acta. 79, 419 (1991). [23] http://www.burkatzki.com/pseudos/index.2.html.
[9] V. Veryazov, P.-O. Widmark, and B. O. Roos, [24] K. Irikura, J. Phys. Chem. Ref. Data 36, 389 (2007).
Theor. Chim. Acta. 111, 345 (2004). [25] Y. R. Luo, Comprehensive Handbook of Chemical Bond
[10] M. Burkatzki, C. Filippi, and M. Dolg, J. Chem. Phys. Energies (CRC Press, 2007).
126, 234105 (2007). [26] L.A.Curtiss,K.Raghavachari,G.W.Trucks,andJ.A.
[11] A.Becke, J. Chem. Phys.98, 5648 (1993). Pople, J. Chem. Phys.94, 7221 (1991).
[12] See Supplementary Material at [27] D. Feller and K. A. Peterson, J. Chem. Phys. 110, 8384
http://dx.doi.org/10.1063/1.3551512 for atom spe- (1999).
cific basis sets, Gaussian fits of Gauss-Slater functions, [28] J.ToulouseandC.Umrigar,J.Chem.Phys.126,084102
and data used for making figures. (2007).
[13] M. W. Schmidt, J. A. Boatz, K. K. Baldridge, S. T. El- [29] J.ToulouseandC.Umrigar,J. Chem. Phys.128,174101
bert, M. S. Gordon, J. H. Jensen, S. Koseki, N. Mat- (2008).
sunaga, K. A.Nguyen,S.Su,T. L.Windus, M. Dupuis, [30] C. Umrigar, J. Toulouse, C. Filippi, S. Sorella, and
andJ.A.Montgomery,J. Comp. Chem.14,1347(1993). R. Hennig, Phys.Rev.Lett. 98, 110201 (2007).
[14] K. P. Huber and G. Herzberg, Constants of Diatomic [31] R. J. Needs, M. D. Towler, N. D. Drummond, and
Molecules, Molecular Spectra and Molecular Structure P. L´opez R´ıos, J. Phys.Condens. Matter 22, 023201
Vol 4 (Van Nostrand Reinhold Company, 1979). (2010).
[15] L. Gurvich, I. V. Veyts, and C. B. Alcock, Thermody- [32] T. Helgaker, P. Jorgensen, and J. Olsen, Molecular
namicPropertiesofIndividualSubstances, FouthEdition Electronic-Structure Theory (John Wiley & Sons LTD,
(HemispherePub. Co., 1989). Chichester, England, 2000).
[16] NIST Chemistry WebBook, NIST Standard Reference [33] J. C. Grossman, J. Chem. Phys. 117, 1434 (2002).
DatabaseNumber69,editedbyP.J.LinstromandW.G. [34] N. Nemec, M. D. Towler, and R. J. Needs,
Mallard (NIST,Gaithersburg, MD, 2005). J. Chem. Phys. 132, 034111 (2010).
[17] NISTComputationalChemistryComparisonandBench-
mark Database, NIST Standard Reference Database