We construct a new modification of correlation consistent effective core potentials (ccECPs) for late 3d elements Cr–Zn with Ne-core that are adapted for efficiency and low energy cut-offs in plane wave calculations. The decrease in accuracy is rather minor, so that the constructions are in the same overall accuracy class as the original ccECPs. The resulting new constructions work with energy cut-offs at or below ≈400 Ry and, thus, make calculations of large systems with transition metals feasible for plane wave codes. We also provide the basic benchmarks for atomic spectra and molecular tests of this modified option that we denote as ccECP-soft.

Key properties of matter, such as cohesion and magnetic or optical responses, can be derived from valence electronic structure calculations. Fortunately, electronic levels in atoms show a significant distinction between core and valence states so that it is possible to introduce valence-only effective Hamiltonians. The fact that core and valence states occupy different ranges both in spatial and in energy domains enables us to partition the atomic states into core and valence subspaces. The theory of these well-known pseudopotentials or effective core potentials (ECPs) has been perfected over the decades, and it includes a number of criteria that model the influence of the core on valence electrons as closely as possible to the original, all-electron atom. This includes concepts such as norm-conservation of one-particle states, consistency of energy differences for atomic excitations, different forms, and beyond.1–16 The advantages of using ECPs are both qualitative and quantitative. On the quantitative side, the valence energy scale and number of degrees of freedom are essentially unchanged across the periodic table with resulting orders of magnitude speed-ups in calculations. On the qualitative front, ECPs can be constructed to mimic true atoms by effectively taking into account not only single-particle pictures but also core–core correlations, impact of core–valence correlations, spin–orbit and other relativistic effects by employing surprisingly simple, nonlocal, one-particle operators with lm- or jlm-projectors. Of course, there is a price for these advantages since there is always some ECP-related bias present. However, over the years, the accuracy of high-quality ECPs has been steadily improving and, at present, the corresponding errors typically do not dominate electronic structure calculations.

Recently, we have introduced a new generation of such effective Hamiltonians, which we call correlation consistent ECPs (ccECPs).17 We have emphasized several key principles in order to significantly improve the fidelity with respect to the original, all-electron Hamiltonian in broad classes of calculations, such as molecular systems and condensed matter materials with periodic or mixed boundary conditions. These principles involve (i) constructions that minimize discrepancies between all-electron and ccECP many-body atomic spectra as well as one-particle properties, such as charge norm conservation; (ii) use of many-body methods that include Coupled Cluster in constructions as well as quantum Monte Carlo (QMC) in testing and benchmarking; (iii) for certain elements, we introduced several sizes of valence and core subspaces including all-electron regularized Coulomb potentials for very light elements (for example, reg-ccECPs for H–Be elements); (iv) easy use with simple parametrizations in Gaussian expansions; and (v) open data website with full access and further adapting to particular types of calculations.

The derived ccECPs indeed turned out to be, in general, more consistently accurate than previous constructions, and they also provide a better balance of accuracy in various settings. In particular, the tests on molecules in non-equilibrium geometries have demonstrated significantly improved transferability apart from equilibrium atomic conformations. Some deviations from chemical accuracy have occurred for early main group elements in 3s3p and 4s4p columns at very short bond lengths of oxide molecules. This is a well-known limitation due to the small number of valence electrons and polarizability of the most shallow core states. For 3s3p elements, this issue has been addressed by providing a [He]-core option that makes the calculations essentially equivalent to the all-electron setting. Being well-defined and tested, ccECPs also enabled us to put a bound on systematic biases in quantum Monte Carlo calculations, for example, in the study of molecular and solid state systems, including very large supercell sizes with hundreds of valence electrons.18,19 Therefore, ccECPs provide both a practical tool and a valuable accuracy standard for the benchmarking of other constructions, and as such, they are being independently probed by the community at large.20 

At the same time, the higher accuracy implies deeper potential functions that make the valence-only electronic states less smooth and more curved in the core region. This is fundamentally correct since it mimics more accurately the distribution of effective valence charges inside the core. This is corroborated by capturing the correct shape of molecular binding curves for hydride and oxide dimers that probe both covalent and polarized bonds. For most elements, the conventional core ccECPs can be used both in plane wave calculations with energy cut-offs below roughly 400 Ry with converged energy to 1 meV/electron or so. Unfortunately, for late 3d transition metals with very deep 3s, 3p semicore and localized 3d states, the cost of plane wave calculations goes up very significantly, needing a cut-off of around 1400 Ry or more.

In order to overcome this limitation, we have identified possible ways to decrease the energy cut-offs while keeping the accuracy at a level comparable to the original ccECPs. Note that there is often a very subtle balance between the parameterization details on the one hand and their impact on ECP properties on the other. As we commented upon before,17 even rather minor changes in ECP parameters could push the accuracy, charge smoothness, and curvature in directions that are counterproductive. Similarly, over constraining of cost (objective) functions and/or over parameterizations could be counterproductive as well since that could lead to linear dependencies, costly undue optimizations, and overall inefficiency. In what follows, we briefly describe the results of these adapted constructions, the corresponding forms, relevant updates of previously introduced methods, results, and testing. The resulting ccECP-soft constructions exhibit much shallower potential functions and offer major efficiency gains in plane wave calculations.

The constructed ccECP-soft should be useful for a broad set of calculations since for many valence properties, they exhibit a high degree of fidelity to the original all-electron Hamiltonians. Computationally, this is particularly relevant for highly developed plane wave codes that employ Density Functional Theory (DFT) and post-DFT methods. We note that DFT orbitals are now routinely used as inputs for trial wave functions in QMC calculations. On the other hand, ccECP-soft should also be attractive for DFT and post-DFT calculations due to the mentioned genuine increase in accuracy and transferability. Of course, it is understood that the DFT biases with ECPs vs with all-electron in the core might differ due to the major change in electron density inside the core. Therefore, one has to pay attention to the consistency of DFT error cancellations for energy differences and other properties. This is a known issue that has been studied previously, for example, in Refs. 21 and 22.

The parameterization of the ECPs is unchanged from23 following a semi-local format,

Vipp=Vloc(ri)+=0maxV(ri)m|mm|,
(1)

where ri indicates the radial contribution of the ith electron, and max is 1 for the elements investigated in this work. Vloc is chosen to cancel out the Coulomb singularity. Vloc(r) has the form

Vloc(r)=Zeffr(1eαr2)+αZeffreβr2+i=12γieδir2,
(2)

where Zeff is the effective charge of the valence space while α, β, δi, and γi are parameters determined by the optimization. This form explicitly cancels the Coulomb singularity and insures smooth behavior at very small radii.8 

The non-local potentials are parameterized as

Vl(r)=j=12βjrnj2eαjr2,
(3)

where βℓj and αℓj are optimized for each non-local channel.

The methods employed in this work are largely adapted and updated from our previous studies;23 therefore, here we recount only the basic points and current modifications. Most of the differences from the previous publication come in the form of more stringent constraints during the optimization. In the cited work, the exponents of the Gaussian expansions that form the pseudopotential were allowed to freely move during optimization; thus, they could increase or decrease as needed within a very wide range of permitted values. Since the original objective was to maximize accuracy and fidelity when compared to the all-electron setting, the resulting pseudopotential functions were in general varied so as to reflect the true all-electron and bare ion interactions experienced by a valence electron.

In order to construct softer ECPs, constraints were placed on the maximum allowed values for the exponents. This has two effects. First, it smooths out the curvature of the potentials and, second, it indirectly restricts the potentials’ amplitudes. The resulting potentials are therefore more shallow and exhibit lower curvatures. This in turn leads to lower cutoffs in plane wave expansions. Finally, in our previous research we included highly ionized states in the objective function to construct each ccECP, while in this work we only ionize to around +4 or +5. This focuses the objective function to only account for likely configurations that could feasibly occur in molecules and solids and makes the optimization easier to perform reliably. Thus, the optimization protocol is as follows:

  1. Calculate high accuracy coupled-cluster with singles, doubles, and perturbative triples [CCSD(T)] all-electron (AE) data for each element that involves a set of excited states within the desired energy window, and generate initial ECP candidates with confirmed cutoffs below the desired threshold.

  2. Applying the same techniques we previously employed,23 incorporate correlation energy into the optimization implicitly by finding the contributions to the energy from the scalar relativistic Dirac–Fock calculation and then correlation contributions from CCSD(T) calculations. Using the AE data, the gaps are shifted as ΔEshifted=ΔEAEΔEECPcorr, where ΔEAE is the AE coupled cluster (CC) gap, and ΔEECPcorr is the ECP’s correlation energy contribution for the same gap. This is viable because the correlation energy between ECPs that share the same valence space tends to remain largely constant, as shown in our previous paper.23 These new gaps will serve as the major component of the objective function.

  3. Considering a many-body spectrum S, the objective function Γ is given by
    Σ(S)=sSws(ΔEECP(s)ΔEshifted(s)),
    (4)
    Γ=Σ(S)+γ(ϵECPϵAE)2,
    (5)
    and it is minimized from an initial guess using the DONLP2 routine.24ϵAE denote the AE one-particle scalar relativistic Hartree–Fock eigenvalues for the semi-core orbitals (3s and 3p in this case) and similarly ϵECP denote the corresponding Hartree–Fock eigenvalues.
  4. Step 3 will keep iterating to minimize the objective function, resulting in the final ECP once a sufficient minimum is reached.

Equation (5) is the same as seen in previous work by the group on the 3d transition metals.23 As these optimizations follow a multi-variate minimization scheme and the objective function landscape exhibits multiple valleys, the quality of the final output depends on the initial guess and is also optimization procedure dependent. In order to overcome this, we have used a number of different initial starting points as well as different sets of weights that enable the routine to explore a much larger part of the parameter space. This has provided well-optimized solutions that fulfil the imposed accuracy requirements.

The overall target for the plane wave energy cut-off was 400 Ry. We found this value to be close to the “sweet spot” with regard to the balance between accuracy vs gains in efficiency. Of course, this value serves only as a guiding parameter since the actual cut-offs in various codes will depend on calculated systems, accuracy criteria, etc.

All the ECP parameters are given in Table I. For all updated ECPs in this work, the core removed is the innermost ten electrons, referred to as a [Ne]-core for simplicity.

TABLE I.

Parameters for the ccECP-soft. For all ECPs, the highest value corresponds to the local channel L. Note that the highest non-local angular momentum channel max is related to it as max = L − 1.

AtomZeffnℓkαℓkβℓkAtomZeffnℓkαℓkβℓk
Cr 14 9.800 322 89.846 846 Mn 15 11.244 397 57.880 958 
8.010 010 18.997 257 11.614 251 92.965 750 
8.785 958 44.926 062 8.702 628 44.447 892 
7.014 726 14.003 861 14.217 018 41.889 380 
3.497 383 14.000 000 4.039 945 15.000 000 
3.611 831 48.963 362 4.200 000 60.599 175 
3.449 201 −56.466 431 4.139 297 −65.806 234 
2.009 794 0.968 440     
Fe 16 13.221 833 153.088 061 Co 17 11.423 427 90.855 286 
7.769 539 11.680 385 9.920 127 25.185 194 
9.100 629 40.685 923 9.811 352 48.270 556 
7.483 933 14.200 485 9.340 854 14.620 602 
3.798 917 16.000 000 3.932 921 17.000 000 
3.576 729 60.782 672 4.547 187 66.859 664 
3.514 698 −66.518 840 4.242 934 −76.154 505 
3.058 692 1.621 670 2.106 360 1.483 219 
Ni 18 10.199 961 41.053 383 Cu 19 12.068 348 78.019 159 
11.552 726 66.727 192 9.360 313 27.107 011 
8.131 870 24.281 961 13.173 488 54.905 280 
11.380 045 36.306 696 6.969 207 14.661 758 
3.641 646 18.000 000 3.806 452 19.000 000 
3.641 643 65.549 624 4.021 416 72.322 595 
3.637 271 −73.527 489 3.885 376 −84.688 200 
3.327 582 −0.856 416 2.626 437 3.393 685 
Zn 20 12.006 960 56.869 394       
9.103 589 34.859 484 
10.245 529 32.153 902 
7.286 335 15.898 530 
3.465 445 20.000 000 
3.528 420 69.308 902 
3.545 575 −83.673 652 
2.234 272 0.840 046 
AtomZeffnℓkαℓkβℓkAtomZeffnℓkαℓkβℓk
Cr 14 9.800 322 89.846 846 Mn 15 11.244 397 57.880 958 
8.010 010 18.997 257 11.614 251 92.965 750 
8.785 958 44.926 062 8.702 628 44.447 892 
7.014 726 14.003 861 14.217 018 41.889 380 
3.497 383 14.000 000 4.039 945 15.000 000 
3.611 831 48.963 362 4.200 000 60.599 175 
3.449 201 −56.466 431 4.139 297 −65.806 234 
2.009 794 0.968 440     
Fe 16 13.221 833 153.088 061 Co 17 11.423 427 90.855 286 
7.769 539 11.680 385 9.920 127 25.185 194 
9.100 629 40.685 923 9.811 352 48.270 556 
7.483 933 14.200 485 9.340 854 14.620 602 
3.798 917 16.000 000 3.932 921 17.000 000 
3.576 729 60.782 672 4.547 187 66.859 664 
3.514 698 −66.518 840 4.242 934 −76.154 505 
3.058 692 1.621 670 2.106 360 1.483 219 
Ni 18 10.199 961 41.053 383 Cu 19 12.068 348 78.019 159 
11.552 726 66.727 192 9.360 313 27.107 011 
8.131 870 24.281 961 13.173 488 54.905 280 
11.380 045 36.306 696 6.969 207 14.661 758 
3.641 646 18.000 000 3.806 452 19.000 000 
3.641 643 65.549 624 4.021 416 72.322 595 
3.637 271 −73.527 489 3.885 376 −84.688 200 
3.327 582 −0.856 416 2.626 437 3.393 685 
Zn 20 12.006 960 56.869 394       
9.103 589 34.859 484 
10.245 529 32.153 902 
7.286 335 15.898 530 
3.465 445 20.000 000 
3.528 420 69.308 902 
3.545 575 −83.673 652 
2.234 272 0.840 046 

The errors of the atomic spectrum for ccECP and ccECP-soft are evaluated by the mean absolute deviation (MAD) of considered atomic excitations that include bounded anions, sd transfers and cations up to third/fourth ionizations

MAD=1NiNΔEiECPΔEiAE.
(6)

For all the elements, the errors are provided at their complete basis limit (CBS). The CBS atomic state energies are estimated from the extrapolation method we used in our previous papers.18,23,25 The detailed data of the AE spectrum and ECP discrepancies for each atom can be found in the supplementary material.

In Fig. 1, we show the summary of spectral errors of ccECP and ccECP-soft compared to AE calculations. ccECP-soft constructions show mildly larger discrepancies compared to ccECPs, but the MAD errors are essentially within chemical accuracy for all the elements. All energies involved are in eV unless specified otherwise.

FIG. 1.

Mean absolute deviations (MADs) of atomic excitations for ccECP and ccECP-soft with the reference represented by the scalar relativistic RCCSD(T) method.

FIG. 1.

Mean absolute deviations (MADs) of atomic excitations for ccECP and ccECP-soft with the reference represented by the scalar relativistic RCCSD(T) method.

Close modal

Cr starts off the set and showcases some of the compromises made to soften the previously made ccECPs. The spectrum performance is of a similar quality for the chosen spectra sets as the standard ccECP but is decidedly worse at higher ionizations. In general, all of the ccECP-soft pseudopotentials struggle to capture the nature of highly ionized states and, as such, require truncated training sets compared to the standard ccECPs. Despite this, Cr does well to demonstrate the benefits that such compromises can lead to. The highest discrepancy observed in the spectral states was 0.042 eV for the [Ar]d6 state, confirming that all states tested lie within chemical accuracy for the spectrum. Similarly, the binding energy curves for CrH and CrO in Fig. 2 show that the discrepancies lie within chemical accuracy across all of the geometries investigated. Thus, with the minor compromise on the overall quality of the ccECP, we were able to specialize the ccECP-soft version with a significant decrease in plane wave cutoff to about 300 Ry.

FIG. 2.

Binding energy discrepancies for (a) CrH and (b) CrO molecules with the reference being the scalar relativistic all-electron CCSD(T) result. The shaded region indicates the band of chemical accuracy. The dashed vertical line represents the equilibrium geometry.

FIG. 2.

Binding energy discrepancies for (a) CrH and (b) CrO molecules with the reference being the scalar relativistic all-electron CCSD(T) result. The shaded region indicates the band of chemical accuracy. The dashed vertical line represents the equilibrium geometry.

Close modal

The new ccECP-soft construction for Mn has a slightly different form from the rest, possessing only three Gaussians in the local channel as opposed to four used for the rest of the ECPs in the series. However, the spectrum shows very good agreement with all electron values and the largest error of ≈0.06 eV for the 4F(3d7) state and very good balance for the rest of the states. Similarly, the molecular curves are within the chemical accuracy (Fig. 3).

FIG. 3.

Binding energy discrepancies for (a) MnH and (b) MnO molecules with the reference being the scalar relativistic all-electron CCSD(T) result. The shaded region indicates the band of chemical accuracy. The dashed vertical line represents the equilibrium geometry.

FIG. 3.

Binding energy discrepancies for (a) MnH and (b) MnO molecules with the reference being the scalar relativistic all-electron CCSD(T) result. The shaded region indicates the band of chemical accuracy. The dashed vertical line represents the equilibrium geometry.

Close modal

1. MnO2

We additionally tested the agreement between the ccECP-soft and the original ccECP for the Mn atom to ensure that softening the potential did not hamper its transferability in chemical environments that result in higher oxidation states. Manganese dioxide (MnO2) was picked as a representative for the 3d transition metals due to its position near the center of the series. We first picked a suitable state and geometry, and then ran a series of bond lengths for the oxygen atoms. We restricted the probe to the Hartree–Fock level since that is typically, as shown previously, the dominant source of overall bias. We chose a state with the highest spin possible, 2S = 3, and we specified the geometry to be unbent, θ = 180°, as can occur in some of the phases of related solid MnO2 material. Finally, we checked the differences in binding energies between ccECP and ccECP-soft at different geometries, as shown in Fig. 4 using the original ccECP as the reference. The original ccECP was optimized using an objective function ionizing the atom all the way to the semi-core electrons, as a result, its accuracy at high oxidation/ionization states is inherent to the ccECP construction.23 

FIG. 4.

HF binding energy discrepancy of Mn ccECP-soft for MnO2 molecule with reference being the ccECP from Annaberdiyev et al.23 The shaded region shows the band of chemical accuracy. The bond length indicates the distance between both oxygen atoms and the manganese atom.

FIG. 4.

HF binding energy discrepancy of Mn ccECP-soft for MnO2 molecule with reference being the ccECP from Annaberdiyev et al.23 The shaded region shows the band of chemical accuracy. The bond length indicates the distance between both oxygen atoms and the manganese atom.

Close modal

For the Fe atom, we achieved very convincing results, our ccECP-soft version is almost fully comparable to the original ccECP in overall accuracy. For comparison purposes, we also show molecular discrepancies from the DFT-derived set of alternative ECPs by Krogel–Santana–Reboredo, which are constructed to have very small cut-offs in general.26 In order to do so, we refitted the original form with Gaussians so that we could analyze its behavior and quantify the energy discrepancies using Gaussian-based codes. For the oxide dimer, we see noticeable overbinding bias, which increases with shortening of the bond length [Fig. 5(b)].

FIG. 5.

Binding energy discrepancies for (a) FeH and (b) FeO molecules with the reference being the scalar relativistic all-electron CCSD(T) result. KSR denotes the soft, DFT-derived ECP by Krogel, Santana, and Reboredo.26 The shaded region indicates the band of chemical accuracy. The dashed vertical line represents the equilibrium geometry.

FIG. 5.

Binding energy discrepancies for (a) FeH and (b) FeO molecules with the reference being the scalar relativistic all-electron CCSD(T) result. KSR denotes the soft, DFT-derived ECP by Krogel, Santana, and Reboredo.26 The shaded region indicates the band of chemical accuracy. The dashed vertical line represents the equilibrium geometry.

Close modal

Similar results as for Fe have been obtained for Co. The atomic spectrum that included ionizations up to the 5D(d6) state came out marginally worse than for the original ccECP construction. However, both molecular binding curves show excellent agreement with the coupled cluster (CC) calculations within the chemical accuracy range (Fig. 6). The reasonable plane wave energy cut-off is ∼360 Ry, with the obvious caveat that this might differ somewhat depending on the system, accuracy thresholds, and the used code.

FIG. 6.

Binding energy discrepancies for (a) CoH and (b) CoO molecules with the reference being the scalar relativistic all-electron CCSD(T) result. The shaded region indicates the band of chemical accuracy. The dashed vertical line represents the equilibrium geometry.

FIG. 6.

Binding energy discrepancies for (a) CoH and (b) CoO molecules with the reference being the scalar relativistic all-electron CCSD(T) result. The shaded region indicates the band of chemical accuracy. The dashed vertical line represents the equilibrium geometry.

Close modal

The results for the Ni atom show a similar pattern apart from a minor accuracy compromise. We observe that the hydride molecule bias is roughly constant along the whole binding curve. On the other hand, we clearly see small bias approaching chemical accuracy at the very shortest bond lengths for the oxide dimer [Fig. 7(b)]. We suspect that the fidelity tuning is complicated by the well-known 3d ↔ 4s degeneracy/instability. Note, however, that around the bond equilibrium we see excellent agreement with the all-electron reference and, therefore, we deem this minor accuracy compromise to be acceptable. Probing such small bond lengths would correspond to extremely high pressures in solid systems, such as NiO crystals. The plane wave cut-off is ∼375 Ry based on atomic criteria in Opium, although the corresponding value in solid state calculations might differ somewhat depending on the system.

FIG. 7.

Binding energy discrepancies for (a) NiH and (b) NiO molecules. The rest of the notations are the same as in the figures above.

FIG. 7.

Binding energy discrepancies for (a) NiH and (b) NiO molecules. The rest of the notations are the same as in the figures above.

Close modal

The quality of ccECP-soft constructions for Cu and Zn atoms is comparable. The requirement of smooth and as shallow as possible local potential is in contradiction with the large number of valence electrons and corresponding Zeff. As a consequence, the spectral errors are a bit larger but still acceptable. However, the molecular data show essentially the same accuracy as the original ccECP construction (Figs. 8 and 9). The approximate plane wave energy cut-off is about 400 Ry, which enables highly accurate calculations of solids and other systems in periodic settings. It is actually remarkable that even with rather deep semicore state 3s, one can achieve such a low cut-off within rather stringent accuracy requirements.

FIG. 8.

Binding energy discrepancies for (a) CuH and (b) CuO molecules. The rest of the notations are the same as in the figures above.

FIG. 8.

Binding energy discrepancies for (a) CuH and (b) CuO molecules. The rest of the notations are the same as in the figures above.

Close modal
FIG. 9.

Binding energy discrepancies for (a) ZnH and (b) ZnO molecules. The rest of the notations are the same as in the figures above.

FIG. 9.

Binding energy discrepancies for (a) ZnH and (b) ZnO molecules. The rest of the notations are the same as in the figures above.

Close modal

In Table II, we present the summary of molecular binding property discrepancies for ccECP and ccECP-soft. Overall, we see comparable quality in both ECPs with all-electron binding parameters. This can be expected since the binding energy discrepancy figures for both ECPs are balanced within chemical accuracy for all the molecules, with the minor exception mentioned above for NiO, where ccECP-soft mildly underbinds at very compressed bond lengths.

TABLE II.

Mean absolute deviations of binding parameters for various core approximations with respect to AE correlated data for related hydride and oxide molecules. All parameters were obtained using Morse potential fit. The parameters shown are dissociation energy De, equilibrium bond length re, vibrational frequency ωe, and binding energy discrepancy at dissociation for compressed bond length Ddiss.

De (eV)re (Å)ωe (cm−1)Ddiss (eV)
ccECP 0.0119(50) 0.0008(17) 2.8(5.5) 0.015(47) 
ccECP-soft 0.0174(50) 0.0012(17) 4.6(5.5) 0.023(47) 
De (eV)re (Å)ωe (cm−1)Ddiss (eV)
ccECP 0.0119(50) 0.0008(17) 2.8(5.5) 0.015(47) 
ccECP-soft 0.0174(50) 0.0012(17) 4.6(5.5) 0.023(47) 

To verify the performance of the ccECP-soft series in planewave codes, we ran several solid state calculations in QUANTUM ESPRESSO at various cutoffs to benchmark the softened ccECPs. We investigated FeO, CuO, and ZnO in the local-density approximation (LDA)+U framework as a representative set for the series to determine the kinetic energy (KE) cutoffs for the ccECP-softs produced in this work. Each system has a neutral unit cell, for FeO, the unit cell is rhombohedral with two iron and two oxygen atoms, and the solid is in an anti-ferromagnetic phase. The unit cell of CuO studied here is in C2/c symmetry,27 with four copper and four oxygen atoms; the solid is in the non-magnetic phase. Finally, the unit cell of ZnO has a P63mc structure,28 with two zinc and two oxygen atoms. The solid is in the non-magnetic phase. First, a reference calculation is run at a KE cutoff of 2000 Ry, at which point the total energy should not meaningfully change if the cutoff is increased further, and beyond this value would become needlessly expensive to compute. Then, a range of cutoffs are run from 300 to 500 Ry in increments of 10 Ry to scan the discrepancy from the reference energy. The system is deemed “converged” with regard to cutoff when the energy discrepancy per valence electron falls below 1 meV. Further discussion of atomic and molecular cutoff can be found in the supplementary material. Finally, we repeated the process described above for the original ccECP parameterizations from our previous work23 and found their cutoffs in a similar way for comparison purposes. The major difference was the amount of time each sub-calculation ran and the cutoff at which the parameterizations converged. For the ccECP-softs, the cutoffs for the transition metal oxides (TMOs) we tested converged around 400 Ry, whereas for the ccECPs, they converged well after 1000 Ry. CuO even has a cutoff of 1500 Ry as seen in Table III. The graphs from which the table values were taken are included in the supplementary material.

TABLE III.

Comparison of cutoffs for several TMOs in the series using the newly created ccECP-soft and ccECP-pseudopotentials. The original ccECPs generally triple the cutoff of the soft versions, leading to a much higher computational cost.

TMOccECP-soft (Ry)ccECP (Ry)
FeO 430 1140 
CuO 390 1500 
ZnO 390 1315 
TMOccECP-soft (Ry)ccECP (Ry)
FeO 430 1140 
CuO 390 1500 
ZnO 390 1315 

We present a new modified set of ccECP-soft for late 3d elements Mn–Zn with [Ne]-core that can be used with plane wave codes. The two key goals we address are accuracy through many-body construction and, at the same time, efficient application to plane wave codes with low energy cut-offs. For these purposes, we use mildly less stringent criteria on spectral fidelity and on tests of molecular binding. We do not include highly ionized states (beyond 6+) into the optimization and we allow smaller exponents in the local channel that proved to be the key parameters for this purpose. This leads to smaller curvatures and more shallow local potentials in the core region.

We obtained very encouraging results. In particular, for Cr–Ni atoms, the spectral errors are small and almost on par with our more stringent original construction. Mildly less accurate spectra for Cu and Zn atoms are, however, still qualified for high accuracy calculations. With all the molecular binding energy discrepancy curves for the hydride and oxide molecules are strictly within the chemical accuracy, showing excellent fidelity also in bonding situations away from equilibrium. Overall, the newly constructed and modified set of ccECP-soft pseudopotentials are basically in the same accuracy class as the original ccECPs. We also expect essentially all the properties, such as correlation energies, calculated previously29 to follow the same pattern with very minor differences, say, at few mHa level or so.

We believe that the constructed set will significantly expand the usefulness of ccECPs and provide tested and consistent data for accurate valence-only electronic structure calculations in many-body approaches.

More figures and tables relating to the benchmarks of the ccECP-soft ECPs are included in the supplementary material for this paper. There are tables showing a comparison of the ccECP and ccECP-soft energy discrepancies in the atomic spectra for each of the elements, a full table of the cutoffs for each new ccECP-soft, and convergence graphs for the cutoffs displayed in Table III.

We are grateful to Paul R. C. Kent and Luke Shulenburger for reading the manuscript and helpful suggestions.

This work was supported by the U.S. Department of Energy, Office of Science, Basic Energy Sciences, Materials Sciences and Engineering Division, as part of the Computational Materials Sciences Program and Center for Predictive Simulation of Functional Materials.

This research used resources of the National Energy Research Scientific Computing Center (NERSC), a U.S. Department of Energy Office of Science User Facility operated under Contract No. DE-AC02-05CH11231. This research used resources of the Argonne Leadership Computing Facility, which is a DOE Office of Science User Facility supported under Contract No. DE-AC02-06CH11357. This research also used resources of the Oak Ridge Leadership Computing Facility, which is a DOE Office of Science User Facility supported under Contract No. DE-AC05-00OR22725.

The authors have no conflicts to disclose.

B.K. and G.W. contributed equally to this work.

Benjamin Kincaid: Conceptualization (equal); Data curation (equal); Investigation (equal); Methodology (equal); Writing – original draft (equal); Writing – review & editing (equal). Guangming Wang: Conceptualization (equal); Investigation (equal); Methodology (equal); Writing – original draft (equal); Writing – review & editing (equal). Haihan Zhou: Investigation (equal); Writing – original draft (equal); Writing – review & editing (equal). Lubos Mitas: Formal analysis (equal); Funding acquisition (equal); Methodology (equal); Project administration (equal); Supervision (equal); Writing – original draft (equal).

All ccECP-soft in various code formats and corresponding basis sets are available at https://pseudopotentiallibrary.org.

1.
L.
Fernandez Pacios
and
P. A.
Christiansen
,
J. Chem. Phys.
82
,
2664
(
1985
).
2.
M. M.
Hurley
,
L. F.
Pacios
,
P. A.
Christiansen
,
R. B.
Ross
, and
W. C.
Ermler
,
J. Chem. Phys.
84
,
6840
(
1986
).
3.
L. A.
LaJohn
,
P. A.
Christiansen
,
R. B.
Ross
,
T.
Atashroo
, and
W. C.
Ermler
,
J. Chem. Phys.
87
,
2812
(
1987
).
4.
R. B.
Ross
,
J. M.
Powers
,
T.
Atashroo
,
W. C.
Ermler
,
L. A.
LaJohn
, and
P. A.
Christiansen
,
J. Chem. Phys.
93
,
6654
(
1990
).
5.
W. J.
Stevens
,
M.
Krauss
,
H.
Basch
, and
P. G.
Jasien
,
Can. J. Chem.
70
,
612
(
1992
).
6.
N.
Troullier
and
J. L.
Martins
,
Phys. Rev. B
43
,
1993
(
1991
).
7.
M.
Burkatzki
,
C.
Filippi
, and
M.
Dolg
,
J. Chem. Phys.
126
,
234105
(
2007
).
8.
M.
Burkatzki
,
C.
Filippi
, and
M.
Dolg
,
J. Chem. Phys.
129
,
164115
(
2008
).
9.
A.
Bergner
,
M.
Dolg
,
W.
Küchle
,
H.
Stoll
, and
H.
Preuß
,
Mol. Phys.
80
,
1431
(
1993
).
10.
M.
Dolg
,
U.
Wedig
,
H.
Stoll
, and
H.
Preuss
,
J. Chem. Phys.
86
,
866
(
1987
).
11.
G. B.
Bachelet
,
D. R.
Hamann
, and
M.
Schlüter
,
Phys. Rev. B
26
,
4199
(
1982
).
12.
J. R.
Trail
and
R. J.
Needs
,
J. Chem. Phys.
146
,
204107
(
2017
).
13.
D.
Vanderbilt
,
Phys. Rev. B
41
,
7892
(
1990
).
14.
M.
Dolg
and
X.
Cao
,
Chem. Rev.
112
,
403
(
2012
).
15.
C.
Hartwigsen
,
S.
Goedecker
, and
J.
Hutter
,
Phys. Rev. B
58
,
3461
(
1998
).
16.
A.
Willand
,
Y. O.
Kvashnin
,
L.
Genovese
,
Á.
Vázquez-Mayagoitia
,
A. K.
Deb
,
A.
Sadeghi
,
T.
Deutsch
, and
S.
Goedecker
,
J. Chem. Phys.
138
,
104109
(
2013
).
17.
M. C.
Bennett
,
C. A.
Melton
,
A.
Annaberdiyev
,
G.
Wang
,
L.
Shulenburger
, and
L.
Mitas
,
J. Chem. Phys.
147
,
224106
(
2017
).
18.
G.
Wang
,
A.
Annaberdiyev
,
C. A.
Melton
,
M. C.
Bennett
,
L.
Shulenburger
, and
L.
Mitas
,
J. Chem. Phys.
151
,
144110
(
2019
).
19.
A.
Annaberdiyev
,
G.
Wang
,
C. A.
Melton
,
M. C.
Bennett
, and
L.
Mitas
,
Phys. Rev. B
103
,
205206
(
2021
).
20.
X.
Zhou
,
H.
Zhao
,
T.
Wang
, and
F.
Wang
,
Phys. Rev. A
100
,
062502
(
2019
).
21.
S. G.
Louie
,
S.
Froyen
, and
M. L.
Cohen
,
Phys. Rev. B
26
,
1738
(
1982
).
22.
J.
Kolorenč
and
L.
Mitas
,
Phys. Rev. B
75
,
235118
(
2007
).
23.
A.
Annaberdiyev
,
G.
Wang
,
C. A.
Melton
,
M.
Chandler Bennett
,
L.
Shulenburger
, and
L.
Mitas
,
J. Chem. Phys.
149
,
134108
(
2018
).
24.
P.
Spellucci
, Donlp2 nonlinear optimization code, https://www2.mathematik.tu-darmstadt.de/fbereiche/numerik/staff/spellucci/DONLP2/,
2009
.
25.
M. C.
Bennett
,
G.
Wang
,
A.
Annaberdiyev
,
C. A.
Melton
,
L.
Shulenburger
, and
L.
Mitas
,
J. Chem. Phys.
149
,
104108
(
2018
).
26.
J. T.
Krogel
,
J. A.
Santana
, and
F. A.
Reboredo
,
Phys. Rev. B
93
,
075143
(
2016
).
27.
S.
Åsbrink
and
L. J.
Norrby
,
Acta Crystallogr., Sect. B: Struct. Crystallogr. Cryst. Chem.
26
,
8
(
1970
).
28.
Y.-N.
Xu
and
W. Y.
Ching
,
Phys. Rev. B
48
,
4335
(
1993
).
29.
A.
Annaberdiyev
,
C. A.
Melton
,
M. C.
Bennett
,
G.
Wang
, and
L.
Mitas
,
J. Chem. Theory Comput.
16
,
1482
(
2020
).

Supplementary Material